Fix sparse LU/QR/Cholesky fill-reduction permutations

[micycle1]

Solves #207.

Problem

Sparse fill-reduction orderings were broken for any non-identity permutation.

LU and Cholesky could return plausible but incorrect results. Sparse QR could crash. This went unnoticed because the only orderings covered by tests were NONE and IDENTITY. FillReducing.RANDOM existed but was not used in the test suite.

There were six issues across LU and QR:

1. ApplyFillReductionPermutation_DSCC

   • Its documentation says it applies P*A*Q.
   • It only applied the row permutation P.
   • It ignored the column permutation Q.
   • It also required a row permutation and threw "No row permutation matrix" otherwise.

2. LuUpLooking_DSCC.performLU()

   • Used the row permutation array (getArrayP()) as the column ordering during factorization.

3. LinearSolverLu_DSCC.solve()

   • Did not apply the fill-reduction row permutation to b.
   • Restored x in the wrong direction: it used permute instead of permuteInv.

4. LinearSolverLu_DSCC.solveSparse()

   • Did not apply fill reduction at all.

5. QrLeftLookingDecomposition_DSCC.decompose()

   • Crashed with ArrayIndexOutOfBoundsException for nontrivial column permutations.
   • The symbolic QR phase analyzed columns in natural order, while the numeric phase processed permuted columns.

6. LinearSolverQrLeftLooking_DSCC.solve() and solveSparse()

   • Restored x using the row permutation in the forward direction.
   • They should use the column permutation in the inverse direction.

While fixing these, there was one more issue:

7. LinearSolverCholesky_DSCC.solveSparse()
   • Passed the fill inverse permutation as pinv to TriangularSolver_DSCC.solve().
   • That argument maps solution indices to factor columns; it is not a right-hand-side row permutation.
   • Sparse Cholesky solves were therefore also wrong under nontrivial orderings. The dense path was already correct.

Cause

Each solver handled fill-reduction permutations differently:

• Cholesky used a symmetric pinv approach.
• LU applied the row permutation to both sides.
• QR mixed physical row permutation with logical column ordering.

ApplyFillReductionPermutation_DSCC was intended to provide a common implementation, but it only handled matrix permutation. Each solver separately handled b and x, which led to inconsistent conventions and direction mistakes.

Fix

The common convention is now:

P * A * Q

The decomposition operates on P*A*Q, where P and Q come from the fill-reduction ordering. For symmetric orderings such as SYMRCM or AMD, Q = Pᵀ.

ApplyFillReductionPermutation_DSCC

It now applies both row and column permutations.

• Missing row or column permutations are treated as identity.
• This allows column-only orderings such as COLAMD.
• Symmetric mode exposes matching row/column permutations through the same interface.

Permutation handling is now centralized in this class:

• rowPermInv(numericRowPinv, work)

  • Combines fill-reduction row permutations with LU numeric pivots or QR structural pivots.
  • Used for both dense and sparse RHS paths.

• undoColumnPermutation(double[] ...)

• undoColumnPermutation(DMatrixSparseCSC ...)

  • Restore the solution using the fill-reduction column permutation.

All LU, QR, and Cholesky solve paths now follow the same pattern:

rowPinv = reduce.rowPermInv(decompositionPivots, work)
apply rowPinv to b
run the solver's triangular solves
undo the column permutation on x

The solver-specific parts are now limited to the decomposition pivot array and the triangular solve implementation.

LU

• LuUpLooking_DSCC now processes column k directly, since the input matrix has already been fully permuted.
• This removes the row-permutation-used-as-column-ordering bug.
• computeDeterminant() now includes the sign of both fill permutations.
• Updated the old CSparse pivoting comment to refer to the correct index.

QR

• QR no longer has separate symbolic and numeric column orders.
• Since the matrix is physically permuted before both phases, they see the same structure.
• This fixes the nontrivial-permutation crash.

API

LU and QR now expose getApplyFillReduction(), matching the existing Cholesky setup.

Tests

Added FillReducing.RANDOM to the permutation test cases for LU, QR, and Cholesky.

Also added or updated tests for:

• LU and QR factor reconstruction against Pfr * A * Qfr, rather than directly against A.
• A deterministic LU regression using two 3-cycle permutations. These are intentionally not self-inverse, so direction mistakes cannot accidentally cancel out.
• Sparse and dense RHS paths.
• QR wide-matrix behavior under RANDOM.

For wide matrices, QR may now return false under a fill-reduction ordering if the resulting structure requires column pivoting. Previously this case crashed.

Verification

• Full org.ejml.sparse dsparse suite: 358/358 passing

Compatibility

There are no public signature changes for the solve APIs.

• NONE is unchanged; the fill-reduction object is null and all new logic is skipped.
• IDENTITY gives the same numerical results.
• Behavior changes only for nontrivial fill-reduction permutations, where results were previously wrong or QR crashed.

One consequence is that LU and QR factors now reconstruct the fill-reduced matrix:

Pfr * A * Qfr

rather than A directly.

Callers reconstructing A from LU or QR factors need the fill-reduction permutation as well as the decomposition's numeric pivot data. It is available through getApplyFillReduction().

Follow-up

This should unblock SYMRCM and similar orderings. A symmetric ordering can return its ordering through both prow and pcol, and it will now work consistently with Cholesky, LU, and QR.