Issue triage: docstring cleanup, test fixes, gvec sync

docs/advanced/index.md

@@ -74,4 +74,5 @@ curved-boundary-conditions
 mesh-adaptation
 troubleshooting
 api-patterns
+SWARM-INTEGRATION-STATISTICS
 ```

docs/api/constitutive_models.md

@@ -1,10 +1,5 @@
 # Constitutive Models
 
-```{eval-rst}
-.. automodule:: underworld3.constitutive_models
-   :no-members:
-```
-
 ## Viscous Flow Models
 
 ### ViscousFlowModel

docs/api/coordinates.md

@@ -1,10 +1,5 @@
 # Coordinate Systems
 
-```{eval-rst}
-.. automodule:: underworld3.coordinates
-   :no-members:
-```
-
 ## CoordinateSystem
 
 The base class for coordinate system handling on meshes.

docs/api/discretisation.md

@@ -1,10 +1,5 @@
 # Discretisation
 
-```{eval-rst}
-.. automodule:: underworld3.discretisation
-   :no-members:
-```
-
 ## Mesh
 
 The computational mesh class that provides the spatial discretisation for finite element calculations.

docs/api/function.md

@@ -1,16 +1,11 @@
 # Function and Expressions
 
-```{eval-rst}
-.. automodule:: underworld3.function
-   :no-members:
-```
-
 ## Expressions
 
 ### UWexpression
 
 ```{eval-rst}
-.. autoclass:: underworld3.function.UWexpression
+.. autoclass:: underworld3.function.expressions.UWexpression
    :members:
    :show-inheritance:
 ```
@@ -75,12 +70,6 @@ Factory function for creating UWQuantity objects with units.
 .. autofunction:: underworld3.function.global_evaluate
 ```
 
-### evalf
-
-```{eval-rst}
-.. autofunction:: underworld3.function.evalf
-```
-
 ### evaluate_gradient
 
 ```{eval-rst}

docs/api/maths.md

@@ -1,11 +1,5 @@
 # Mathematical Operations
 
-```{eval-rst}
-.. automodule:: underworld3.maths
-   :members:
-   :show-inheritance:
-```
-
 ## Integral Operations
 
 ```{eval-rst}

docs/api/solvers.md

@@ -1,10 +1,5 @@
 # Solvers
 
-```{eval-rst}
-.. automodule:: underworld3.systems.solvers
-   :no-members:
-```
-
 ## Stokes Flow
 
 ### SNES_Stokes

docs/api/swarm.md

@@ -1,10 +1,5 @@
 # Particle Swarms
 
-```{eval-rst}
-.. automodule:: underworld3.swarm
-   :no-members:
-```
-
 ## Swarm
 
 ```{eval-rst}
@@ -21,14 +16,6 @@
    :show-inheritance:
 ```
 
