Route constant solver parameters through PetscDS constants[] array

[lmoresi]

Summary

• Constant UWexpression solver parameters (viscosity, diffusivity, dt_elastic, BC amplitudes) are now automatically routed through PETSc's constants[] array instead of being baked as C literals
• Changing a constant value between time steps no longer triggers JIT recompilation — only PetscDSSetConstants() is called (~0ms vs ~10s)
• Fixed pre-existing bug: _build() re-added Null_Boundary every call (string membership check on list of namedtuples), forcing unnecessary recompilation for any solver with natural BCs

How it works

1. _extract_constants() scans expression trees for UWexpression atoms that resolve to pure numbers
2. These are replaced with _JITConstant symbols that render as constants[i] in generated C code
3. Cache hash uses structural form (constants as placeholders) — value changes don't cause cache misses
4. _update_constants() re-packs values and calls PetscDSSetConstants() before every snes.solve()
5. Detection works at any nesting depth (user expression -> constitutive model parameter -> solver template)

Test plan

• All level_1 tests pass (372 passed, 10 skipped, 1 xfailed)
• All level_2 tests pass (425 passed, 42 skipped, 14 xfailed)
• New functional tests: Poisson constant diffusivity, change-without-recompile, constant source, Stokes viscosity, essential BC amplitude, Stokes natural BC traction
• VE Maxwell shear box benchmark produces correct results (verified on combined branch)
• CI pipeline

Underworld development team with AI support from Claude Code

[Copilot]

Pull request overview

This PR updates Underworld3’s JIT compilation pipeline to route spatially-constant UWexpression parameters through PETSc’s constants[] array (via PetscDSSetConstants()), avoiding expensive recompilation when only parameter values change between solves.

Changes:

• Add constant extraction/substitution so constant UWexpression atoms compile as constants[i] (structural hashing becomes value-independent).
• Update solvers to set/update PETSc DS constants before each solve (including MG coarse levels) and fix Null_Boundary duplication in natural BC setup.
• Add functional tests and expand developer/user docs explaining the new constants mechanism and performance implications.

Reviewed changes

Copilot reviewed 10 out of 10 changed files in this pull request and generated 3 comments.

Show a summary per file

File	Description
tests/test_1001_poisson_constants.py	New functional tests covering Poisson/Stokes constants, no-recompile behavior, and the Null_Boundary fix.
tests/test_0004_pointwise_fns.py	Update call sites for new getext() return type.
src/underworld3/utilities/_jitextension.py	Implement _JITConstant, constant extraction/packing, structural hashing, and return constants_manifest from getext().
src/underworld3/cython/petsc_maths.pyx	Update getext() usage to handle namedtuple return.
src/underworld3/cython/petsc_generic_snes_solvers.pyx	Add constants propagation via PetscDSSetConstants(), call updates before solve, and fix Null_Boundary membership check.
src/underworld3/cython/petsc_extras.pxi	Declare PetscDSSetConstants() for Cython use.
docs/developer/subsystems/expressions-functions.md	Document expressions + the constants mechanism and how it affects JIT performance.
docs/developer/guides/HOW-TO-WRITE-UW3-SCRIPTS.md	Add troubleshooting guidance for slow recompiling time-stepping loops and constant unwrapping symptoms.
docs/developer/UW3_Developers_MathematicalObjects.md	Document two-phase unwrap and constants-based updates without recompilation.
docs/beginner/parameters.md	Beginner-facing guidance to wrap solver parameters with uw.expression() for efficient updates.

---

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

[Copilot]

_extract_constants() imports is_constant_expr and extract_expressions but neither is used, and _collect_constant_atoms() takes is_constant_expr as an argument without using it. This adds dead code/indirection and makes it harder to understand the constant detection logic. Consider removing the unused imports and the unused parameter from _collect_constant_atoms() (and its call sites).

[Copilot]

Constants are sorted and named using str(expr), but UWexpression.__str__ is based on the mutable display name (see UWexpression.rename()), not the symbol identity. This can make constants ordering (and thus constants[i] mapping + cache hash) change due to display-only renames, leading to unnecessary recompiles or mismatched packing. Prefer using the stable expr.name (and/or _uw_id where applicable) for deterministic ordering and placeholder naming.

Suggested change

	# Sort by name for deterministic MPI-consistent ordering
	sorted_constants = sorted(constant_exprs, key=lambda e: str(e))
	manifest = []
	subs_map = {}
	for i, expr in enumerate(sorted_constants):
	jit_const = _JITConstant(i, name=f"_jit_const_{str(expr)}")
	# Sort by a stable identifier (prefer UWexpression.name) for deterministic MPI-consistent ordering
	sorted_constants = sorted(constant_exprs, key=lambda e: getattr(e, "name", str(e)))
	manifest = []
	subs_map = {}
	for i, expr in enumerate(sorted_constants):
	base_name = getattr(expr, "name", str(expr))
	jit_const = _JITConstant(i, name=f"_jit_const_{base_name}")

[Copilot]

_pack_constants() silently falls back to 0.0 if a constant cannot be converted to float, which can mask real configuration/errors and produce incorrect physics with no signal. It would be safer to raise a clear exception (or at least warn and include the expression name/value) so callers know constants could not be packed.

Suggested change

	except (TypeError, ValueError):
	# Fallback: try .data property
	try:
	values[idx] = float(uw_expr.data)
	except Exception:
	values[idx] = 0.0
	except (TypeError, ValueError) as exc:
	# Fallback: try .data property
	try:
	values[idx] = float(getattr(uw_expr, "data"))
	except Exception as data_exc:
	raise RuntimeError(
	"Failed to pack JIT constant at index "
	f"{idx}: could not convert UWexpression {uw_expr!r} to float. "
	"Both unwrap_expression(..., mode='nondimensional') and the "
	".data attribute conversion failed."
	) from data_exc