PETSc branch review

Reviewed jose/ksp-multishift at 3a6ae79d28d129b5c0eb68716e1ef891effcdd70 against origin/main: 21 files changed, 1,800 insertions, 17 deletions.

Four MEDIUM findings and three Style findings.

Findings

MEDIUM — Unpreconditioned convergence uses an M⁻¹-scaled residual

Locations: src/ksp/ksp/impls/eksm/eksm.c:276–303, :177, :205, :819–820.

EKSM normalizes M⁻¹b and passes that norm to the convergence test. Its projected residuals likewise correspond to the transformed system, even when the user selects KSP_NORM_UNPRECONDITIONED.

For a scalar problem with K=1, M=10¹², sigma=0, b=1, exact inner solves, and atol=10⁻⁸, the initial reported residual is 10⁻¹². EKSM therefore declares convergence at iteration zero and returns x=0, although the unpreconditioned residual is 1.

Keep the transformed norm for the internal algorithm, but supply convergence and monitoring norms consistent with the selected KSPNormType. In particular, the advertised unpreconditioned mode must measure the original system’s residual.

MEDIUM — Matrix mutations leave multishift metadata stale

Locations: src/mat/impls/nest/multishift.c:140–146, :179–180; src/ksp/ksp/impls/eksm/eksm.c:49–65, :273–275.

The factory attaches the original K, M, and shifts to an ordinary MATNEST. Subsequent supported matrix operations do not update or invalidate that description, while EKSM reconstructs its operator exclusively from it.

For example, create an explicit nest with K=1, M=NULL, and sigma=0, then call MatScale(A,2) before solving with b=1. MatScale_Nest() changes the diagonal block to 2, but the stored K remains 1. EKSM solves the original system and returns 1 instead of ½.

The metadata needs a validity check against the matrix/component states, or maintenance through supported mutations. The current presence-only check must not accept a description that no longer represents A.

MEDIUM — Solution and residual callbacks see the zero initial guess

Locations: src/ksp/ksp/impls/eksm/eksm.c:375–392, :829–834.

The solution is assembled only after the iteration loop, and EKSM retains KSPBuildSolutionDefault(). Consequently, calls to KSPBuildSolution() or KSPBuildResidual() during convergence tests and monitors use the zero initial guess rather than the current approximation.

For K=1, M=NULL, sigma=0, and b=1, the first iteration has the exact solution, but the final in-loop true-residual monitor still sees residual 1. A custom convergence test using KSPBuildResidual() can reject that exact solution and cause KSP_DIVERGED_BREAKDOWN.

Provide an EKSM solution-building operation using the current basis and projected coefficients. Also ensure the final solution is materialized before callbacks for which ksp->reason is already nonzero, since KSPBuildSolution() then bypasses the implementation-specific operation.

MEDIUM — Zero happy-breakdown tolerance misses exact breakdown

Locations: src/ksp/ksp/impls/eksm/eksm.c:335, :362.

KSPEKSMSetHapTol() accepts zero, but both breakdown checks use strict <. An exactly zero new direction therefore fails the check.

With K=1, M=NULL, sigma=1, b=1, and haptol=0, the first direction is exactly dependent. Continuing produces the projected matrix

[1  2]
[0  0]
[0  0]

A standalone call to system LAPACK confirmed DGELS INFO=2, which PetscCallLAPACKInfo() turns into an error. Taking the exact-breakdown branch instead gives x=½ and zero residual.

Suggested patch, not applied or tested in PETSc:

--- a/src/ksp/ksp/impls/eksm/eksm.c
+++ b/src/ksp/ksp/impls/eksm/eksm.c
@@
-    if (nrm < eksm->haptol) {
+    if (nrm <= eksm->haptol) {
@@
-      if (nrm < eksm->haptol) {
+      if (nrm <= eksm->haptol) {

Style — KSP consumes Mat’s private metadata layout directly

Locations: src/ksp/ksp/impls/eksm/eksmimpl.h:7; src/ksp/ksp/impls/eksm/eksm.c:50–55, :274–275.

EKSM includes the nest implementation header and extracts Mat_MultiShift from a composed container. This conflicts with AGENTS.md, Core Working Rules: use PETSc accessors for object metadata and extend the owning API when metadata access is missing.

Provide a Mat-owned accessor for the matrices, shifts, and complex-pair information, documenting array lengths, ownership, and lifetime. Use it from EKSM’s setup and solve paths.

Style — Macro arguments need parentheses

Locations: src/ksp/ksp/impls/eksm/eksmimpl.h:49; src/mat/impls/nest/matnestimpl.h:39–40.

These arguments appear next to operators without parentheses, contrary to doc/developers/style.md, C Usage, item 12.

Suggested patch, untested:

--- a/src/ksp/ksp/impls/eksm/eksmimpl.h
+++ b/src/ksp/ksp/impls/eksm/eksmimpl.h
@@
-#define VEC_VV(i)  eksm->vecs[VEC_OFFSET + i]
+#define VEC_VV(i)  eksm->vecs[VEC_OFFSET + (i)]
--- a/src/mat/impls/nest/matnestimpl.h
+++ b/src/mat/impls/nest/matnestimpl.h
@@
-    PetscCall(PetscObjectQuery((PetscObject)A, "MatMultiShift", (PetscObject *)&container)); \
-    PetscCheck(container, PetscObjectComm((PetscObject)A), PETSC_ERR_ARG_WRONG, "The Mat is not a multi-shift matrix"); \
+    PetscCall(PetscObjectQuery((PetscObject)(A), "MatMultiShift", (PetscObject *)&container)); \
+    PetscCheck(container, PetscObjectComm((PetscObject)(A)), PETSC_ERR_ARG_WRONG, "The Mat is not a multi-shift matrix"); \

Style — Tutorial introduces decorative separator comments

Locations: src/ksp/ksp/tutorials/ex90.c:115–117, :219–221.

The dashed comment banners conflict with doc/developers/style.md, C Formatting, item 13, which prohibits decorative section separators.

Suggested patch, untested:

--- a/src/ksp/ksp/tutorials/ex90.c
+++ b/src/ksp/ksp/tutorials/ex90.c
@@
-  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-              Build the matrices and right-hand side
-     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+  // Build the matrices and right-hand side
@@
-  /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-                      Solve the shifted linear systems
-     - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
+  // Solve the shifted linear systems

Coverage and verification