-## Population Control
-
-```{eval-rst}
-.. autoclass:: underworld3.swarm.PopulationControl
-   :members:
-   :show-inheritance:
-```
-
 ## Index Swarm Variable
 
 ```{eval-rst}

docs/api/systems_ddt.md

@@ -1,10 +1,5 @@
 # Time Derivatives
 
-```{eval-rst}
-.. automodule:: underworld3.systems.ddt
-   :no-members:
-```
-
 Time derivative operators approximate $D\phi/Dt$ or $DF/Dt$ for transient
 solvers.  All operators share a common interface: ``update_pre_solve(dt)``
 before each timestep, ``bdf()`` for the BDF approximation in the weak form,

docs/beginner/quickstart.md

@@ -67,7 +67,7 @@ $$
 
 where $t(N_{ref})$ is the run-time for a reference number of processors, $N_{ref}$, and $t(N)$ is the run-time when $N$ processors are used. In the ideal case, $N$ additional processors should contribute all of its resources in solving the problem and reduce the compute time by a factor of $N$ relative to the reference run time. For example, using $2 N_{ref}$ processors will ideally halve the run-time resulting to a speed-up = 2.
 
-```{figure} media/UW3-StrongScalingSolvers.png
+```{figure} ../media/UW3-StrongScalingSolvers.png
 :name: fig-strong-scaling
 
 Strong parallel-scaling tests run on Australia's peak computing system, [GADI, at the National Computational Infrastructure](https://nci.org.au/our-systems/hpc-systems?ref=underworldcode.org). This is a typical High Performance Computing facility with large numbers of dedicated, identical CPUs and fast communication links.

docs/developer/UW3_Developers_MathematicalObjects.md

@@ -574,8 +574,7 @@ def test_mathematical_operations():
 ## Planned Enhancements
 
 ### 1. Operator Overloading
-```python
-
+```text
 # Future possibilities
 laplacian = velocity @ nabla**2  # Operator syntax
 curl = nabla × velocity          # Cross product notation

docs/developer/container/README.md

@@ -1,3 +1,7 @@
+---
+orphan: true
+---
+
 The `Containerfile` is used to build the container and can be used with docker or podman. 
 At present only amd64 architecture is built, because vtk-osmesa isn't available for arm by default. This may change in future, or we could build vtk-osmesa from source (see https://docs.pyvista.org/extras/building_vtk.html)
 

docs/developer/design/COORDINATE_MIGRATION_GUIDE.md

@@ -7,18 +7,18 @@
 
 ### Pattern Migration Cheat Sheet
 
-```python
-# ❌ DEPRECATED          →  ✅ RECOMMENDED
-mesh.data                 →  mesh.X.coords  or  mesh.points
-mesh.data.shape           →  mesh.X.coords.shape
-mesh.data[:, 0]           →  mesh.X.coords[:, 0]
-mesh.data.min()           →  mesh.X.coords.min()
-mesh.data.copy()          →  mesh.X.coords.copy()
-
-# ✅ UNCHANGED (these work exactly as before)
-mesh.X[0]                 →  mesh.X[0]  (symbolic x-coordinate)
-x, y = mesh.X             →  x, y = mesh.X  (symbolic unpacking)
-mesh.units                →  mesh.units  (coordinate units)
+```text
+DEPRECATED                   RECOMMENDED
+mesh.data                 -> mesh.X.coords  or  mesh.points
+mesh.data.shape           -> mesh.X.coords.shape
+mesh.data[:, 0]           -> mesh.X.coords[:, 0]
+mesh.data.min()           -> mesh.X.coords.min()
+mesh.data.copy()          -> mesh.X.coords.copy()
+
+UNCHANGED (these work exactly as before)
+mesh.X[0]                 -> mesh.X[0]  (symbolic x-coordinate)
+x, y = mesh.X             -> x, y = mesh.X  (symbolic unpacking)
+mesh.units                -> mesh.units  (coordinate units)
 ```
 
 ## Migration Examples by Use Case

docs/developer/design/WHY_UNITS_NOT_DIMENSIONALITY.md

@@ -20,7 +20,7 @@
 - **Dimensionality**: ⚠️ Abstract physics concept ([length], [time], [mass])
 
 ### What the API Shows:
-```python
+```text
 # User-friendly - clear what's happening
 var.units = "m"
 var.to("km")
@@ -345,7 +345,7 @@ from underworld3.utilities.units_mixin import UnitAwareMixin
 - What does that even mean?
 
 ### 3. API Clarity ✅
-```python
+```text
 # Clear and intuitive
 var.units = "m"
 var.to("km")

docs/developer/guides/HOW-TO-WRITE-UW3-SCRIPTS.md

@@ -612,7 +612,7 @@ def test_swarm_functionality():
 
 If you see generated C code with symbolic expressions instead of numbers:
 
-```c
+```text
 // ERROR symptom in generated code:
 out[0] = 1.0/{ \eta \hspace{ 0.0006pt } };  // Should be numeric!
 ```

docs/developer/subsystems/constitutive-models-anisotropy.md

@@ -45,7 +45,7 @@ In `underworld` tensor rotation is provided for rank 2 and rank 4 tensors by `uw
 
 If we apply constraints to the deformation, we expect to reduce the number of independent material constants. Incompressibility should reduce the number of *independent* material constants. In an isotropic medium (or a medium with cubic symmetry), incompressibility eliminates the one volumetric material modulus (e.g. the bulk modulus). In general anisotropic media, it is not the case that changes in pressure result in uniform expansion or contraction, and an incompressibility constraint reduces the number of *independent* material constants. In the transversely isotropic case, there are five independent materia constants in general, reducing to four when the material is incompressible. 
 
-It is not a given that the stiffness matrix is trivial to construct / meaningful for incompressible anisotropy and there is some discussion here: https://rastgaragah.wordpress.com/2013/03/12/incompressibility-of-linearly-elastic-material/ (identifies the issue) and this is explained in more detail by {cite}`destradeIncompressibleLimitLinear2002`
+It is not a given that the stiffness matrix is trivial to construct / meaningful for incompressible anisotropy and there is some discussion here: https://rastgaragah.wordpress.com/2013/03/12/incompressibility-of-linearly-elastic-material/ (identifies the issue) and this is explained in more detail by Destrade et al. (2002)
 
 
 ### Example
@@ -192,7 +192,7 @@ display(C_IJm_MM_iso)
 
 ## Han & Wahr, 1997 (full transverse isotropic tensor)
 
-In the {cite}`hanAnalysisAnisotropicMantle1997` paper, the expression for incompressible transverse-isotropy is as follows
+In the Han & Wahr (1997) paper, the expression for incompressible transverse-isotropy is as follows
 
 $$\left[\begin{matrix}2 \eta_{0} + \mu_{0} & \mu_{0} & 0 & 0 & 0 & 0\\\mu_{0} & 2 \eta_{0} + \mu_{0} & 0 & 0 & 0 & 0\\0 & 0 & - 2 \Delta\eta + 2 \eta_{0} + \mu_{1} & 0 & 0 & 0\\0 & 0 & 0 & - 2 \Delta\eta + 2 \eta_{0} & 0 & 0\\0 & 0 & 0 & 0 & - 2 \Delta\eta + 2 \eta_{0} & 0\\0 & 0 & 0 & 0 & 0 & 2 \eta_{0}\end{matrix}\right]$$
 
@@ -249,7 +249,7 @@ display(C_IJm_HW_R_s2)
 
 ## Orthotropic medium
 
-**Note** all the caveats above regarding incompressibility. The {cite}`browaeysDecompositionElasticTensor2004` elastic tensors have a bulk modulus term, so it is not completely obvious how to square the assumptions in the first two implementations with this set. 
+**Note** all the caveats above regarding incompressibility. The Browaeys & Chevrot (2004) elastic tensors have a bulk modulus term, so it is not completely obvious how to square the assumptions in the first two implementations with this set. 
 
 The full formulation should look like this:
 

src/underworld3/discretisation/__init__.py

@@ -4,20 +4,13 @@
 This module provides the core discretisation infrastructure for finite
 element computations in Underworld3.
 
-Classes
--------
-Mesh : class
-    Unstructured mesh with PETSc DMPlex backend. Supports various cell
-    types and coordinate systems.
-MeshVariable : class
-    Field variable defined on mesh (nodal or cell-based storage).
+**Mesh** -- Unstructured mesh with PETSc DMPlex backend.
 
-Functions
----------
-checkpoint_xdmf : function
-    Save mesh and variables to XDMF format for visualization.
-meshVariable_lookup_by_symbol : function
-    Find mesh variable by its symbolic representation.
+**MeshVariable** -- Field variable defined on mesh (nodal or cell-based storage).
+
+**checkpoint_xdmf** -- Save mesh and variables to XDMF format.
+
+**meshVariable_lookup_by_symbol** -- Find mesh variable by its symbol.
 
 See Also
 --------

src/underworld3/discretisation/discretisation_mesh.py

@@ -1869,11 +1869,11 @@ def write_checkpoint(
         index: Optional[int] = 0,
         unique_id: Optional[bool] = False,
     ):
-        """Write data in a format that can be restored for restarting the simulation
+        """Write data in a format that can be restored for restarting the simulation.
+
         The difference between this and the visualisation is 1) the parallel section needs
         to be stored to reload the data correctly, and 2) the visualisation information (vertex form of fields)
-        is not stored. This routines uses dmplex *VectorView and *VectorLoad functionality
-
+        is not stored. This routine uses dmplex VectorView and VectorLoad functionality.
         """
 
         # The mesh checkpoint is the same as the one required for visualisation

src/underworld3/function/__init__.py

@@ -5,21 +5,16 @@
 on meshes and swarms. It bridges SymPy symbolic mathematics with PETSc
 numerical evaluation.
 
-Key Components
---------------
-expression : class
-    User-facing symbolic expression wrapper (UWexpression).
-evaluate : function
-    Evaluate expressions at mesh/swarm points.
-global_evaluate : function
-    Parallel-safe evaluation gathering results across MPI ranks.
-UnderworldFunction : class
-    Core function evaluation machinery.
-
-Unit Conversion
----------------
-convert_quantity_units, make_dimensionless, add_units
-    Utilities for handling physical units in expressions.
+**expression** -- User-facing symbolic expression wrapper (UWexpression).
+
+**evaluate** -- Evaluate expressions at mesh/swarm points.
+
+**global_evaluate** -- Parallel-safe evaluation gathering results across MPI ranks.
+
+**UnderworldFunction** -- Core function evaluation machinery.
+
+**convert_quantity_units**, **make_dimensionless**, **add_units** --
+Utilities for handling physical units in expressions.
 
 The function module integrates with pint for optional unit-aware
 computations when the mesh has an associated unit registry.

src/underworld3/function/_function.pyx

@@ -133,18 +133,16 @@ class UnderworldAppliedFunctionDeriv(UnderworldAppliedFunction):
 
 class UnderworldFunction(sympy.Function):
     """
-    This is a metaclass, so it returns programmatic class objects rather
-    than instances. This basically follows the pattern of the `sympy.Function`
-    metaclass, with two key differences:
-    1. We set `UnderworldAppliedFunction` as the base class. This is really just a
-       dummy class (see its definition) which allows us to do calls such
-       as `isinstance(someobj, UnderworldAppliedFunction)` to test if a `sympy`
-       object is one we've defined.
-    2. We grab a weakref of the owning meshvariable onto the *class* itself. Note
-       that it's important that it's recorded onto the class (instead of the instance),
-       as Sympy internally sometimes uses calls such as `type(obj)(obj.args)` to
-       replace objects with cloned instances, and therefore 'extra' info must be
-       recorded onto the class so that the clones are _complete_.
+    Metaclass that returns programmatic class objects rather than instances.
+
+    This basically follows the pattern of the ``sympy.Function``
+    metaclass, with two key differences. First, we set
+    ``UnderworldAppliedFunction`` as the base class, which allows
+    ``isinstance(someobj, UnderworldAppliedFunction)`` checks.
+    Second, we grab a weakref of the owning meshvariable onto the
+    class itself (not the instance), because SymPy internally uses
+    ``type(obj)(obj.args)`` to clone instances and extra info must
+    be on the class so that clones are complete.
 
     Consider the calling pattern
 
@@ -157,16 +155,16 @@ class UnderworldFunction(sympy.Function):
 
     Parameters
     ----------
-    name:
+    name : str
         The name of the function.
-    meshvar:
+    meshvar : MeshVariable
         The mesh variable corresponding to this function.
-    vtype:
-        The variable type (scalar,vector,etc).
-    component:
+    vtype : VarType
+        The variable type (scalar, vector, etc).
+    component : int or tuple
         For vector functions, this is the component of the vector.
-        For example, component `1` might correspond to `v_y`.
-        For tensors, the component is a tuple
+        For example, component ``1`` might correspond to ``v_y``.
+        For tensors, the component is a tuple.
         For scalars, this value is ignored.
     """
     def __new__(cls,