Internals

Workgroup size of the assembly kernels.

AlgorithmEnum

Union of the enums set through "default"/"algoN" strings (or the integer N).

BUNCH_KAUFMAN_ALPHA

α = (1 + √17)/8, the Bunch–Kaufman constant that minimizes the element growth bound of 1×1/2×2 pivoting.

CharEnum

Union of the enums spelled with characters (convert(E, ::AbstractChar), convert(Char, x)).

DEVICE_MAPS

Names of the index vectors of a Symbolic that live on the device (VI), in the order adapt moves them.

FGMRES_PROVIDER

Ref to the FGMRES driver of the Krylov.jl extension (SparseDirectSolverKrylovExt), set in its __init__; nothing while Krylov.jl is not loaded (ir_mode = "fgmres" then raises NotSupportedError). Called as provider(cache, A, P, R, len; atol, itmax) -> (D, iterations): FGMRES for A D = R[1:len] with right preconditioner P (mul!(y, P, x)), absolute residual tolerance atol, at most itmax iterations, storage kept in cache (a Ref{Any}).

Workgroup size of the fused regime-B kernel.

Default 1-D workgroup size of the single-workgroup kernels (potrf, getrf, trsm).

Smallest row class of a regime-B bin that the LDLᵀ/LDLᴴ factorization runs on the blocked regime-C path (ldlt_blocked_path) instead of the fused front kernel (together with LDLT_BLOCKED_MIN_WCLASS).

Smallest width class of a regime-B bin that the LDLᵀ/LDLᴴ factorization runs on the blocked regime-C path (ldlt_blocked_path) instead of the fused front kernel (together with LDLT_BLOCKED_MIN_FCLASS).

Pivot columns per block of the regime-C LDLᵀ/LDLᴴ path (panel_ldlt_kernel!).

Workgroup size of the regime-C LDLᵀ/LDLᴴ panel kernel.

Workgroup size of the LDLᵀ/LDLᴴ front kernel (regimes B and C).

Refinement steps of the solvers created by cholesky, ldlt and lu (PLAN §3.1; the handle layer uses 0).

Workgroup size of the LU front kernel (regimes B and C).

ORDERING_LEVEL_FLOPS

Default depth weight of ordering_cost: flops charged per level of the supernodal schedule, 1e8. See ordering_cost for the model.

Workgroup size of the permutation kernels.

PIVOT_PAIRS_MAX_ROUNDS

Maximum number of zero_pivot_pairs! rounds (each followed by a re-ordering) of the pivot_pairs = "default" analysis (8).

PIVOT_PAIR_TOLERANCE

Default relative tolerance of the 2×2 candidate test of pivot_pairs: row i is a candidate when |aᵢᵢ| ≤ τ · maxⱼ≠ᵢ |aᵢⱼ|, τ = the parameter "pivot_pair_tolerance" (default 1e-6).

Workgroup size of the residual-norm reduction (one workgroup per right-hand side).

REGIME_A, REGIME_B, REGIME_C

Regime codes (Int8) of Schedule.regime: fused subtree kernels, fused per-front level kernels, vendor dense calls on large fronts (PLAN §2.2).

Largest front width with a fused regime-B kernel; wider bins (raised regime_c_width) take the regime-C path.

Width classes with a fused regime-B kernel (@localmem holds the packed W×W lower triangle).

Workgroup size of the Schur kernels (metadata, LU assembly, export).

Workgroup size of the solve kernels.

Workgroup size of the per-phase reduction of the front statistics.

Bytes of a regime-A kernel's local memory kept for its control words and pivot (not for fronts).

Bytes of the regime-A LDLᵀ/LDLᴴ kernel's local memory kept for its control words, pivot block and the reduction slots of the cooperative pivot search (structures "S"/"H"; at most 944 bytes for ComplexF64 with Int64 maps).

@localmem sizes (bytes) of the regime-A kernel instances; a budget uses the largest one it holds. @localmem is static shared memory, so the usable sizes depend on the backend (max_local_bytes): 48 KiB on CUDA (ptxas rejects more static shared data), 64 KiB on AMDGPU (LDS per workgroup); the analysis keeps the budgets within the backend's cap.

Default subtree_parallelism: a regime-A subtree does at most this fraction (inverse) of the factorization flops.

Workgroup size of the fused regime-A kernel.

StringEnum

Union of the enums spelled with strings (convert(E, ::AbstractString), convert(String, x)).

SparseMatrixCSC(F::FullPatternMap, nzval) -> SparseMatrixCSC

The full matrix M described by F with the values taken from the user's nzval (host copy; for tests and debugging).

SparseMatrixCSC(P::SymmetricPattern) -> SparseMatrixCSC{Bool,Int}

The pattern as a Boolean sparse matrix (no diagonal), for tests and debugging.

SparseMatrixCSC(A::CSR[, k = 1])

Host copy of batch member k of the stored CSR matrix (size(A), index type INT, base rebased to one). The transposed flag is not applied: for B = csr_of_transpose(M), SparseMatrixCSC(B) == transpose(M).

AmalgamationParams

Type of the "amalgamation" parameter: (max_width, zero_fraction, min_width), the relaxed-amalgamation limits of PLAN §2.3 step 4 (maximum width of a merged panel; a fundamental supernode wider than that stays whole, maximum fraction of explicit zeros added to the factor, target minimum width).

BatchMap{V <: Union{Nothing, AbstractVector{Int32}}}

Batch part of a kernel launch (an argument of every numeric and solve kernel): members[1:nact] are the active batch members (1-based, increasing; members = nothing: a single matrix, member 1), nbatch the batch size, first the launch's first entry in its node list (the former first argument of the kernels) and nrhs the right-hand sides per member of a solve (1 in the numeric phase). Adapted to the device with its members vector.

ColumnPartition

A partition of the columns of an elimination tree into supernodes (contiguous column ranges after a renumbering), without row structure:

• order::Vector{Int}: column k of the supernodal numbering is column order[k] of the etree numbering (a topological order of the etree);
• super_ptr::Vector{Int}: supernode s owns columns super_ptr[s]:(super_ptr[s+1]-1) of the supernodal numbering;
• snparent::Vector{Int}: parent supernode, 0 for a root.

Returned by fundamental_supernodes and amalgamate.

DenseCapabilities

Result of capabilities: one Bool per probe.

• generic_mul, generic_trsm (triangular ldiv!/rdiv!), generic_cholesky (cholesky!), generic_lu (lu!): LinearAlgebra entry points.
• vendor_gemm, vendor_syrk, vendor_herk, vendor_trsm, vendor_potrf, vendor_getrf, vendor_sytrf, vendor_gemm_strided_batched, vendor_trsm_batched, vendor_potrf_batched, vendor_getrf_batched: the vendor_* bindings of the backend.
• reshape_view_mul: mul! accepts a panel view reshape(view(buf, r), f, w).
• atomic_add: Atomix.@atomic x[i] += v works on T in a KA kernel.

FactorPreconditioner(ws, S, N, nrhs, ncols, conjugate, deterministic, interrupt, transpose = false,
                     scaling = nothing)

op(A)⁻¹ through the factorization in N (permutation, forward, diagonal and backward sweeps in the solve workspace ws; conj(M⁻¹ conj(x)) when conjugate; the sweeps of Mᵀ when transpose, LU only) on the stacked compact columns (nrhs per active member, ncols in all) of the RefinementOperator; with the matching permutations and scalings of a SolveScaling scaling (solve_io). mul!(y, P, x) polls interrupt first (InterruptedError); asynchronous otherwise. The right preconditioner of FGMRES-IR.

FactorizationAlg

factorization_alg values (cuDSS cudssFactorizationAlg_t, reinterpreted, PLAN §1.3): "default" automatic, "algo1" very-sparse-factor path (regimes A/B only, no vendor calls), "algo2" vendor dense calls above the subtree regime.

FullPatternMap

Host description of the full matrix M behind a user CSR input (PLAN §2.3 step 1, used by the residual SpMV of refinement): M as a 1-based CSR (rowptr, colval, sorted columns, duplicates merged) plus, for every entry e of M, the user nzval entries that make it up:

M.nzval[e] = Σ_{s ∈ srcptr[e]:srcptr[e+1]-1}  conjflag[s] ? conj(nzval[src[s]]) : nzval[src[s]]

src holds 1-based positions in the user's nzval (first batch member; member k adds (k - 1) * nnz). For symmetric structures the triangle given by the view is mirrored; for "H"/"HPD" the mirrored entries are conjugated. The map describes the stored CSR matrix; the transposed flag of a CSR is applied by the solver (through solve_mode).

IRMode

ir_mode values (PLAN §1.7): "ir" plain iterative refinement, "fgmres" FGMRES with the factorization as preconditioner.

IndexBase

Index base of rowPtr/colVal: 'Z' zero-based, 'O' one-based (cuDSS cudssIndexBase_t).

Layout

Static memory layout of the numeric phase (PLAN §2.3 steps 6 and 8), element offsets, 1-based:

• panel_ptr (length ns + 1): panel of supernode s at panel_ptr[s]:(panel_ptr[s+1]-1), column-major f×w; factor_len entries in all;
• d_ptr (length ns + 1): D entries of supernode s at d_ptr[s]:(d_ptr[s+1]-1) (its columns); the subdiagonal of 2×2 pivots at n more; d_len = 2n;
• cb_ptr, cb_len: contribution block of s on the update stack (0/0 if it has none there); cb_first, cb_last: the steps it lives in (step 0 is regime A, see Schedule);
• step_top[t + 1]: last update-stack entry in use during step t (t = 0:nsteps); stack_len = maximum(step_top) is the high-water mark;
• work_len: entries of the regime-C syrk workspace (largest m^2 of a front on the regime-C path, takes_c_path, with a block on the stack; LDLᵀ/LDLᴴ: the largest sum of ldlt_c_work_len over a chunk of concurrent fronts on the blocked path, ldlt_blocked_path and ldlt_c_chunk_end);
• local_front, local_cb (regime A, 0 elsewhere): local-memory offset of the packed front of s and of its contribution block after the move (0 for a subtree root, whose block goes to the update stack); local_len[t]: entries of the serial stack of subtree t (its peak).

Matching

Result of compute_matching (host):

• alg::MatchingAlg: the job that ran (MATCHING_AUTO resolved to job 5);
• perm::Vector{Int}: the matching permutation, row i matched to column perm[i] (the "perm_matching" data parameter; for a structurally singular matrix the unmatched rows get the unmatched columns in increasing order);
• matched::BitVector: row i has a matched entry (false only for a structurally singular matrix);
• rscale, cscale::Vector{Float64}: row and column scaling factors ("scale_row", "scale_col"; ones for jobs 1–4, equal for the symmetric structures);
• symmetric::Bool: symmetric scaling only (structures other than "G"): the factored matrix is D A D, not permuted.

MatchingAlg

matching_alg values (cuDSS cudssMatchingAlg_t): "default" no matching, "algo1"–"algo5" MC64 jobs 1–5, "algo6" automatic (job 5).

MatchingState{VI, VR, VT}

Device state of the matching and scaling of a DirectSolver (matching_alg ≠ "default"), built by "symbolic_factorization" from the host Matching: io (the SolveScaling of the solve phases), weights (rᵢ cⱼ of every stored entry, entry_scaling) and values (the scaled values of every batch member, the input of the numeric phase, scale_values!).

MatrixView

Which part of the matrix is read: 'F' full, 'L' lower, 'U' upper triangle (cuDSS cudssMatrixViewType_t).

MemberPanels{V}

The panel pointers of member k of an nb-member batch in the interleaved factor layout: P[s] = nb (front_ptr[s] - 1) + (k - 1)(front_ptr[s+1] - front_ptr[s]) + 1 (the first entry of the f×w panel of supernode s of member k). Passed to the kernels in place of front_ptr for the factor offsets; P.ptr is the analysis's front_ptr. Note that P[s + 1] - P[s] is not a panel length.

Numeric{T, VT <: AbstractVector{T}, VS <: AbstractVector{Int64}, VI <: AbstractVector{Int32},
        VK <: AbstractVector{Int8}}

Numeric storage of a factorization (PLAN §3.2), laid out by the Layout of a Symbolic:

• factor (layout.factor_len entries): the panel of supernode s is the column-major f×w block factor[panel_ptr[s]:(panel_ptr[s+1]-1)] (leading dimension f, rows snrows(s), the upper triangle of its diagonal block is unused and kept zero);
• d (layout.d_len = 2n entries): D of LDLᵀ/LDLᴴ (unused by Cholesky): d[k] is the diagonal of D at factor column k, d[n + k] the subdiagonal entry D[k+1, k] of a 2×2 block starting at k (zero otherwise);
• piv (n Int32, LDLᵀ/LDLᴴ only): the local pivot order, piv[k] is the column of P A Pᵀ (supernodal numbering before pivoting) eliminated at factor column k; pivoting stays inside a supernode, so piv permutes each sncols(s) (see ref_ldlt!);
• pivot_kind (n Int8, LDLᵀ/LDLᴴ only): PIVOT_KIND_1X1, PIVOT_KIND_PERTURBED, PIVOT_KIND_2X2_FIRST, PIVOT_KIND_2X2_SECOND (0: not factored);
• stack (layout.stack_len entries): the update stack of the device path (contribution block of s at cb_ptr[s], packed lower triangle of the m×m block, m(m+1)/2 entries);
• work (layout.work_len entries): the full m×m syrk/herk result of a regime-C front before it is packed into its block;
• stats (FRONT_STATS_FIELDS × ns Int64, column s = front s): (npos, nneg, nzero, nperturbed, n2x2, info), info = the local column of the first failed pivot of the front (0 = none);
• totals (FRONT_STATS_FIELDS Int64): the sums of stats over the fronts (reduce_stats!; info field: number of failed fronts);
• psign (n Int8, LDLᵀ/LDLᴴ): the pivot_sign request per original row (0 = none), copied from the options at every factorization;
• aux (1 entry): max |aᵢⱼ| of the values of the last LDLᵀ/LDLᴴ factorization with pivot_epsilon_alg = "algo1" (scaled perturbation);
• info (ns + 1 Int32): device status of the numeric phase, the potrf info of front s at s, the reduced result (smallest failed factor column, 0 = none) at ns + 1;
• plan: the host NumericPlan;
• members (nbatch Int32): the active batch members of the last phase (first plan.nact[] entries, see BatchMap);
• nbatch: the number of matrices of the uniform batch (PLAN §3.5);
• ufactor, ustack (structure "G" only, empty otherwise): the Uᵀ structure of LU (src/numeric/lu.jl), laid out exactly as factor and stack: the panel of s holds U[1:w, :]ᵀ (unit diagonal stored), the packed contribution block the transposed strict upper triangle of the m×m block (its diagonal stays zero). For LU, factor holds the unit-lower L, d the diagonal D of L D U and piv the local row pivot order.

A uniform batch of nbatch members multiplies every length above by nbatch, in the layouts of src/numeric/batch.jl: interleaved per front for factor (panel_offset), interleaved per entry for info, member after member for the others (the regime-C workspace holds the m×m blocks of the members of one front).

NumericPlan

Host launch plan of the level driver (factorize!), derived from the Schedule once at allocation: for every step t of regimes B/C, step_first[t]:step_last[t] is its range in schedule.group_nodes (empty when step_first[t] > step_last[t]) and step_maxchild[t] the largest number of children of one of its fronts (the trip count of the owner-pull extend-add). The launch groups of regimes B/C, in execution order: group k covers group_first[k]:group_last[k] of group_nodes, its fronts have at most group_maxchild[k] children, and group_width[k] is the width class W of its fused regime-B kernel (factorize_fronts_b!), or 0 when the group takes the regime-C path (regime C, or a regime-B bin wider than REGIME_B_MAX_WIDTH). The regime-A launch groups (one per budget class, run first): group k covers the subtree ids sub_first[k]:sub_last[k] of group_nodes and its kernel has sub_local[k] bytes of local memory (subtree_local_bytes of its budget). info_host is the host staging buffer of the one info read per phase (one entry per batch member). Batch state (PLAN §3.5): members_host[1:nact[]] are the active members of the last phase (mirrored in numeric.members), runs the starts of their runs of consecutive members (member_runs!, the strided batches of regime C). vendor_ptrs[s] (regime-C front s of a batch on a backend with batched vendor potrf/trsm, else nothing): the device vectors of member pointers (F11, F21) of its nbatch panels, built once at allocation (vendor_batch_pointers). Schur complement mode (PLAN §3.6): the launch group of the Schur root (schedule.schur) is not one of the B/C groups; schur_first is its position in group_nodes (0 = no Schur root) and schur_maxchild its number of children (assemble_schur!).

Ordering

Result of compute_ordering:

• perm, iperm: the fill-reducing permutation (perm[k] is the original index of the k-th pivot, so the factorized matrix is A[perm, perm]) and its inverse;
• alg_used: :natural, :amd, :mmd, :nd or :user;
• stats::NamedTuple: nnz_L, flops, nlevels (column-etree depth), sdepth (supernodal schedule depth), nsupernodes, cost of the chosen ordering, candidates::Vector{OrderingCandidate} (every ordering evaluated, the chosen one included), auto::Bool (whether the automatic choice ran), level_flops (the depth weight of ordering_cost used);
• pairs::Vector{Tuple{Int,Int}}: the 2×2 pivot candidate pairs (partner, candidate) the ordering was computed with (empty for no pair): perm puts each candidate right after its partner, and the symbolic factorization runs on factor_pattern(P, ordering). stats.pair_rounds is the number of zero_pivot_pairs! rounds that added pairs for the chosen ordering (0 without that pairing).

OrderingCandidate

Evaluation of one ordering by the cost model of compute_ordering (evaluate_ordering): fields alg::Symbol, nnz_L::Int, flops::Float64, nlevels::Int (column-etree depth), sdepth::Int (supernodal schedule depth), nsupernodes::Int (fundamental supernodes), cost::Float64 (ordering_cost).

PairCandidates

2×2 pivot candidates of a matrix, from pivot_candidates (host arrays):

• candidate::BitVector: row i has an absent diagonal or |aᵢᵢ| ≤ τ maxⱼ≠ᵢ |aᵢⱼ|, τ = PIVOT_PAIR_TOLERANCE;
• w::Vector{R}: |aᵢⱼ| (duplicates summed first), aligned with P.rowval of the pattern P it was computed on (R = real(T));
• amax::Vector{R}: maxⱼ≠ᵢ |aᵢⱼ| per row (0 for a row without off-diagonal entries).

Phase

Execution phases with the cuDSS bit values (cudssPhase_t). Composite phases are the bitwise OR of their parts: "analysis" = reordering | symbolic factorization, "solve" = all six solve sub-phases, and the CUDSS.jl shorthands "solve_fwd_schur" = fwdperm | fwd and `"solvebwdschur"` = bwd | bwdperm. See phase_includes.

PivotEpsilonAlg

pivot_epsilon_alg values (cuDSS cudssPivotEpsilonAlg_t): "default", "algo1" scaled perturbation, "algo2" static perturbation.

PivotPairsMode

pivot_pairs values (beyond cuDSS), for the analysis of "S"/"H" matrices: "default" (2×2 pivot pairs for the candidates whose pivot is structurally zero in the ordering, see zero_pivot_pairs!), "all" (every candidate paired, see pivot_pairs) and "none".

PivotType

pivot_type values (cuDSS cudssPivotType_t): 'A' auto, 'N' none, 'C'/'R' global column/row pivoting (not supported, PLAN §3.3), 'D' diagonal, 'L' local block, 'B' Bunch–Kaufman 1×1/2×2.

RefinementOperator(W, nzval, n, ncols, bm, conjugate, roff = 0)

op(A) (the full matrix M of the RefinementWorkspace W with the values nzval, conj(M) when conjugate; roff = n: Mᵀ from the rows of a "G" map) acting on the n × ncols compact columns stacked into vectors of length n ncols (column r uses the values of the member of compact column r of the BatchMap bm). Supports size, eltype and mul!(y, op, x) (one launch, asynchronous): the matrix operator of FGMRES-IR.

RefinementWorkspace{T, R, VI, MT, VR}

Device storage of iterative refinement, allocated at the first solve that refines (allocate_refinement) and grown with the number of right-hand sides:

• rowptr, colval, src (INT vectors): the full matrix M as a CSR over the contributions of the user's nzval (src[c] > 0: nzval[src[c]], src[c] < 0: conj(nzval[-src[c]])), from refinement_map; for "G" the rows n+1:2n of rowptr are the rows of Mᵀ;
• R (n × nrhs): the residual, then the correction;
• Bc (n × nrhs): a copy of the right-hand side when X and B alias;
• norms (2 nrhs, real): ‖Rₖ‖₂² and ‖Bₖ‖₂² per right-hand side, and norms_host, its host copy (the early-exit test of ir_tol);
• krylov: the storage of FGMRES-IR (ir_mode = "fgmres"), owned by the Krylov.jl extension (nothing until its first solve).

ReorderingAlg

reordering_alg values (cuDSS cudssReorderingAlg_t): "default" automatic AMD/ND choice by flops plus a charge per level of the supernodal schedule (the same on every backend), "algo1" BTF+COLAMD and "algo2" COLAMD (both fall back to the symmetric-pattern path), "algo3" AMD, "algo4" nested dissection, "algo5" natural ordering.

Schedule

Regime assignment and launch order of the supernodes of a SupernodePartition for element type T (budgets are bytes):

• width[s], rows[s]: front width w and number of rows f;
• level[s]: height of supernode s above the leaves of the supernodal tree (leaves 1, level[parent] > level[child]), nlevels the maximum;
• regime[s] ∈ (REGIME_A, REGIME_B, REGIME_C);
• slevel[s]: schedule level of a B/C front (height counting only B/C fronts), 0 in regime A; nslevels the maximum;
• step[s]: launch step of a B/C front (1:nsteps, one per (level, chunk)), 0 in regime A; step_level[t], step_chunk[t] describe step t;
• bin[s]: regime-B bin id ((wi - 1) * length(fclasses) + fi for width class wclasses[wi] and row class fclasses[fi]), 0 otherwise;
• regime A subtrees: subtree[s] (0 outside A), subtree_ptr/subtree_nodes (the nodes of subtree t, in processing order: a postorder whose children are visited by decreasing peak - cb, which minimizes the stack), subtree_root, subtree_peak (bytes of the serial stack of packed fronts and contribution blocks), subtree_class (index of the smallest budget of budgets whose subtree_capacity holds the peak);
• groups and group_nodes: the launch groups in execution order (regime A classes first, then step by step: B bins in increasing bin id, then C);
• wclasses, fclasses, budgets, regime_c_width, regime_c_rows, memory_budget, elsize (bytes per packed entry, sizeof(T); 2 sizeof(T) for LU), vendor_c (whether regime C uses vendor calls; factorization_alg = "algo1" turns them off);
• schur: the Schur complement root (schur_mode = 1, PLAN §3.6), a regime-C front alone in the last step that the numeric phase assembles but does not factor, or 0.

ScheduleGroup

One launch group of a Schedule: step (0 for regime A), schedule level and chunk (0 for regime A), regime, class (regime A: index into budgets; regime B: bin id; regime C: 0) and the range first:last of Schedule.group_nodes it covers. For regime A the range holds subtree ids, for B and C supernode ids.

ScheduleKind

schedule values (PLAN §1.7): "auto", "subtree+level" (portable baseline), "syncfree" (CUDA/ROCm only).

SchurState

Host state of the Schur complement mode of a DirectSolver (PLAN §3.6), set by "reordering" when schur_mode = 1: flags (the Schur rows and columns of user_schur_indices), indices (them, increasing: row k of S is row indices[k] of A), rowptr/colval (the symbolic pattern of S, schur_pattern: 1-based CSR, both triangles, diagonal included) and dest, dest_view (the destination registered with setparam!(solver, "schur_matrix", …), or nothing).

SolveAlg

solve_alg values (cuDSS cudssSolveAlg_t, reinterpreted, PLAN §1.3): "default" level-batched sweeps, "algo1" partitioned-inverse diagonal blocks.

SolvePlan

Host launch plan of the sweeps, in forward order (the backward sweep runs it in reverse): launch k has kind kind[k] (SOLVE_SUBTREES: the subtree ids group_nodes[first[k]:last[k]], all regime-A subtrees in one launch; SOLVE_FRONTS: the regime-B fronts of one step; SOLVE_DENSE: the regime-C-path fronts of one step, solved one after the other with dense calls; SOLVE_SCHUR: the unfactored Schur root of schur_mode = 1, skipped by the sweeps except for the owner-pull of its children's updates in the deterministic forward sweep). maxm is the largest m = f - w of a regime-C-path front.

SolveScaling{VI, VR}

Device permutations and scalings of the solve with matching (matching_alg ≠ "default", Matching): the factored matrix is M = Dr A Dc Q ("G"; Q the matching column permutation) or D A D (symmetric structures), so A x = b is solved as M y = Dr b, x = Dc Q y:

• rperm: factor row k is original row rperm[k] (the symbolic perm);
• cperm: factor column k is original column cperm[k] (q[perm[k]] for "G", perm otherwise);
• rscale, cscale: the row and column scalings Dr, Dc (original numbering, real(T) vectors).

SolveWorkspace{T, MT <: AbstractMatrix{T}}

Device storage of the solve phase for up to nrhs right-hand sides, allocated once by allocate_solve:

• Y (n × nrhs): the permuted right-hand side, overwritten by the solution;
• U (length(rowval) × nrhs): the per-front update buffers of the deterministic forward sweep (front s uses rows rowptr[s] + w : rowptr[s+1] - 1);
• tmp (maxm × nrhs): the gemm result of a regime-C-path front;
• atomic: whether the backend has an atomic add for T (capabilities); without one the forward sweep is always the deterministic variant;
• plan: the host SolvePlan.

Structure

Matrix structure, "G", "S", "H", "SPD" or "HPD" (cuDSS cudssMatrixType_t).

SupernodePartition

Supernodes of the factor with their row structure (PLAN §2.3 steps 3–4), in the supernodal numbering: column k of the factor is original column perm[k].

• n, perm, iperm: size, final fill-reducing permutation (the T05 ordering composed with the supernodal renumbering) and its inverse;
• parent, counts: column elimination tree and true column counts (diagonal included) in the supernodal numbering;
• super_ptr: supernode s owns columns super_ptr[s]:(super_ptr[s+1]-1); col2sn[j] is the supernode of column j;
• snparent (0 = root) and snpost, a postorder of the supernodal tree;
• rowptr, rowval: the rows of supernode s are rowval[rowptr[s]:(rowptr[s+1]-1)] (see snrows), sorted, so its own columns come first; the panel is length(rows) × width;
• nnz_stored: entries of the stored lower-trapezoidal panels, explicit zeros of amalgamation included; nnz_L: true nonzeros of L; flops: Cholesky flop count on the stored panels (Σ over panel columns of rows², real arithmetic, the cholesky_flops convention);
• amalgamated: whether amalgamate ran.

SupernodePartition(P::SymmetricPattern, perm, parent, counts, cp::ColumnPartition; amalgamated)

Supernodal symbolic factorization: the rows of every supernode of cp are its own columns plus the union of the rows of A[perm, perm] in its columns and of the rows of its children, restricted to rows after its last column. perm is the T05 ordering and parent, counts the etree and column counts in its numbering.

Symbolic{INT, VI <: AbstractVector{INT}}

Result of the analysis of one sparsity pattern (PLAN §2.3, §3.2): host tree data, schedule and layout, plus the device maps as VI vectors.

Host fields: n, nnz (stored entries of the user's CSR), structure, view, index, elsize (sizeof(T) the schedule was built for), partition::SupernodePartition, schedule::Schedule, layout::Layout.

Device maps (DEVICE_MAPS; supernodal numbering, 1-based, element offsets):

• perm, iperm (length n): the fill-reducing permutation of the factor and its inverse (column k of L is original column perm[k]);
• super_ptr, snparent: supernode column ranges and parents (0 = root);
• solve gather lists rowptr, rowval: the (permuted) rows of supernode s are rowval[rowptr[s]:(rowptr[s+1]-1)], its own columns first;
• per-front descriptors: front_ptr (= layout.panel_ptr, length ns + 1), front_nrows (f), front_ncols (w), cb_ptr (update-stack offset of the contribution block, 0 if it does not go through the stack), snparent;
• child_ptr, child_list: the children of s in the fixed order their contribution blocks are pulled (increasing id), for owner-pull extend-add;
• relind_ptr, relind: for a child c with m = f_c - w_c > 0, relind[relind_ptr[c]:(relind_ptr[c+1]-1)] are the positions (1-based) of its rows w_c+1:f_c among the rows of its parent;
• amap (length nnz): for every entry p of the user's nzval, the factor-buffer offset it is added to; 0 when the view ignores the entry, and -offset when the conjugate is added (a "H"/"HPD" entry that lands in the upper triangle of P A Pᵀ and is mirrored; for "G": an upper-triangle entry that goes to the Uᵀ structure). Duplicated entries map to the same offset and are summed;
• amap_ptr, amap_src: the same map grouped by supernode for owner-pull assembly: supernode s receives the entries amap_src[amap_ptr[s]:(amap_ptr[s+1]-1)] (positions in nzval, sorted by destination, then by position);
• subtree_ptr, subtree_nodes: regime-A subtree descriptors (the supernodes of subtree t in processing order, see Schedule);
• local_front, local_cb: local-memory offsets of the packed front and of the contribution block of a regime-A supernode (see Layout);
• group_ptr, group_nodes: the launch groups of the schedule (group g covers group_nodes[group_ptr[g]:(group_ptr[g+1]-1)]: subtree ids for regime A, supernode ids otherwise).

Symbolic(sp::SupernodePartition, schedule, layout, rowptr, colval, n, structure;
         view = 'F', index = 'O') -> Symbolic{Int, Vector{Int}}

Assemble the host Symbolic from the supernodes, schedule and layout of a pattern and its user CSR arrays (only the pattern is read).

SymmetricPattern

Adjacency structure of an n × n symmetric sparsity pattern on the host: colptr (length n + 1) and rowval hold, for every column j, the sorted row indices i ≠ j with (i, j) in the pattern, 1-based, without duplicates and without the diagonal, both triangles stored ((i, j) present iff (j, i) present).

SymmetricPattern(n, colptr, rowval)
SymmetricPattern(A::CSR, structure; view = 'F')
SymmetricPattern(rowptr, colval, n, structure; view = 'F', index = 'O', device = nothing)

The first form wraps arrays that already satisfy the invariants (checked). The others build the pattern from a CSR pattern, following PLAN §2.3 step 1:

• structure "S", "H", "SPD", "HPD" (or a Structure): only the triangle selected by view is read ('L': col ≤ row, 'U': col ≥ row; 'F' reads the lower triangle, as cuDSS does) and mirrored;
• "G": the pattern of A + Aᵀ; the view is ignored (every stored entry is read, as in cuDSS);
• duplicates are dropped, index ('O'/'Z') is the base of rowptr/colval, non-square matrices and out-of-range column indices raise InvalidValueError.

abs_max!(aux, nzval[, bm]) -> aux

aux[k] = max |nzval of member k| for the active members k of the BatchMap bm (default: a single matrix, aux[1]), one workgroup per member, @localmem tree reduction; the scale of pivot_epsilon_alg = "algo1". Asynchronous.

active_members(nbatch, ubatch_index, ubatch_mask) -> Vector{Int}

The batch members a phase processes (1-based, increasing): all of them for ubatch_index = -1 and no mask, member ubatch_index + 1 (ubatch_index is 0-based, as in cuDSS), or the members whose ubatch_mask flag is 1 (the mask has nbatch entries; with both set, their intersection). Raises InvalidValueError for an index or mask that does not fit the batch.

add_correction!(X, Y, perm; nrhs, transposed = false, conjugate = false, bm = nothing, scale = nothing) -> X

X[perm[k], r] += Y[k, r] (conj(Y[k, r]) when conjugate; Y[k, r] times scale[perm[k]] first with a matching scale vector): the solution of the correction system, in factor order in the workspace Y, added to the user array X of nrhs right-hand sides (uniform batch: the compact columns of Y to the active members' columns of X, BatchMap bm). One launch, conflict-free (perm is a permutation).

allocate_numeric(symbolic, T, backend = CPU(); nbatch = 1) -> Numeric{T}

Allocate (zero-filled) the factor panels, D, update stack, regime-C workspace, per-front statistics and their totals, status vector, pivot order, pivot kinds, pivot sign requests, aux (and for "G" the Uᵀ panels and update stack) of symbolic's Layout for element type T on the KernelAbstractions backend, for a uniform batch of nbatch matrices (every length times nbatch, all members active), and build its NumericPlan. This is the only allocation of the numeric phase. Raises InvalidValueError when the panels or the regime-C workspace of the batch do not fit the index type INT of symbolic, or when a regime-A kernel of the analysis needs more local memory than backend has (max_local_bytes).

allocate_refinement(map, T, INT, backend, nrhs; n = rows of map) -> RefinementWorkspace

Move the host refinement_map map to backend as INT vectors (overflow-checked) and allocate the residual storage for nrhs right-hand sides of element type T and n rows (pass n for a "G" map, which holds 2n rows).

allocate_solve(symbolic, T, backend = CPU(), nrhs = 1) -> SolveWorkspace{T}

Allocate (zero-filled) the solve storage of SolveWorkspace for nrhs right-hand sides of element type T on backend, build the SolvePlan and look up the atomic capability once. This is the only allocation of the solve phase.

amalgamate(sn::ColumnPartition, parent, counts, params) -> ColumnPartition

GPU-tuned relaxed amalgamation (Ashcraft–Grimes) of the exact supernodes sn (fundamental supernodes, see fundamental_supernodes); parent and counts are in the etree numbering. params = (max_width, zero_fraction, min_width) is opts.amalgamation:

1. the units are the fundamental supernodes themselves. A supernode wider than max_width stays whole: it is a dense separator that goes to regime C as a single front; splitting it into a chain of max_width panels multiplied the update-stack footprint by the chain length and added tree levels (issue #48);
2. bottom-up over the supernodal tree, the children c of each supernode p are tried in decreasing order of their off-diagonal row count (the order of increasing explicit zeros per column, which merging does not change) and merged into p when the merged width is ≤ max_width, the explicit zeros of the whole factor stay ≤ zero_fraction · nnz(L), and either the merged panel has at most zero_fraction explicit zeros per true nonzero or its width is still ≤ min_width (tiny panels are merged more eagerly).

Merging c (width w_c, f_c rows) into p (f_p rows) gives a panel of width w_c + w_p and w_c + f_p rows, adding w_c (w_c + f_p - f_c) explicit zeros. The result is renumbered so that every merged supernode is contiguous (postorder of the merged tree, columns in their previous order).

analysis_pairs(P::SymmetricPattern, rowptr, colval, nzval, n, structure, opts; view, index)
    -> (; pairs, candidates)

The 2×2 pivot pair input of compute_ordering for the analysis: nothing when pairs_enabled is false; pivot_candidates for pivot_pairs = "default" (pairs for the structurally zero pivots); the fixed pivot_pairs for "all". Reads nzval on the host once. Pairs are decided from the values present at "analysis" (the first batch member): an all-zero nzval gives no pairs and an undefined one arbitrary pairs, so run "analysis" after the first assembly of the matrix (MadNLP: after the first KKT assembly); "all" reads the values as well.

assemble_schur!(numeric, symbolic, nzval) -> numeric

Schur complement mode (PLAN §3.6): assemble the Schur root front (symbolic.schedule.schur, run by factorize! after every other front) without factoring it: zero its panel, scatter the values nzval of A into it and extend-add its children's contribution blocks (zero_fronts!, scatter_A!, extend_add!; for LU one kernel doing the same on the L and Uᵀ structures). The panel then holds the Schur complement S (schur_matrix!). The root's factor columns get the identity local pivot order, d = 1 (unit 1×1 pivots, so diagonal_sweep! leaves the Schur block unchanged), zero statistics and status 0. A no-op without a Schur root. Single matrices only (Schur mode refuses uniform batches). Asynchronous.

assembly_map(sp, layout, rowptr, colval, n, structure; view = 'F', index = 'O') -> amap

The amap of Symbolic: for every stored entry p of the user's CSR pattern (rowptr, colval, base index), its offset in the factor buffer of layout, -offset when its conjugate is stored instead, 0 when view ignores it. Entry (r, c) lands at (i, j) = (iperm[r], iperm[c]) of P A Pᵀ, or at (j, i) when i < j. Structure "G" (LU, two lower-triangular structures L and Uᵀ with the same layout, see src/numeric/lu.jl): an entry with i ≥ j goes to offset off of the L structure, an entry with i < j to the offset of (j, i) in the Uᵀ structure, encoded as -off.

backward_sweep!(ws, symbolic, numeric; nrhs = max_rhs(ws), impl = :auto, transpose = false) -> ws

Backward sweep Lᴴ X = W (Lᵀ for complex symmetric LDLᵀ) in place on the first nrhs columns of ws.Y, the forward plan in reverse: per step the regime-C-path fronts (gather kernel, gemm with L21ᴴ, trsm with L11ᴴ) and one launch over the regime-B fronts, then one launch over all regime-A subtrees (supernodes in reverse processing order). Gather-based and conflict-free: no atomics, deterministic. Uniform batch as in forward_sweep!. LU ("G"): U X = W, or Lᵀ X = W followed by the inverse local row order with transpose = true. Schur complement mode: the rows of the Schur root hold its solution and are only read. Asynchronous.

batch_map(numeric; first = 1, nrhs = 1) -> BatchMap

The BatchMap of the active members of numeric for a launch.

build_layout(sp::SupernodePartition, schedule::Schedule; ldlt = false) -> Layout

Panel offsets in supernode order, D offsets, and the update-stack offsets of the contribution blocks that leave their front through global memory (B/C fronts and regime-A subtree roots that have a parent; packed lower triangles). Blocks whose lifetimes [step(s), step(parent)] overlap never share entries; the offsets are the placement with the lowest high-water mark among a step-by-step first fit and two offline placements (lowest free offset, largest blocks first and largest size × lifetime first; issue #48). Also the regime-C workspace (with ldlt, structures "S"/"H", the larger one of the blocked LDLᵀ/LDLᴴ path, ldlt_c_work_len) and the local-memory offsets of the regime-A subtrees (subtree_local_layout).

build_schedule(sp::SupernodePartition, opts::Options = Options(), ::Type{T} = Float64;
               elsize = sizeof(T), reserve = SUBTREE_LOCAL_RESERVE, max_local = typemax(Int)) -> Schedule

Regime assignment, binning, levels and launch groups (PLAN §2.3 step 5) for the supernodes of sp, with byte budgets for the element type T (elsize bytes per packed entry: 2 sizeof(T) for LU, whose fronts and contribution blocks are two packed triangles, see src/numeric/lu.jl):

1. regime C: fronts with w > opts.regime_c_width or f > opts.regime_c_rows (every non-A front with factorization_alg = "algo2");
2. regime A: a front is eligible when it is not C, its children are eligible and its subtree's serial stack peak (packed lower-triangular f×f fronts plus the packed m×m contribution blocks waiting for their parent, m = f - w) fits the subtree_capacity of the largest of opts.subtree_budgets (local memory minus reserve bytes, subtree_local_reserve of the structure; the budgets are those of resolve_subtree_budgets for the backend's max_local bytes of local memory, max_local_bytes), and its subtree's flops (front_flops) are at most total ÷ opts.subtree_parallelism (0: no flop limit), so that a large tree of small fronts is split into enough subtrees to fill the device instead of running on one workgroup; the maximal eligible subtrees are the regime-A subtrees, each with the smallest budget class that holds it (empty subtree_budgets disables regime A);
3. regime B: the rest, binned by the smallest width class in 8, 16, 32, 64, … ≥ w and row class in 64, 128, 256, 512, … ≥ f (up to the regime-C thresholds);
4. schedule levels of the B/C fronts, split into chunks whose produced (packed) contribution-block bytes stay ≤ opts.memory_budget (a single front larger than the budget gets its own chunk; a negative budget means no chunking).

schur (a root supernode, 0 = none): the Schur complement root of schur_mode = 1 (PLAN §3.6), forced into regime C and onto a schedule level of its own above every other front, so that it is the only front of the last launch group (assembled, never factored, see assemble_schur!).

capabilities(backend, T) -> DenseCapabilities

Dense-layer capabilities of backend (a KernelAbstractions backend) for element type T, probed on first use with small try/catch calls and cached per (backend, T). A probe counts as available only if it neither throws nor returns a wrong result.

cholesky_flops(counts) -> Float64

Floating-point operation count of the Cholesky factorization with column counts counts (diagonal included), Σⱼ cⱼ² as reported by CHOLMOD (one division and cⱼ - 1 multiply–adds per entry of the update, real arithmetic). Multiply by 4 for complex element types.

cholesky_stats!(numeric, symbolic) -> numeric

Fill the per-front statistics of a Cholesky factorization from the per-front potrf status numeric.info[1:ns] (npos = w, or the failed local column minus one; info = the failed local column; the unfactored Schur root counts nothing) and reduce the smallest failed factor column into numeric.info[ns + 1] (0 = none), for every active batch member. One workgroup per member, no atomics. Asynchronous.

colcounts(P::SymmetricPattern, perm, parent, post) -> counts::Vector{Int}

Number of nonzeros in each column of the Cholesky factor L of A[perm, perm] (diagonal included), by the Gilbert–Ng–Peyton algorithm in O(nnz(A) α). parent comes from etree and post from postorder.

compressed_pattern(P::SymmetricPattern, pairs) -> (Pc::SymmetricPattern, ptr, members)

The pattern compressed by the 2×2 pairs (partner, candidate): a pair is one vertex whose adjacency is the union of both rows (self-loops removed), every other column is one vertex. Compressed vertex g stands for the columns members[ptr[g]:(ptr[g+1]-1)] (partner first); vertices are numbered by increasing first column (the partner for a pair).

compute_matching(rowptr, colval, nzval, n, structure, opts::Options; view = 'F', index = 'O')
    -> Union{Nothing, Matching}

The Matching of the analysis for opts.matching_alg (nothing for "default"), from the CSR matrix (rowptr, colval, nzval; the first batch member when nzval is a matrix) read with the view and index rules of SymmetricPattern:

• "G": mc64 on the stored rows (rscale, cscale from job 5);
• the symmetric structures: mc64 on the full matrix (both triangles), then the symmetric scaling dᵢ = sqrt(rᵢ cᵢ) of Duff & Pralet (2005), so |dᵢ aᵢⱼ dⱼ| ≤ 1 (rscale = cscale = d). An index left unmatched by a structurally singular matrix gets dᵢ = 1 / maxⱼ |aᵢⱼ| dⱼ over its matched neighbours (1 without one).

Reads nzval on the host once.

compute_ordering(P::SymmetricPattern, opts::Options; T = Float64, alg = nothing, pairs = [], candidates = nothing,
                 level_flops = ORDERING_LEVEL_FLOPS) -> Ordering

Fill-reducing ordering of the pattern P (PLAN §2.3 step 2):

• opts.user_perm set: that permutation, validated, accepted 0- or 1-based and returned 1-based (alg_used = :user);
• reordering_alg = "algo5": natural ordering; "algo3": AMD (CliqueTrees.AMD()); "algo4": nested dissection (nd_algorithm): METIS_NodeND with ufactor = nd_ubfactor, nseps = nd_nseps and seed = nd_seed (-1: METIS defaults; on fill-bound problems nd_nseps = 4 trades a slower ordering for less fill); nd_nlevels is the cuDSS minimum number of dissection levels, which METIS' full recursion meets on every graph large enough to be split that often; "algo1"/"algo2": AMD on the symmetric pattern;
• "default": AMD and ND are both evaluated with evaluate_ordering and the one with the lower ordering_cost flops + level_flops × sdepth (supernodal schedule depth) is chosen; stats.candidates lists both with their schedule depths.

The keyword alg (:natural, :amd, :mmd, :nd, :auto) overrides reordering_alg (MMD has no cuDSS spelling); T scales the flop counts in stats (complex: 4×) and does not change the choice. level_flops is the depth weight of ordering_cost; the schedule depth is that of the amalgamation of opts (opts.amalgamation, none with use_superpanels = 0).

2×2 pivot pairs ("S"/"H", issues #64 and #66), from analysis_pairs in the analysis when pairs_enabled:

• pairs: disjoint (partner, candidate) pairs fixed in advance (pivot_pairs = "all", pivot_pairs);
• candidates: a PairCandidates (pivot_pairs = "default"): for every ordering algorithm evaluated, the plain ordering is computed, the candidates whose pivot is structurally zero in it are paired (zero_pivot_pairs!, partners accepted down to opts.pivot_threshold of the row maximum) and the graph compressed by the pairs is re-ordered, until a round adds no pair (at most PIVOT_PAIRS_MAX_ROUNDS). Each algorithm keeps its own pairs. The automatic choice takes the ordering with the fewest structurally zero pivots left (structural_zero_pivots) and compares the cost only among those: a zero pivot is perturbed in the factorization, an accuracy loss that no saving in flops or depth pays for (on the slack KKT generators ND's search stops with zero pivots left where AMD's has none).

Both come from the values of the matrix at analysis: an all-zero nzval gives no pairs and an undefined one arbitrary pairs, so the analysis must run after the first assembly of the values (MadNLP: after the first KKT assembly), for "default" and "all" alike.

AMD, MMD and ND then order the graph compressed by the pairs (compressed_pattern: a pair is one vertex with the union adjacency), and the order is expanded with each partner right before its candidate (expand_permutation). Every ordering is evaluated, and the factor is built, on pair_pattern, where both columns of a pair have the union structure, so the pair lies in one fundamental supernode (never split by amalgamate) and in-front Bunch–Kaufman can take it as a 2×2 pivot. With user_perm (and with the natural ordering) the pairs are not applied, the user owns the order: ordering.pairs is empty. Without pairs (no candidate, or no structurally zero pivot) the result is the one of the plain ordering, bitwise.

compute_schur_ordering(P::SymmetricPattern, opts::Options, schur::AbstractVector{Bool}; T = Float64,
                       level_flops = ORDERING_LEVEL_FLOPS) -> Ordering

Schur-constrained ordering (PLAN §3.6): the rows and columns flagged in schur (schur_flags) are ordered last, in increasing original order, and the others first. The fill-reducing ordering of compute_ordering (reordering_alg) runs on the pattern induced by the other vertices (induced_pattern); since the Schur vertices come last they are never intermediate vertices of a fill path, so this is the ordering of the factored part. A user_perm keeps its relative order of the other vertices. The 2×2 pivot pairs of "S"/"H" are not used (pairs empty). stats describe the whole ordering (the Schur block included).

copy_rhs!(C, B; nrhs, transposed = false) -> C

C[:, 1:nrhs] = the right-hand sides of the user array B (column-major copy; the refinement keeps B when the solution overwrites it). One launch.

dense_impls(op::Symbol, backend, T) -> Vector{Symbol}

The implementations of dense op op (:gemm, :syrk, :herk, :trsm, :potrf, :getrf, :laswp, :gemm_strided_batched, :trsm_strided_batched, :potrf_batched) available for element type T on backend according to capabilities, in :auto preference order (:vendor, :generic, :ka). :ka is always present.

device_assembly_map(backend, sp, layout, rowptr, colval, n, structure; view, index)
    -> amap::Vector{Int}

assembly_map on backend for the symmetric structures with a triangular stored view; identical output. The caller checks device_maps_supported.

device_group_amap(backend, amap, layout, ns) -> (amap_ptr, order)

_group_amap (unsigned) on backend: identical order and per-owner counts. The signed ("G") grouping stays on the host.

device_map_bytes(symbolic, INT = eltype of its maps) -> Int64

Bytes of the device maps (DEVICE_MAPS) as INT vectors.

device_maps_supported(backend, structure, view, n, nnz, factor_len) -> Bool

Whether the device symbolic-map routes apply: a GPU backend, a symmetric structure ("SPD", "HPD", "S", "H"), a triangular stored view, and sizes that fit the packed 64-bit sort keys.

device_symmetric_pattern(backend, rowptr, colval, n, view, index) -> (colptr, rowval)

The SymmetricPattern arrays (1-based, sorted, deduplicated, no diagonal) of the stored entries, built on backend. rowptr/colval are host or device vectors with base index.

diagonal_sweep!(ws, symbolic, numeric; nrhs = max_rhs(ws)) -> ws

Diagonal step D W = Z between the sweeps, in place on the first nrhs columns of ws.Y: the identity for Cholesky ("SPD"/"HPD"); for LDLᵀ/LDLᴴ one launch, one work item per (factor column, right-hand side), dividing by the 1×1 pivots of numeric.d and solving the 2×2 blocks (d[k], d[n + k], d[k + 1]; Hermitian or complex symmetric as the structure; LU: the pivots of L D U, transpose is accepted for symmetry). Uniform batch: nrhs per active member, each with its member's D. Asynchronous.

entry_scaling(m::Matching, rowptr, colval, n; index = 'O') -> Vector{Float64}

The factor rᵢ cⱼ of every stored entry (i, j) of the CSR pattern (original column numbering), aligned with colval: the factored values are rᵢ aᵢⱼ cⱼ (scale_values!).

enum_spellings(E) -> Tuple{Vararg{Pair}}

The CUDSS.jl spellings of every instance of the enum type E, as spelling => instance pairs. Spellings are Strings, except for MatrixView, IndexBase and PivotType, which use Chars.

etree(P::SymmetricPattern, perm) -> parent::Vector{Int}

Elimination tree of A[perm, perm] for the symmetric pattern P (Liu's algorithm with path compression). parent[k] is the parent of column k of the permuted matrix, 0 for a root; parent[k] > k for every non-root.

evaluate_ordering(P::SymmetricPattern, perm; T = Float64, level_flops = ORDERING_LEVEL_FLOPS,
                  amalgamation = DEFAULT_AMALGAMATION)
    -> (; nnz_L, flops, nlevels, sdepth, nsupernodes, cost, parent, post, counts)

Elimination tree, column counts, column-etree height (nlevels) and supernodal schedule depth and supernode count (sdepth, nsupernodes, schedule_depth with amalgamation) of A[perm, perm], condensed into the cost of ordering_cost. flops counts real operations (4× for complex T); cost uses the real-arithmetic count whatever T, so the element type does not change the choice. The returned parent, post and counts are the elimination tree, its postorder and the column counts the evaluation computed, so the caller can reuse them (compute_ordering carries the chosen candidate's triple for supernode_partition).

expand_permutation(cperm, ptr, members) -> Vector{Int}

The ordering of the original columns given by an ordering cperm of the compressed vertices of compressed_pattern: the columns of each vertex in turn, partner before candidate.

extend_add!(numeric, symbolic, first, count, maxchild) -> numeric

Owner-pull extend-add: the workgroup of each front s of symbolic.group_nodes[first:(first + count - 1)] adds its children's packed contribution blocks (update stack, cb_ptr) to its panel and its own packed contribution block through relind, children in child_list order (maxchild ≥ the largest child count of these fronts). No atomics; deterministic. Asynchronous.

extract_L(symbolic, numeric) -> SparseMatrixCSC

The Cholesky factor L of P A Pᵀ = L Lᴴ (P = symbolic.partition.perm, supernodal numbering) as a sparse lower-triangular matrix with every stored panel entry, explicit zeros of amalgamation included (nnz == nnz_stored). A device numeric.factor is copied to the host once.

extract_ldlt(symbolic, numeric) -> (L, D, p)

The factors of an LDLᵀ/LDLᴴ factorization (ref_ldlt! or a device numeric, copied to the host once) as sparse matrices with A[p, p] ≈ L D Lᴴ (Lᵀ for complex symmetric "S"): p is the analysis permutation composed with the local pivot orders, L unit lower triangular with every stored panel entry, D block diagonal with 1×1 and 2×2 blocks.

extract_lu(symbolic, numeric) -> (L, D, U, p, q)

The LU factor in numeric (host or device, batch member 1) as host sparse matrices in factor order: L unit lower and U unit upper triangular, D the pivots, p and q the row and column permutations (original indices, p[k] = perm[piv[k]], q = perm), so that A[p, q] ≈ L * Diagonal(D) * U (with the perturbed pivots, A + E). Rows of L below a diagonal block are the rows of later supernodes, which their own local row interchanges moved to ipiv[row].

factor_pattern(P::SymmetricPattern, ord::Ordering) -> SymmetricPattern

The pattern of the symbolic factorization for the ordering ord: pair_pattern(P, ord.pairs).

factorize_cholesky!(numeric, symbolic, nzval; impl = :auto, opts = Options()) -> info::Int

The Cholesky path of factorize! (structures "SPD" and "HPD"), callable without compiling the LDLᵀ/LDLᴴ path (factorize_ldlt!).

factorize_fronts_b!(numeric, symbolic, nzval, first, count, maxchild, W) -> numeric

Regime-B launch: assemble and factor the count fronts symbolic.group_nodes[first:(first + count - 1)] (widths ≤ W, W ∈ (8, 16, 32, 64), at most maxchild children each) of every active batch member with one front_cholesky_kernel! launch, one workgroup per (front, member); info[s] receives each front's status. Asynchronous.

factorize_ldlt!(numeric, symbolic, nzval; impl = :auto, opts = Options(), nb = LDLT_C_NB) -> info::Int

Multifrontal P A Pᵀ = L D Lᴴ (structure "H", or "S" with real T) or L D Lᵀ (complex symmetric "S") on the device, with the in-front pivoting and the static perturbation of the reference ref_ldlt! (same pivot choice, same storage): opts.pivot_type, pivot_threshold, pivot_epsilon, pivot_epsilon_alg and pivot_sign (copied to numeric.psign; a length other than n raises InvalidValueError). Regime-A groups are one subtree_ldlt_kernel! launch per budget class; a regime-B launch group is one front_ldlt_kernel! launch (fused assembly and factorization, one workgroup per front, F₁₁ in local memory above the narrow width classes); regime-C groups and the wide regime-B bins of ldlt_blocked_path take the blocked path of src/numeric/ldlt_c.jl (panel_ldlt_kernel! per block of nb pivot columns, then GEMMs through the dense interface: impl = :auto means vendor GEMMs when the analysis chose vendor calls for regime C, else :ka; other values as in select_impl; nb < LDLT_C_NB is for tests). Then reduce_stats!. opts.user_host_interrupt is polled before every launch group (InterruptedError). With pivot_epsilon_alg = "algo1" abs_max! computes the scale first. Fills numeric.factor (unit-lower panels), d, piv, pivot_kind, stats and totals.

Returns info: an exactly zero 1×1 pivot is always perturbed, and with an effective ε = 0 it stays zero; the factorization then fails there (the following divisions are left to IEEE arithmetic) and info is the original column of the failed pivot with the smallest factor column, as for ref_ldlt!. With ε > 0 the factorization always completes (info = 0). The phase allocates nothing on the device and is deterministic; it synchronizes with the host only to read info when the factorization can fail (pivot_epsilon = 0, or pivot_epsilon_alg = "algo1"), otherwise never.

factorize_lu!(numeric, symbolic, nzval; opts = Options()) -> info::Int

Multifrontal P_r P A Pᵀ = L D U of a general matrix (structure "G") on the device, with the in-front row pivoting and the static perturbation of the reference ref_lu! (same pivot choice, same storage): opts.pivot_type (lu_pivoting), pivot_threshold, pivot_epsilon and pivot_epsilon_alg. Regime-A groups are one subtree_lu_kernel! launch per budget class; every regime-B and regime-C launch group is one front_lu_kernel! launch (fused assembly and factorization, one workgroup per front; no vendor calls); then reduce_stats!. opts.user_host_interrupt is polled before every launch group (InterruptedError). With pivot_epsilon_alg = "algo1" abs_max! computes the scale first. Fills numeric.factor (L), numeric.ufactor (Uᵀ), d, piv (local row order), pivot_kind, stats and totals. Returns info as factorize_ldlt!: an exactly zero pivot stays zero with an effective ε = 0 and the factorization fails there (info = the original column, as ref_lu!); with ε > 0 it always completes (info = 0). The phase allocates nothing on the device, is deterministic, and synchronizes with the host only to read info when the factorization can fail (pivot_epsilon = 0, or pivot_epsilon_alg = "algo1"). Uniform batch as factorize!.

factorize_subtrees!(numeric, symbolic, nzval, first, count, local_bytes) -> numeric

Regime-A launch: assemble and factor the count subtrees symbolic.group_nodes[first:(first + count - 1)] (subtree ids of one budget class) of every active batch member with one subtree_cholesky_kernel! launch, one workgroup per (subtree, member), local_bytes (subtree_local_bytes of the class) of local memory per workgroup. Panels and statuses of every supernode of the subtrees are written, and each subtree root's contribution block goes to the update stack. Asynchronous.

fgmres_available() -> Bool

Whether the Krylov.jl extension is loaded (ir_mode = "fgmres" works).

fgmres_refine!(X, B, W, ws, symbolic, numeric, nzval; nsteps, tol = 0, transposed = false,
               b_transposed = transposed, conjugate = false, deterministic = false,
               interrupt = nothing, progress = Ref(0), transpose_matrix = false, scaling = nothing)
    -> iterations

FGMRES-IR (ir_mode = "fgmres", needs Krylov.jl) of the solution X of op(A) X = B, with the arguments of refine!:

R₀ = B - op(A) X                         (residual!)
stop if tol > 0 and maxₖ ‖R₀ₖ‖₂ / ‖Bₖ‖₂ ≤ tol
D = FGMRES(op(A), R₀; right preconditioner op(A)⁻¹ by the factor,
           at most nsteps iterations, ‖R₀ - op(A) D‖₂ ≤ tol minₖ ‖Bₖ‖₂)
X += D                                   (add_correction!)

on the stacked compact columns (RefinementOperator, FactorPreconditioner). Each iteration costs one SpMV and one solve, as a step of plain refinement. tol = 0 runs nsteps iterations (FGMRES stops earlier only on an exact residual). The joint criterion bounds every column's relative residual by tol. Returns the iterations performed (also in progress[]); an interrupt (polled before every preconditioner application) leaves X unchanged and progress[] = 0. FGMRES synchronizes with the host every iteration (Krylov.jl's dot products and norms).

finish_ldlt_kernel!(backend, WG)(factor, d, pivot_kind, info, nodes, qa, super_ptr, front_ptr, k, nbatch,
                                 front_nrows, front_ncols, Val(HERM), Val(WG); ndrange = WG * count)

After the blocks of the regime-C fronts nodes[qa:(qa + count - 1)], one workgroup each: scale the pivot columns into L (_lt_finalize!), clear the upper triangle of F₁₁, set the unit diagonal and the front's status.

forward_sweep!(ws, symbolic, numeric; nrhs = max_rhs(ws), deterministic = false, impl = :auto,
               transpose = false) -> ws

Forward sweep L Z = Y in place on the first nrhs columns of ws.Y (the permuted right-hand side, permute_rhs!), with the factor of factorize! in numeric and the device maps of symbolic: one launch over all regime-A subtrees, then step by step one launch over the regime-B fronts and, for each regime-C-path front, trsm + gemm through the dense interface (impl as in factorize!) and a scatter kernel. The atomic variant runs unless deterministic is set or the backend lacks an atomic add for T (complex T, issue #36); the deterministic variant goes through the per-front buffers ws.U and an owner-pull of the children (one extra launch per regime-C step). Uniform batch: nrhs right-hand sides per active member of numeric (nrhs × nact columns of ws.Y, each solved with its member's factor; regime-C-path fronts use strided-batched trsm/gemm over runs of consecutive members). LU ("G"): L Z = Y with the local row order, or Uᵀ Z = Y with transpose = true (the solve of Aᵀ). Schur complement mode (schedule.schur > 0): the Schur root is not solved, its rows end as the condensed right-hand side Y₂ − L₂₁ Z₁. Asynchronous; allocates nothing on the device.

front_cholesky!(factor, stack, info, nodes, first, count, front_ptr, front_nrows, front_ncols, cb_ptr;
                width = 64, workgroup = Val(FRONT_WORKGROUP)) -> info

Dense part of the regime-B kernel on already assembled fronts (benchmarks and tests): for each front s = nodes[first + k - 1], k = 1:count, with panel factor[front_ptr[s]:(front_ptr[s + 1] - 1)] (f×w, f = front_nrows[s], w = front_ncols[s] ≤ width, width ∈ (8, 16, 32, 64)) and contribution block stack[cb_ptr[s]:(cb_ptr[s] + m(m+1)/2 - 1)] (m = f - w, packed lower triangle, skipped when cb_ptr[s] == 0): F11 = L11 L11ᴴ, F21 ← F21 L11⁻ᴴ, F22 ← F22 − F21 F21ᴴ, info[s] = status. Asynchronous.

front_cholesky_kernel!(backend, WG)(factor, stack, info, nzval, amap, amap_ptr, amap_src, nodes, bm,
                                    front_ptr, front_nrows, front_ncols, cb_ptr, child_ptr, child_list,
                                    relind_ptr, relind, maxchild, Val(ASM), Val(W), Val(NL), Val(WG);
                                    ndrange = WG * count * bm.nact)

Fused regime-B kernel: workgroup G takes front s = nodes[bm.first + g - 1] (width w ≤ W) of batch member k (g, k from the BatchMap bm), and, when ASM, zeroes and assembles it (A through the amap, then the maxchild children's contribution blocks in child_list order), then factors F11 = L11 L11ᴴ in @localmem (NL = W(W+1)/2 entries), solves F21 ← F21 L11⁻ᴴ, updates its contribution block on the update stack (cb_ptr[s] > 0, packed lower triangle) and sets info[s] (0, or the local column of the first non-positive pivot; the front is then left partially factored). Without ASM the panels and contribution blocks must already be assembled.

front_flops(f, w) -> Int64

Multiply-add count of the partial dense factorization of an f×f front with w pivot columns: Σ_{j=f-w+1}^{f} j² (the elimination of column k updates the trailing (f-k)×(f-k) block). The work measure of the regime-A split.

front_ldlt_kernel!(backend, WG)(factor, stack, info, stats, d, piv, pivot_kind, psign, perm, aux, nzval, amap,
                                amap_ptr, amap_src, nodes, bm, super_ptr, front_ptr, front_nrows,
                                front_ncols, cb_ptr, child_ptr, child_list, relind_ptr, relind, maxchild, prm,
                                Val(W), Val(WG); ndrange = WG * count * bm.nact)

LDLᵀ/LDLᴴ kernel of regime B: workgroup G takes front s = nodes[bm.first + g - 1] of batch member k (g, k from the BatchMap bm; the per-member arrays are the member's), zeroes and assembles it (A through the amap, then the maxchild children's packed contribution blocks in child_list order), factors its w fully-summed columns in place with in-block Bunch–Kaufman pivoting (W = 0, the narrow width classes ≤ _LT_GLOBAL_MAX_W: in the panel in global memory; width class W ≥ w above that: F₁₁ staged as a packed W×W lower triangle in @localmem, the rows below it in the panel; the pivot search is a workgroup reduction), threshold acceptance and perturbation (prm, a _LDLTDevice{R, HERM}; HERM: Hermitian or real symmetric, else complex symmetric), writes D, the local pivot order piv, the pivot kinds and the front's statistics, updates its packed contribution block on the update stack (cb_ptr[s] > 0) with F₂₂ − (L₂₁ D) L₂₁ᴴ and leaves a unit-lower panel.

front_lu_kernel!(backend, WG)(factor, ufactor, stack, ustack, info, stats, d, piv, pivot_kind, aux, nzval, amap,
                              amap_ptr, amap_src, nodes, bm, super_ptr, front_ptr, front_nrows, front_ncols,
                              cb_ptr, child_ptr, child_list, relind_ptr, relind, maxchild, prm, Val(WG);
                              ndrange = WG * count * bm.nact)

LU kernel of regimes B and C: workgroup G takes front s = nodes[bm.first + g - 1] of batch member k (g, k from the BatchMap bm; the per-member arrays are the member's), zeroes and assembles its L and Uᵀ structures (A through the amap, then the maxchild children's two packed contribution blocks in child_list order), factors its w fully-summed columns in place in the two panels with in-block row pivoting, threshold test and perturbation (prm, an _LUDevice), writes D, the local row order piv, the pivot kinds and the front's statistics, updates its two packed contribution blocks (cb_ptr[s] > 0) with F₂₂ − (L₂₁ D) U₁₂ and leaves unit-lower panels of L and Uᵀ.

full_pattern_map(A::CSR, structure; view = 'F') -> FullPatternMap
full_pattern_map(rowptr, colval, n, structure; view = 'F', index = 'O') -> FullPatternMap

Build the FullPatternMap of a CSR input with the same structure and view rules as SymmetricPattern (for "G" the matrix is taken as is). Duplicated user entries are summed.

full_values(F::FullPatternMap, nzval::AbstractVector) -> Vector

Host evaluation of the map: the values of the full matrix M in the order of F.colval (see FullPatternMap).

fundamental_supernodes(parent, post, counts) -> ColumnPartition

Fundamental supernodes of the elimination tree parent with column counts counts (both in the etree numbering, from etree and colcounts) and its postorder post (postorder). The columns are renumbered by post (order == post); column k + 1 joins the supernode of column k when k is its only child and counts[k] == counts[k+1] + 1, so the row structure of every supernode is exactly that of its first column.

gemm!(C, A, B, α = 1, β = 0; transA = :N, transB = :N, impl = :auto) -> C

C ← α op(A) op(B) + β C with op from :N, :T, :C.

gemm_strided_batched!(C, A, B, α = 1, β = 0; transA = :N, transB = :N, impl = :auto) -> C

Uniform strided batch of gemm!: C[:, :, i] ← α op(A[:, :, i]) op(B[:, :, i]) + β C[:, :, i] for every member i of the 3-D arrays (see strided_batch). Implementations: :vendor, :ka.

getrf!(A, ipiv; impl = :auto) -> info::Int

LU factorization with partial pivoting of the m × n matrix A in place, P A = L U (unit L), LAPACK getrf: ipiv[1:min(m, n)] receives the row interchanges (1-based, see laswp!); info is 0 or the first column with an exactly zero pivot. ipiv lives on the backend of A. Reading info synchronizes with the device.

herk!(C, A, α, β; uplo = :L, impl = :auto) -> C

Hermitian rank-k update of the uplo triangle of C: C ← α A Aᴴ + β C with real α, β (BLAS herk; the diagonal of the result is real). For real T this is syrk!. Triangle semantics as in syrk!.

host_numeric(numeric) -> Numeric{T, Vector{T}}

Host copy of numeric (panels, D, update stack, workspace, statistics and totals, status, pivot order, kinds and sign requests, aux, the Uᵀ panels and stack of LU; the plan is shared), usable by ref_solve! and extract_L. A host numeric is returned as is.

induced_pattern(P::SymmetricPattern, vertices) -> SymmetricPattern

The pattern of P restricted to vertices (sorted, distinct), renumbered 1:length(vertices) in their order.

inertia(numeric) -> (npos, nneg)

Inertia of D of an LDLᵀ/LDLᴴ factorization (PLAN §1.4 inertia). With perturbed pivots this is the inertia of A + E; read it together with npivots.

ka_chol_check_info!(info, offset, L::AbstractArray{T,3}) -> info

Batched form: the check of ka_chol_check_info!(info, offset + b, L[:, :, b]) for every member b of the 3-D strided batch L, one launch.

ka_chol_check_info!(info, idx, L) -> info

Device-side validation of the info[idx] a vendor or generic Cholesky left for the computed factor L (a matrix): when info[idx] == 0, it is replaced by the first j whose diagonal entry L[j, j] is not a finite positive real (still 0 if there is none); a nonzero info[idx] is kept. Asynchronous (see potrf_info!).

ka_chol_diag_info!(info, L) -> info

Pivot check of a computed Cholesky factor L (matrix or 3-D strided batch): info[b] is set to the first j whose diagonal entry L[j, j] is not a finite positive real, 0 if there is none. Used to validate the info = 0 of vendor and generic factorizations (see potrf!). Asynchronous.

ka_gemm!(C, A, B, α, β; transA = 'N', transB = 'N', tile = Val(16) or Val(32)) -> C

KernelAbstractions fallback for C ← α op(A) op(B) + β C (op from 'N', 'T', 'C'). Tiled, one TILE × TILE tile of C per 1-D workgroup of TILE² work items; tile defaults to Val(32) when both dimensions of C are at least 64. C is not read when β == 0. Asynchronous: no host synchronization.

ka_gemm_strided_batched!(C, A, B, α, β; transA = 'N', transB = 'N', tile) -> C

Strided batched ka_gemm!: A, B, C are 3-D arrays (for example strided_batch views of a flat buffer) whose third dimension is the batch; one launch, the batch index is the second ndrange dimension.

ka_getrf!(A, ipiv, info; workgroup = Val(128)) -> info

KernelAbstractions fallback LU with partial pivoting of the m × n matrix A in place (P A = L U, unit L below the diagonal, U on and above), one 1-D workgroup per matrix. A is a matrix (ipiv a vector of length min(m, n)) or a 3-D strided batch (ipiv a min(m, n) × count matrix). info is a device Int32 vector with one entry per matrix: 0, or the first column with an exactly zero pivot, as LAPACK getrf. Asynchronous.

ka_laswp!(A, ipiv; npiv = size(ipiv, 1), reverse = false, workgroup = 128) -> A

KernelAbstractions row interchanges: for i = 1:npiv (in reverse order if reverse), swap rows i and ipiv[i] of A, as LAPACK laswp with k1 = 1, k2 = npiv; entries of ipiv beyond npiv are not read. Applied after ka_getrf! this computes P A; reverse = true applies Pᵀ. One work item per column; A may be a 3-D strided batch with ipiv a matrix (one column per member).

ka_potrf!(uplo, A, info; workgroup = Val(128), offset = 0) -> info

KernelAbstractions fallback Cholesky of the uplo triangle of A in place (A = LLᴴ for 'L', UᴴU for 'U'), unblocked, one 1-D workgroup per matrix. A is a matrix or a 3-D strided batch; info is a device vector of Int32 with one entry per matrix, set to 0 on success or to the first column j (1-based) whose pivot is not positive, as LAPACK potrf; matrix b of the batch writes info[offset + b]. Asynchronous: the caller reads info at a phase boundary.

ka_syrk!(C, A, α, β; uplo = 'L', conjugate = false, tile) -> C

KernelAbstractions fallback for the uplo triangle of C ← α A Aᵀ + β C (conjugate = true: α A Aᴴ + β C with real α, β and a real diagonal, as BLAS herk). The other triangle is not touched.

ka_trsm!(side, uplo, trans, diag, α, A, B; workgroup = Val(128)) -> B

KernelAbstractions fallback for BLAS trsm: B ← α op(A)⁻¹ B (side = 'L') or B ← α B op(A)⁻¹ (side = 'R'), A triangular (uplo), unit diagonal if diag = 'U', op from trans ∈ ('N', 'T', 'C'). One 1-D workgroup per right-hand side column (row for side = 'R') runs the forward/back substitution. Asynchronous.

ka_trsm_strided_batched!(side, uplo, trans, diag, α, A, B; workgroup = Val(128)) -> B

Strided batched ka_trsm! on 3-D arrays (batch = third dimension, the second ndrange dimension of the launch).

laswp!(A, ipiv; npiv = length(ipiv), reverse = false, impl = :auto) -> A

Apply the row interchanges of ipiv (as returned by getrf!) to A: for i = 1:npiv swap rows i and ipiv[i], so that laswp!(A, ipiv) computes P A; reverse = true applies them backwards (Pᵀ A). npiv is LAPACK's k2 (with k1 = 1): pass npiv = min(m, n) when ipiv is an oversized buffer filled by getrf!. Only :ka exists (no generic entry point; vendor bindings may be added per backend).

ldlt_blocked_path(sc, s) -> Bool

Front s takes the blocked LDLᵀ/LDLᴴ path (panel_ldlt_kernel! and GEMMs, src/numeric/ldlt_c.jl): it is on the regime-C path (takes_c_path), or in a regime-B bin of width class ≥ LDLT_BLOCKED_MIN_WCLASS and row class ≥ LDLT_BLOCKED_MIN_FCLASS, where the GEMMs beat the fused kernel's in-kernel updates.

ldlt_c_chunk_end(sc, cb_len, q, last, cap) -> Int

The fronts sc.group_nodes[q:e] of a regime-C launch group (e ≤ last) that the blocked LDLᵀ/LDLᴴ path factors concurrently, one workgroup each: the longest run from q whose ldlt_c_work_len slices, laid out one after another, fit cap entries (at least one front). cb_len[s] > 0 when front s has a contribution block on the update stack.

ldlt_c_work_len(f, w, cb::Bool, nb = LDLT_C_NB) -> Int

Workspace entries of a regime-C LDLᵀ/LDLᴴ front with f rows and w fully-summed columns: the m×m contribution block (m = f - w, when cb), Lb and Wb (f × (nb + 1) each), the saved column, the next pivot column of pass 1 and the rejection flags of the fallback scan (f each).

lu_pivoting(pivot_type::PivotType) -> Bool

Whether the LU of structure "G" interchanges rows inside the fully-summed block of a front: 'A' (default), 'L' and 'B' mean in-block threshold partial pivoting; 'N' and 'D' keep the diagonal pivots (perturbed when tiny). Global pivoting ('C', 'R') is not planned (PLAN §3.3).

matched_colval(m::Matching, colval, index) -> Vector

The column indices of the CSR pattern of the matched matrix M = (Dr A Dc)[:, q] of "G" (q = m.perm): column c of the stored matrix becomes column q⁻¹[c], in the index base index and the element type of colval.

matching_job(alg::MatchingAlg) -> Int

The MC64 job of matching_alg (PLAN §1.3): "algo1"–"algo5" are jobs 1–5, "algo6" (auto) is job 5, "default" is 0 (no matching).

matching_pairs(m::Matching, rowptr, colval, nzval, n, structure, opts; view = 'F', index = 'O')
    -> Vector{Tuple{Int,Int}}

The 2×2 pivot pairs (partner, candidate) of "S"/"H" derived from the symmetric matching m (issue #67, Duff & Pralet 2005), sorted by candidate, for compute_ordering(…; pairs):

• the matching permutation is split into cycles (only cycles of matched indices); a 2-cycle (i, j) is a pair; a longer cycle is split into pairs of consecutive indices (adjacent in the matrix, since aᵢ,perm[i] ≠ 0): for an even length the alternative with the larger product of scaled |dᵢ aᵢⱼ dⱼ|, for an odd length leaving out the index with the largest scaled diagonal;
• opts.pivot_pairs = "default" keeps a pair only when one of its indices is a 2×2 candidate of the scaled matrix: |dᵢ² aᵢᵢ| ≤ τ maxⱼ≠ᵢ |dᵢ aᵢⱼ dⱼ|, τ = opts.pivot_pair_tolerance; "all" keeps every pair; "none" none;
• the candidate of a pair is its index with the smaller relative scaled diagonal, the other one the partner (eliminated first).

max_local_bytes(backend) -> Int

Largest static @localmem size (bytes) one workgroup of a KernelAbstractions kernel may declare on backend (issue #60). The regime-A subtree kernels pick their local-memory class from SUBTREE_LOCAL_SIZES up to this cap: the analysis clamps the default subtree_budgets to it and rejects an explicit budget above it with an InvalidValueError.

49152 (48 KiB) for a backend without a method, the static shared-memory limit of CUDA; the CPU backend has no limit (typemax(Int)). The backend extensions add their own methods (CUDA: 48 KiB; AMDGPU: the device's LDS size per workgroup, 64 KiB on CDNA).

max_rhs(ws::SolveWorkspace) -> Int

Number of right-hand sides ws was allocated for.

max_threads_per_multiprocessor(backend) -> Int

Resident work items per compute unit of the device backend runs on, 0 when the backend does not report it. The CUDA extension reads the device attribute (2048 on sm70/80/90, 1536 on sm86/89/120); the fused pipeline's per-architecture tuning (Fused.default_tuning) keys on it.

mc64(ptr, adj, w, n, job) -> (lmate, rscale, cscale)

MC64-style matching of the rows 1:n of a square matrix (CSR ptr/adj, absolute values w) to its columns:

• job 1: maximum cardinality (structural; every stored entry counts);
• jobs 2 and 3: maximize the smallest matched |aᵢⱼ| (bottleneck);
• job 4: maximize the sum of the matched |aᵢⱼ|;
• job 5: maximize their product, with the scalings rᵢ = exp(uᵢ) / maxⱼ |aᵢⱼ|, cⱼ = exp(vⱼ) from the duals of the assignment on the costs log maxₖ |aᵢₖ| - log |aᵢⱼ|, so that |rᵢ aᵢⱼ cⱼ| ≤ 1, with equality on the matched entries.

Jobs 2–5 use the nonzero entries only. lmate[i] is the column matched to row i (0 for an unmatched row of a structurally singular matrix); the scalings are ones for jobs 1–4.

member_info!(info, numeric, symbolic) -> info

Write the info of the last factorize! of every active batch member k into info[k] (the original column of the member's first failed pivot: non-positive for Cholesky, zero with pivot_epsilon = 0 for LDLᵀ/LDLᴴ and LU; 0 = success; inactive members keep their entries). Host data only: the statuses were read by factorize!.

member_numeric(numeric, symbolic, k) -> Numeric{T, Vector{T}}

Host copy of batch member k of numeric as a single-matrix storage (the layout of allocate_numeric(symbolic, T)), for the reference comparisons and extract_L.

memory_estimates(symbolic, T[, INT]) -> Vector{Int64}

The "memory_estimates" data parameter (16 entries, bytes) for factors of element type T and device indices INT (default: the index type of symbolic's maps, Int for a host analysis):

Slots 1–6 follow cuDSS's CUDSS_DATA_MEMORY_ESTIMATES. getparam(solver, "memory_estimates") adds to slots 1, 2 and 5 the device buffers of the matching when there is one (matching_alg): the scaled values (nnz · nbatch of T), the entry weights (nnz) and the two scale vectors (n each) of real(T), and for "G" the column permutation (n of INT).

min_cost_matching(ptr, adj, cost, n) -> (lmate, u, v)

Minimum-cost matching of maximum cardinality between the left vertices 1:n (adjacency adj[ptr[i]:ptr[i+1]-1], costs cost ≥ 0 aligned with adj) and the right vertices 1:n, by successive shortest augmenting paths (Dijkstra on the reduced costs cost - u[i] - v[j] ≥ 0, a binary heap; the sparse Hungarian method of MC64). Returns lmate[i] (the right vertex matched to i, 0 if none) and the dual variables u (left) and v (right): u[i] + v[j] ≤ cost on every edge, with equality on the matched ones.

nd_algorithm(nd_nlevels, nd_ubfactor, nd_nseps, nd_seed) -> CliqueTrees.METIS

The CliqueTrees nested-dissection algorithm object of the nd_* parameters: METIS_NodeND with ufactor = nd_ubfactor, nseps = nd_nseps and seed = nd_seed (-1: METIS defaults; nseps < 1 and seed < 0 fall back to the defaults). nd_nlevels is cuDSS' minimum number of dissection levels; METIS_NodeND dissects recursively until the parts are small, so the minimum holds whenever the graph is large enough to be split that often. CliqueTrees' level-capped ND{S} (AMD below the cap) was 4–13× slower and gave 15–80% more fill on Laplacians (T05 Report).

neighbors(P::SymmetricPattern, j)

Sorted indices i ≠ j adjacent to j (a view into P.rowval).

nlaunches(schedule) -> Int

Number of kernel and vendor launches of one numeric factorization following schedule: one per regime-A budget class, one per regime-B group (fused assemble + factor + update kernel), and per regime-C group one batched assembly/extend-add launch plus, per front, potrf and, when the front has a contribution block, trsm, syrk/herk and the pack_add! of the workspace product into the packed contribution block.

nnz_L(counts) -> Int

Number of nonzeros of the Cholesky factor L (diagonal included): sum(counts).

npivots(numeric) -> Int64

Number of perturbed pivots of an LDLᵀ/LDLᴴ factorization (PLAN §1.4 npivots).

nsubtrees(schedule) -> Int

Number of regime-A subtrees.

nsupernodes(cp) -> Int

Number of supernodes of a ColumnPartition or SupernodePartition.

nsuperpanels(sp::SupernodePartition) -> Int

Number of supernodes after amalgamation (the "nsuperpanels" data parameter).

ordering_cost(flops, sdepth, level_flops = ORDERING_LEVEL_FLOPS) -> Float64

Cost model of the automatic ordering choice (PLAN §2.3 step 2, issue #108): the time of one factorization and its solves, in flops,

cost = flops + level_flops × sdepth.

Terms. flops is the real-arithmetic Cholesky operation count of the candidate (cholesky_flops), the throughput-bound part. sdepth is the depth of its supernodal schedule (schedule_depth): the numeric phase launches the fronts level by level and both solve sweeps walk the same levels, so every level costs a latency (launches, synchronization, a tail of few fronts) whatever its size, and the solve, nearly all latency, scales with the depth alone. level_flops is that latency times the throughput.

Why the supernodal schedule and not the column etree. The column-etree depth (nlevels) counts every column of a chain, and chains collapse into supernodes, so it does not predict the schedule depth (issue #108, 78k-bus condensed KKT: cuDSS's ordering 1279 column levels → 30 schedule levels, METIS 1216 → 29, AMD a comparable column depth → 65). The fundamental partition does not predict it either: on the random KKT generators it ranks AMD 33 vs ND 63 where the amalgamated schedules both have 12 levels, and it ranks AMD ahead of ND on apache2 (48 vs 57) where the amalgamated schedule ranks ND ahead (35 vs 39). So sdepth is the height of the amalgamated supernodal tree, the nlevels of the Schedule the analysis builds.

Weight. The term is absolute, not relative to n (the previous model, flops × (1 + nlevels / n), weighed depth at 0.2% at n ≈ 7e5 and was a flop contest). level_flops = 1e8 (ORDERING_LEVEL_FLOPS) is ~100 µs per level at ~1 Tflop/s: one level of the level-batched factorization (several launches) plus a forward and a backward sweep for each of a few solves per factorization (refinement), as measured on the 78k-bus KKT in PR #107 (ND's 36 fewer levels: solve −36%, factorization −11 to −17% despite AMD's lower flop count). Consequences: matrices below ~1e10 flops choose the shallower schedule, unless the depths are within a few levels; 3-D problems, where ND saves 10¹⁰ flops and more, choose the fewer flops. One weight serves every backend, so that the CPU and GPU backends factor the same ordering. Ties keep the first candidate (AMD).

pack_add!(stack, c0, work, m) -> stack
pack_add!(stack, c0, work, m, members, j0, count, nbatch) -> stack

Add the lower triangle of the column-major m×m matrix work[1:m^2] to the packed contribution block stack[c0:(c0 + m(m+1)/2 - 1)] (the regime-C step after the vendor syrk/herk into the workspace). Batched form: the count blocks work[(b-1)m² + 1 : b m²] go to the stacks of batch members members[j0 + b] (device vector) of an nbatch-member update stack. One launch, one work item per entry. Asynchronous.

packed_length(m) -> Int

Entries of the packed lower triangle of an m×m matrix, m(m+1)/2: the storage of a contribution block on the update stack and of a regime-A front in local memory (column-major packed, see _packed).

pair_pattern(P::SymmetricPattern, pairs) -> SymmetricPattern

The pattern the symbolic factorization runs on when pairs is not empty (P itself otherwise): the expansion of compressed_pattern, where both columns of a pair have the union structure of the pair (explicit structural zeros of the factor) and are adjacent to each other. With the partner ordered right before its candidate, the candidate is the only child of the partner in the elimination tree and their columns have the same structure below the pair, so the pair is (part of) one fundamental supernode; amalgamation only merges supernodes, so the pair is never split.

pairs_enabled(structure, opts::Options) -> Bool

Whether the analysis looks for 2×2 pivot candidate pairs: opts.pivot_pairs is "default" or "all", structure "S" or "H", no user_perm and not the natural ordering (reordering_alg = "algo5").

panel(symbolic, numeric, s, k = 1) -> Matrix

Host copy of the f×w panel of supernode s (rows snrows(s), the factor columns of s; the strict upper triangle of its first w rows is zero) of batch member k.

panel_ldlt_kernel!(backend, WG)(factor, lw, d, piv, pivot_kind, psign, perm, aux, stats, nodes, qa, b, nb,
                                super_ptr, front_ptr, k, nbatch, front_nrows, front_ncols, cb_ptr, prm, Val(WG);
                                ndrange = WG * count)

Block b (columns (b - 1) nb + 1 : min(b nb, w)) of the regime-C LDLᵀ/LDLᴴ pivot steps of the count fronts nodes[qa:(qa + count - 1)], one workgroup each (panel of batch member k of nbatch at member_panels(front_ptr, k, nbatch)[s], workspace slice in lw after those of the fronts before it: Lb, Wb, the saved column, the pass-1 column): the reference's pivot sequence on the lazily updated panel (see the top of src/numeric/ldlt_c.jl), D, piv, the pivot kinds and the front's statistics (added to those of the earlier blocks). Fronts with fewer blocks do nothing.

panel_offset(panel_ptr, s, k, nb) -> Int

First factor index of the panel of supernode s of batch member k in the interleaved layout (host side, panel_ptr of the Layout).

parameter_spec(name) -> ParameterSpec

Look up a parameter name; unknown names raise ArgumentError.

permute_rhs!(Y, B, perm; transposed = false, conjugate = false, bm = nothing, scale = nothing) -> Y

Y[k, r] = B[perm[k], r] (conjugated when conjugate, for the solves with conj(A) of solve_mode; times scale[perm[k]] when a real scale vector is given, the row or column scaling of matching, SolveScaling) for the nrhs right-hand sides of the user array B (rhs_count: vector, matrix, strided vector; row-major when transposed) into the first nrhs columns of the n × ≥ nrhs device matrix Y. With a BatchMap bm (uniform batch), B holds bm.nrhs right-hand sides for each of the bm.nbatch members and only those of the active members are copied, member slot j to the columns (j - 1) bm.nrhs + 1 : j bm.nrhs of Y. One launch; asynchronous.

phase_includes(phase::Phase, part::Phase) -> Bool

Whether the bits of part are all set in phase, e.g. phase_includes(PHASE_SOLVE, PHASE_SOLVE_DIAG) == true.

pivot_candidates(P::SymmetricPattern, rowptr, colval, nzval, n, structure; view = 'F', index = 'O')
    -> PairCandidates
pivot_candidates(P::SymmetricPattern, A::CSR, structure; view = 'F') -> PairCandidates

The 2×2 pivot candidates (PairCandidates) of the CSR matrix (rowptr, colval, nzval; the first batch member when nzval is a matrix) whose symmetric pattern is P, read with the view and index rules of SymmetricPattern (duplicates summed), with the candidate tolerance τ = tolerance (the analysis passes opts.pivot_pair_tolerance). The only read of the values at analysis.

pivot_pairs(P::SymmetricPattern, rowptr, colval, nzval, n, structure; view = 'F', index = 'O')
    -> Vector{Tuple{Int,Int}}
pivot_pairs(P::SymmetricPattern, A::CSR, structure; view = 'F')
pivot_pairs(P::SymmetricPattern, cands::PairCandidates)

The pairs of pivot_pairs = "all": every 2×2 pivot candidate of pivot_candidates matched to a neighbour, as (partner, candidate), sorted by candidate.

• Pairing: a greedy matching of candidates to non-candidate neighbours with aᵢⱼ ≠ 0. Candidates are taken by increasing number of such neighbours (then by index); each takes the free neighbour with the largest |aᵢⱼ|, ties by smaller degree in P, then smaller index. Each non-candidate is matched at most once; unmatched candidates stay single columns.

On KKT systems about half of the rows are candidates, and ordering with all of them paired roughly doubles nnz(L) (issue #66); "default" pairs only the candidates that need it, see zero_pivot_pairs!. Does not look at structure beyond the view rules; pairs_enabled decides whether the analysis calls it.

pivot_stats(numeric) -> (npos, nneg, nzero, nperturbed, n2x2)

Sum of the per-front statistics of an LDLᵀ/LDLᴴ factorization (PLAN §1.7 pivot_stats): positive and negative eigenvalues of D (after perturbation; a 2×2 block with negative determinant counts one of each), pivots that were exactly zero, perturbed pivots, 2×2 blocks. A device numeric.stats is copied to the host. npos/nneg are 0 for complex symmetric ("S") matrices.

pivot_totals(numeric, k = 1) -> (npos, nneg, nzero, nperturbed, n2x2)

The statistics of batch member k reduced on the device by reduce_stats! (one copy of numeric.totals to the host: a synchronization, done at the phase boundary by getparam).

postorder(parent) -> post::Vector{Int}

A postorder of the forest parent (0 = root): post[k] is the node visited k-th; children are visited in increasing index order and every node comes after all its descendants. Non-recursive.

potrf!(uplo, A; impl = :auto) -> info::Int

Cholesky factorization of the uplo triangle of the square matrix A in place (A = LLᴴ or UᴴU). Returns 0 on success, or the first column j (1-based) whose leading minor is not positive definite, as LAPACK potrf. Reading info synchronizes with the device.

potrf_batched_info!(uplo, A, info, offset = 0; impl = :auto) -> info

Uniform strided batch of potrf_info!: Cholesky of the uplo triangle of every member A[:, :, b] of the 3-D array A in place, its status (0, or the first non-positive pivot column) into info[offset + b]. Implementations: :vendor (vendor_potrf_batched!, validated on the device as potrf_info!) and :ka (ka_potrf!, one workgroup per member). Asynchronous.

potrf_info!(uplo, A, info, idx = 1; impl = :auto) -> info

Asynchronous potrf!: Cholesky of the uplo triangle of A in place, with the LAPACK info (0, or the first non-positive pivot column) written to info[idx] of the device Int32 vector info instead of being returned, so the numeric phase can factor many fronts and read their status once (PLAN §3.9). :vendor (vendor_potrf_info!) and :generic results are validated on the device by ka_chol_check_info! (cuSOLVER can miss a non-positive last pivot, see potrf!); :ka is ka_potrf!. The :vendor and :ka paths do not synchronize with the host; :generic does when the backend's cholesky! reads its info (CUDA).

print_capabilities([io = stdout], backend; eltypes = (Float32, Float64, ComplexF32, ComplexF64))

Print the capabilities table of backend: one row per probe, one column per element type.

reduce_stats!(numeric, symbolic) -> numeric

Sum the per-front statistics numeric.stats into numeric.totals (npos, nneg, nzero, nperturbed, n2x2, and the number of fronts with a failed pivot) of every active batch member, and the smallest failed factor column (super_ptr[s] + info - 1 over the fronts with a nonzero info statistic, 0 = none) into numeric.info[ns + 1]: one workgroup per member, @localmem tree reduction, no atomics. Asynchronous; pivot_totals reads the result.

ref_factorize!(numeric::Numeric{T, Vector{T}}, symbolic, nzval) -> info

Reference multifrontal Cholesky P A Pᵀ = L Lᴴ of the matrix whose stored values are nzval (same CSR pattern, view and index base as the analysis; a host vector or anything Array converts) into the host numeric. The supernodes are processed in snpost order: zero the f×f front, add A through amap_ptr/amap_src/amap (conjugating negative offsets), extend-add the children's contribution blocks through relind in child_list order, then potrf on the w×w block, trsm on the (f-w)×w block and syrk/herk on the contribution block. The front's first w columns are stored in its panel.

Returns info = 0 on success, else the first non-positive pivot (the smallest failed column of the factor, as the column of the original matrix); the factorization stops there and later panels stay zero. Requires structure "SPD" (real T) or "HPD".

ref_ldlt!(numeric::Numeric{T, Vector{T}}, symbolic, nzval; opts = Options()) -> info

Reference multifrontal LDLᵀ (real T, or complex T with structure "S", complex symmetric) or LDLᴴ (structure "H") of the matrix whose stored values are nzval (same CSR pattern, view and index base as the analysis), into the host numeric. Assembly and extend-add are those of ref_factorize!; each front is then factored densely with pivots chosen inside its fully-summed block (its first w columns), and the contribution block is updated once with F₂₂ ← F₂₂ − (L₂₁ D) L₂₁ᴴ (ᵀ for complex symmetric).

Pivoting (opts.pivot_type; 'A', 'B' and 'L' are Bunch–Kaufman, the default):

• 'B': Bunch–Kaufman with α = (1+√17)/8 on the remaining block, choosing a 1×1 pivot at the current column, a 1×1 pivot at the row r of the largest block entry of that column, or the 2×2 pivot on both. The choice is acceptable when it is not tiny (|d| ≥ ε, or for a 2×2 block a nonzero determinant and an entry ≥ ε) and passes the threshold test with u = opts.pivot_threshold against the whole remaining column of the front, fully-summed or not (|d| ≥ u maxᵢ |aᵢc|; for a 2×2 block, every entry of the resulting L columns bounded by 1/u through |D⁻¹|). Otherwise the acceptable 1×1 pivot of the block with the largest |aⱼⱼ| is taken; if there is none, the Bunch–Kaufman choice is used anyway when it is not tiny (stable inside the block, only the bound on the rows outside fails); else the column is perturbed.
• 'D': 1×1 pivots only: the current column if acceptable, else the acceptable diagonal with the largest |aⱼⱼ|, else the largest diagonal if not tiny, else the current column, perturbed.
• 'N': no search: the current column, perturbed when tiny.

Perturbation: a 1×1 pivot with |d| < ε becomes d ← sign · ε, where ε is pivot_epsilon (default_pivot_epsilon if unset), multiplied by max |aᵢⱼ| of nzval when pivot_epsilon_alg = "algo1" (scaled; "default" and "algo2" are static), and sign is opts.pivot_sign[row] for the original row of that pivot when given and nonzero, else the sign of d (+1 for d == 0; for complex symmetric T, d/|d|).

Fills numeric.factor (unit-lower panels), numeric.d, numeric.piv, numeric.pivot_kind and the per-front numeric.stats (npos, nneg, nzero, nperturbed, n2x2, info; nzero counts pivots that were exactly zero before their perturbation; npos/nneg count the signs of D after perturbation, a 2×2 block with negative determinant counting one of each, and stay 0 for complex symmetric matrices, which have no inertia). An exactly zero 1×1 pivot is always tiny; with ε = 0 its perturbation leaves it zero and the factorization fails there: the front's info statistic is the failed local column (assembled order, the first in pivot order), the following divisions are left to IEEE arithmetic (the factor is unusable), and the return value is the original column of the failed pivot with the smallest factor column, as for the Cholesky oracle. Otherwise returns info = 0: with static perturbation ε > 0 the factorization always completes.

ref_lu!(numeric::Numeric{T, Vector{T}}, symbolic, nzval; opts = Options()) -> info

Reference multifrontal L D U of the general matrix whose stored values are nzval (structure "G", same CSR pattern and index base as the analysis) into the host numeric. The supernodes of the symmetric pattern of A + Aᵀ are processed in postorder: the full f×f front is assembled (A through the amap, whose negative offsets are the entries above the diagonal; the children's full contribution blocks through relind), its w fully-summed columns are factored with row interchanges inside the fully-summed rows, and the contribution block is updated once with F₂₂ ← F₂₂ − (L₂₁ D) U₁₂.

Pivot choice at step k (opts.pivot_type, see lu_pivoting): with γ = maxᵢ |F[i, k]| over all remaining rows of the front (fully-summed or not), the diagonal is kept when |F[k, k]| ≥ u γ (u = pivot_threshold), otherwise the fully-summed row r ≥ k with the largest |F[r, k]| (the first one) is swapped into row k; it is taken even when it fails the threshold (static pivoting: no delayed pivots). A pivot with |d| < ε becomes d ← sign(d) ε (+ε for d = 0; d/|d| for complex T), with ε as in ref_ldlt! (pivot_epsilon, scaled by max |aᵢⱼ| with pivot_epsilon_alg = "algo1"); pivot_sign is not used by LU.

Fills numeric.factor (unit-lower L panels), numeric.ufactor (unit-lower Uᵀ panels), numeric.d, numeric.piv (local row order), numeric.pivot_kind (PIVOT_KIND_1X1 or PIVOT_KIND_PERTURBED) and the per-front numeric.stats (nzero, nperturbed; npos, nneg and n2x2 stay 0: a general matrix has no inertia). An exactly zero pivot is always tiny; with ε = 0 it stays zero and the factorization fails there: the front's info statistic is the failed local column and the return value the original column of the failed pivot with the smallest factor column (as ref_ldlt!). Otherwise returns 0: with static perturbation ε > 0 the factorization always completes.

ref_solve!(X, symbolic, numeric, B) -> X

Solve A X = B with the reference factor P A Pᵀ = L Lᴴ of ref_factorize!: Y = B[perm, :], forward sweep L Z = Y over the supernodes in snpost order (trsm on the diagonal block, gemm update of the rows below), backward sweep Lᴴ W = Z in reverse order, X[perm, :] = W. B and X are host vectors or n × nrhs matrices (X === B is allowed). For structures "S"/"H" this is the LDLᵀ/LDLᴴ solve of ref_solve_ldlt! with a factor of ref_ldlt!, for "G" the LU solve of ref_solve_lu! with a factor of ref_lu!.

ref_solve_ldlt!(X, symbolic, numeric, B) -> X

Solve A X = B with the reference factor of ref_ldlt!: Y = B[perm, :]; forward sweep over the supernodes in snpost order (reorder the supernode's slice of Y by its local pivot order piv, unit-lower trsm, update of the rows below); diagonal step Y ← D⁻¹ Y (1×1 and 2×2 blocks); backward sweep in reverse order (update from the rows below, trsm with L₁₁ᴴ (ᵀ for complex symmetric), undo the local pivot order); X[perm, :] = Y. Called by ref_solve! for structures "S"/"H". X === B is allowed.

ref_solve_lu!(X, symbolic, numeric, B; transpose = false) -> X

Host solve of A X = B (Aᵀ X = B when transpose) with the LU factor in numeric through extract_lu (sparse triangular solves; tests and debugging). B, X: vectors or matrices.

refine!(X, B, W, ws, symbolic, numeric, nzval; nsteps, tol = 0, transposed = false,
        b_transposed = transposed, conjugate = false, deterministic = false,
        interrupt = nothing, progress = Ref(0), transpose_matrix = false, scaling = nothing) -> steps

Plain iterative refinement of the solution X of op(A) X = B (op(A) = M, or conj(M) when conjugate; Mᵀ or conj(Mᵀ) with transpose_matrix, LU only), at most nsteps steps of

R = B - op(A) X            (residual!)
stop if tol > 0 and maxₖ ‖Rₖ‖₂ / ‖Bₖ‖₂ ≤ tol   (residual_norms!, one host synchronization)
X += op(A)⁻¹ R             (permute, forward / diagonal / backward sweeps, add_correction!)

with the factor in numeric, the solve workspace ws and the refinement workspace W (with the permutations and scalings of matching when scaling is a SolveScaling, solve_io; the residual is always the one of the original matrix values nzval). tol = 0 (the default of ir_tol) never stops early and never synchronizes. interrupt (a Threads.Atomic{Bool} or nothing) is polled before every step (InterruptedError; X then holds the last completed iterate). Returns the number of corrections applied, which is also kept in progress[] after every step (so it survives an interrupt).

refinement_map(F::FullPatternMap; transpose = false) -> (rowptr, colval, src)

Host CSR of the full matrix of F over its contributions (duplicates are not merged; their products are summed by the SpMV): row i holds the entries rowptr[i]:rowptr[i+1]-1, column colval[c], value nzval[src[c]], conjugated when src[c] < 0. 1-based. With transpose = true (structure "G") the rows of Mᵀ follow as rows n+1:2n (rowptr has 2n + 1 entries).

relative_indices(sp::SupernodePartition) -> (relind_ptr, relind)

Extend-add maps: for every supernode c with a parent p the positions of snrows(sp, c)[w_c+1:end] inside snrows(sp, p) (they are a subset, both sorted), concatenated; empty for roots and fronts without a contribution block.

residual!(W, nzval, X, B; nrhs, transposed = false, b_transposed = transposed, conjugate = false,
          bm = nothing, transpose_matrix = false) -> W.R

W.R[:, 1:nrhs] = B - M X (conj(M) when conjugate; Mᵀ when transpose_matrix, which needs a "G" map) with the full matrix M of the map of W and the user values nzval; X, B in the layouts of rhs_count (row-major when transposed/b_transposed) with nrhs right-hand sides. Uniform batch (BatchMap bm): nzval holds the values of the bm.nbatch members, and the residuals of the active members go to the compact columns 1:(bm.nrhs bm.nact) of W.R. One launch, one work item per (row, right-hand side), gather-based (no atomics). Asynchronous.

residual_norms!(W, B; nrhs, transposed = false, bm = nothing) -> W.norms_host

‖W.R[:, k]‖₂² and ‖B[:, k]‖₂² for the compact columns k (1:nrhs, or the active members' columns of the BatchMap bm) reduced on the device (one workgroup per right-hand side, @localmem tree reduction) into W.norms, then copied to W.norms_host (one host synchronization).

resolve_subtree_budgets(budgets, max_local) -> Vector{Int}

The regime-A budgets an analysis uses on a backend with max_local bytes of static local memory per workgroup (max_local_bytes, issue #60), sorted. The default budgets (DEFAULT_SUBTREE_BUDGETS, 16, 32 and 48 KiB) are clamped to max_local (duplicates dropped), so the defaults work on every backend; an explicit budget above max_local raises InvalidValueError, since its kernel could not be launched there.

resolved_pivot_epsilon(opts::Options, T) -> Float64

opts.pivot_epsilon, or default_pivot_epsilon(T) when it was not set.

rhs_count(B, n; transposed = false) -> nrhs

Number of right-hand sides of the user array B for a system of size n: 1 for an n vector, size(B, 2) for an n × nrhs matrix (size(B, 1) for a transposed, i.e. row-major, nrhs × n one), length(B) ÷ n for a strided vector. Raises DimensionMismatch when the shape does not fit n.

scale_values!(out, weights, nzval) -> out

out[k + (m - 1) nnz] = weights[k] nzval[k + (m - 1) nnz] for the nnz = length(weights) stored entries of every batch member m (one launch, asynchronous): the values rᵢ aᵢⱼ cⱼ the numeric phase factors with matching.

scatter_A!(numeric, symbolic, nzval, first, count) -> numeric

Assemble the stored values nzval (on the backend of numeric) into the zeroed panels of the fronts symbolic.group_nodes[first:(first + count - 1)] through amap_ptr/amap_src/amap (one workgroup per front, one work item per nonzero; duplicated entries are summed by one work item in nzval order, negative offsets add the conjugate). Asynchronous.

schedule_depth(parent, post, counts, amalgamation = DEFAULT_AMALGAMATION) -> (sdepth::Int, nsupernodes::Int)

Depth of the supernodal schedule of an ordering and its number of supernodes: the fundamental_supernodes of the etree parent with postorder post and column counts counts (etree numbering, from etree, postorder, colcounts), relaxed by amalgamate with the amalgamation limits (nothing: no amalgamation, use_superpanels = 0), and the number of levels of their tree (tree_levels). It equals the nlevels of the Schedule of the analysis of that ordering (outside Schur complement mode).

schedule_elsize(structure, T) -> Int

Bytes per packed entry of the schedule's local-memory and update-stack budgets: sizeof(T), or 2 sizeof(T) for structure "G" (LU fronts and contribution blocks are two packed triangles, L and Uᵀ).

schur_flags(indices, n) -> BitVector

The Schur set of user_schur_indices (n flags 0/1, 1 = the row and column belong to the Schur complement) as a BitVector; raises InvalidValueError when indices is missing, has the wrong length or selects no index.

schur_matrix!(dest, numeric, symbolic; transposed = false) -> dest

Write the Schur complement assembled by assemble_schur! into dest: an ns × ns dense matrix on the backend of numeric (every entry), or a CSR whose rowptr and colval already hold the pattern to export (the solver writes them, see getparam(solver, "schur_matrix")): only its values are written. transposed: the solver's matrix is a transposed CSR (a CSC input), whose stored Schur complement is the transpose of the user's. One launch; asynchronous.

schur_pattern(P::SymmetricPattern, schur::AbstractVector{Bool}) -> (rowptr, colval)

Symbolic pattern of the Schur complement S = A₂₂ − A₂₁ A₁₁⁻¹ A₁₂ of the flagged rows and columns (A₂₂), in their increasing original order, as 1-based CSR arrays with sorted columns, both triangles and the whole diagonal: (i, j) is in the pattern when i == j, when (i, j) is in P, or when both are adjacent to one connected component of the graph of the other vertices (a fill path through eliminated vertices only). P is the symmetric pattern of the analysis (A + Aᵀ for "G"), so for LU the pattern may hold entries that are numerically zero by structure (explicit zeros).

schur_supernode_partition(P::SymmetricPattern, perm, ns, opts::Options = Options()) -> SupernodePartition

supernode_partition for Schur complement mode (PLAN §3.6): the last ns columns of perm (the Schur block, compute_schur_ordering) form one dense root supernode, the last one, with its columns in the order of perm; its etree is the chain of those columns and its column counts are the dense ones (ns, ns - 1, …, 1). The other columns are partitioned as usual (fundamental supernodes and amalgamate on their elimination forest, so nothing is merged into the Schur root); a supernode whose last column has its etree parent in the Schur block becomes a child of the Schur root.

select_impl(op::Symbol, X::AbstractArray, impl::Symbol) -> Symbol

Resolve impl for dense op op on the backend and element type of X: :auto becomes the first entry of dense_impls; an explicit :generic or :vendor the capability table does not list raises NotSupportedError; anything else raises InvalidValueError.

SparseDirectSolver.set_log_level!(level) -> Int

Set the verbosity of the solver's log messages: 0/"none" (default; the messages are @debug records), 1/"info" (phase summaries: analysis, factorization, refinement steps performed, as @info), 2/"debug" (also the relative residual of every refinement step). The initial value comes from the environment variable SDS_LOG_LEVEL. Returns the previous level.

set_members!(numeric, members) -> numeric

Make members (1-based, increasing, a subset of 1:numeric.nbatch) the active batch members of the next phases: copied to numeric.members (a host-to-device copy, only when they change) and split into runs of consecutive members. Asynchronous.

sncols(sp, s) -> UnitRange{Int}

Columns of supernode s (supernodal numbering).

snrows(sp, s)

Sorted rows of supernode s (a view), its own columns first.

snwidth(sp, s) -> Int

Number of columns w of supernode s; its panel is length(snrows(sp, s)) × w.

solve_bwd_kernel!(backend, WG)(Y, U, factor, piv, list, bm, count, subtree_ptr, subtree_nodes, super_ptr,
                               rowptr, rowval, front_ptr, front_nrows, front_ncols, child_ptr, Val(SUB),
                               Val(LDL), Val(CONJ), Val(WG); ndrange = WG * count * ncols)

Backward sweep Lᴴ x = z (Lᵀ when CONJ = false) of right-hand side r (of batch member _bm_cmember(bm, r)) by workgroup G = g + count (r - 1) over its fronts (the subtree list[bm.first + g - 1] in reverse processing order, SUB = true, or one front): gather the finished rows below, GEMV with L21ᴴ, TRSV with L11ᴴ, then (LDL) restore the order before the front's local pivoting.

solve_conjugated(structure, T, csc_input::Bool, solve_mode) -> Bool

Whether op(A) (solve_mode 0: A, 1: Aᵀ, 2: Aᴴ) of the user matrix is the conjugate of the stored, factorized matrix M (A = M, or A = Mᵀ for a CSC input). Real T: never. Complex symmetric "S" (M = Mᵀ): for Aᴴ. Hermitian "H"/"HPD" (Mᵀ = conj(M)): for Aᵀ of a CSR input, and for A and Aᴴ of a CSC input. General "G": for Aᴴ (op(A) = conj(Mᵀ) or conj(M), see solve_transposed).

solve_fwd_kernel!(backend, WG)(Y, U, factor, piv, list, bm, count, subtree_ptr, subtree_nodes, super_ptr,
                               rowptr, rowval, front_ptr, front_nrows, front_ncols, child_ptr, child_list,
                               relind_ptr, relind, Val(SUB), Val(DET), Val(LDL), Val(WG);
                               ndrange = WG * count * ncols)

Forward sweep L z = y of right-hand side r (a column of Y, of batch member _bm_cmember(bm, r)) by workgroup G = g + count (r - 1) over its fronts: the subtree list[bm.first + g - 1], its supernodes in processing order (SUB = true, regime A), or the single front list[bm.first + g - 1] (regime B). Per front: (deterministic, DET) zero its update buffer and pull the children's; (LDL) reorder its columns of Y by its local pivot order piv; TRSV with L11 on its columns of Y; GEMV with L21 into its update buffer or (atomic variant) into the rows of Y below.

solve_io(S::Symbolic, scaling, transpose::Bool) -> (inperm, inscale, outperm, outscale)

The permutations and scalings of the right-hand side (permute_rhs!) and of the solution (unpermute_solution!, add_correction!) of a solve with the factored matrix M (with Mᵀ when transpose): S.perm both ways and no scaling without matching (scaling === nothing); with a SolveScaling (rperm, rscale, cperm, cscale) for M and (cperm, cscale, rperm, rscale) for Mᵀ = Qᵀ Dc Aᵀ Dr.

solve_memory(symbolic, T, nrhs = 1) -> Int64

Bytes of the SolveWorkspace of nrhs right-hand sides.

solve_transposed(structure, csc_input::Bool, solve_mode) -> Bool

Whether the solve of op(A) (solve_mode 0: A, 1: Aᵀ, 2: Aᴴ) uses the transpose Mᵀ of the stored, factorized matrix M (A = M, or A = Mᵀ for a CSC input), conjugated or not (solve_conjugated). Only for structure "G" (LU): Mᵀ for Aᵀ/Aᴴ of a CSR input and for A of a CSC input. Symmetric and Hermitian structures never need it.

strided_batch(buf, offset, m, n, stride, count) -> AbstractArray{T,3}

The uniform batch of count column-major m × n matrices stored in the flat vector buf, the first starting after offset entries and each stride ≥ m n entries after the previous one, as an m × n × count view sharing buf (no copy). stride must be a multiple of m, and when stride > m n the buffer must extend over the full stride of the last member too (offset + stride * count ≤ length(buf)).

structural_zero_pivots(cands::PairCandidates, P::SymmetricPattern, perm, pairs = []) -> Vector{Int}

The candidates of cands, not in pairs, whose pivot is structurally zero when A[perm, perm] is factored column by column, in elimination order.

With generic values, the k-th pivot is nonzero exactly when the leading block A[perm[1:k], perm[1:k]] is structurally nonsingular (has a perfect matching on its pattern) given that the previous one is; the diagonal of a candidate counts as absent. The test runs incrementally: a perfect matching of the prefix is kept, a non-candidate is matched to its own diagonal, and a candidate looks for an augmenting path from its row to its column through earlier columns (|aᵢⱼ| > 0 edges, and the diagonal of a non-candidate row that an earlier path moved off it). A candidate without one is reported and matched to itself, as if a partner fixed it. A pair (a, c) with c right after a is matched crosswise (a 2×2 block). An elimination-tree leaf candidate is always reported; so is a candidate whose earlier neighbours are all used up by earlier candidates (two dual rows coupled to the same single primal).

subtree_capacity(budget, elsize, reserve = SUBTREE_LOCAL_RESERVE) -> Int

Entries of elsize bytes a regime-A subtree may keep in local memory under budget: subtree_local_bytes minus reserve (the kernel's other local memory, subtree_local_reserve), divided by elsize (0 when the budget is too small).

subtree_cholesky_kernel!(backend, WG)(factor, stack, info, nzval, amap, amap_ptr, amap_src, trees, bm,
                                      subtree_ptr, subtree_nodes, front_ptr, front_nrows, front_ncols, cb_ptr,
                                      local_front, local_cb, child_ptr, child_list, relind_ptr, relind,
                                      Val(NE), Val(WG); ndrange = WG * count * bm.nact)

Fused regime-A kernel: workgroup G takes subtree t = trees[bm.first + g - 1] of batch member k (g, k from the BatchMap bm) and processes its supernodes subtree_nodes[subtree_ptr[t]:(subtree_ptr[t+1]-1)] in order with the serial stack of packed fronts and contribution blocks in a @localmem buffer of NE entries (offsets local_front, local_cb): zero, scatter A, extend-add the children, factor the front's first w columns (Cholesky of F11, L21, Schur complement), write the panel and info[v], then move the contribution block down, or, for the subtree root (cb_ptr > 0), write it to the update stack (packed lower triangle).

subtree_ldlt_kernel!(backend, WG)(factor, stack, info, stats, d, piv, pivot_kind, psign, perm, aux, nzval, amap,
                                  amap_ptr, amap_src, trees, bm, subtree_ptr, subtree_nodes, super_ptr,
                                  front_ptr, front_nrows, front_ncols, cb_ptr, local_front, local_cb, child_ptr,
                                  child_list, relind_ptr, relind, prm, Val(NE), Val(WG);
                                  ndrange = WG * count * bm.nact)

LDLᵀ/LDLᴴ kernel of regime A: workgroup G takes subtree t = trees[bm.first + g - 1] of batch member k and processes its supernodes in order with the serial stack of packed fronts and contribution blocks in a @localmem buffer of NE entries (as subtree_cholesky_kernel!): zero, scatter A, extend-add the children, factor the front's first w columns with the pivoted LDLᵀ/LDLᴴ of front_ldlt_kernel!, update the trailing contribution block, write the unit-lower panel, D, piv, the pivot kinds and the statistics, then move the block down or (subtree root) write it to the update stack.

subtree_local_bytes(budget) -> Int

@localmem bytes of the regime-A kernel of a budget class: the largest of SUBTREE_LOCAL_SIZES (8, 16, 32, 48, 64 KiB) that is ≤ budget, 0 below 8 KiB. The size is a Val parameter of the kernel, so a short list of instances serves every budget.

subtree_local_layout(sp, schedule) -> (local_front, local_cb, local_len)

Local-memory layout (1-based entry offsets) of the regime-A subtrees: a serial stack in each subtree's processing order. The packed front of v (f(f+1)/2 entries) starts at local_front[v], right above the contribution blocks of its children (which are the top of the stack, contiguous); after the factorization its packed m(m+1)/2 contribution block (the trailing triangle of the front) moves down to local_cb[v], the start of the first child's block, or of the front when it has no children. Subtree roots keep local_cb = 0 (their block goes to the update stack). local_len[t] is the peak of subtree t in entries (schedule.subtree_peak[t] / elsize).

subtree_local_reserve(structure) -> Int

SUBTREE_LOCAL_RESERVE_LDLT for the LDLᵀ/LDLᴴ structures "S"/"H", else SUBTREE_LOCAL_RESERVE.

subtree_lu_kernel!(backend, WG)(factor, ufactor, stack, ustack, info, stats, d, piv, pivot_kind, aux, nzval, amap,
                                amap_ptr, amap_src, trees, bm, subtree_ptr, subtree_nodes, super_ptr, front_ptr,
                                front_nrows, front_ncols, cb_ptr, local_front, local_cb, child_ptr, child_list,
                                relind_ptr, relind, prm, Val(NE), Val(WG); ndrange = WG * count * bm.nact)

LU kernel of regime A: workgroup G takes subtree t = trees[bm.first + g - 1] of batch member k and processes its supernodes in order with the serial stacks of the packed L and Uᵀ structures in two @localmem buffers of NE entries each (the layout of subtree_local_layout, shared by both): zero, scatter A, extend-add the children, factor the front's first w columns with the pivoted L D U of front_lu_kernel!, update the trailing contribution blocks, write the unit-lower panels, D, piv, the pivot kinds and the statistics, then move both blocks down or (subtree root) write them to the update stacks.

supernode_partition(P::SymmetricPattern, perm, opts::Options = Options();
                    etree_hint = nothing) -> SupernodePartition

The whole supernode step for the ordering perm (e.g. compute_ordering(P, opts).perm): etree, postorder, colcounts, fundamental_supernodes, amalgamate with opts.amalgamation unless opts.use_superpanels == 0, and the supernodal symbolic factorization.

etree_hint is the (parent, post, counts) triple of A[perm, perm] when the caller already computed it on the same pattern — compute_ordering carries the chosen candidate's triple on the Ordering — so the elimination tree, postorder and column counts are not recomputed. The triple must match (P, perm); the length is checked, the content is trusted.

sweep_solve!(X, ws, symbolic, numeric, B; transposed = false, deterministic = false, impl = :auto,
             transpose = false) -> X

Solve A X = B with the factor P A Pᵀ = L Lᴴ (or the LDLᵀ/LDLᴴ factor with its local pivot orders) of factorize!: permute_rhs! into ws.Y, forward_sweep!, diagonal_sweep!, backward_sweep!, unpermute_solution! into X. B and X are device arrays in the layouts of rhs_count (vector, n × nrhs matrix, strided vector; row-major with transposed = true), with at most max_rhs(ws) right-hand sides; X === B is allowed. Uniform batch: B and X hold the right-hand sides of all numeric.nbatch members (member after member), only the active members are solved. LU ("G"): transpose = true solves Aᵀ X = B. Asynchronous; allocates nothing on the device.

symbolic_analysis(A::CSR, structure, view = 'F'; opts = Options(), T = eltype(A), device = nothing, max_local = typemax(Int))
    -> Symbolic{Int, Vector{Int}}

The whole host analysis of PLAN §2.3: SymmetricPattern, compute_ordering (with the 2×2 pivot pair candidates or fixed pairs of analysis_pairs for "S"/"H"), supernode_partition on factor_pattern, build_schedule for element type T, build_layout and the maps of Symbolic. A.rowptr/A.colval are copied to the host once; A.nzval too when pairs are looked for (pairs_enabled), and only then. max_local: the bytes of local memory per workgroup of the backend the analysis is for (max_local_bytes), which bound the regime-A budgets (resolve_subtree_budgets); the default leaves them unbounded (host analyses, the CPU backend).

syrk!(C, A, α, β; uplo = :L, impl = :auto) -> C

Symmetric rank-k update of the uplo triangle of C: C ← α A Aᵀ + β C (plain transpose, also for complex T; see herk! for A Aᴴ). :vendor and :ka leave the other triangle untouched; :generic (mul!) overwrites it with the same symmetric update.

takes_c_path(schedule, s) -> Bool

Whether supernode s runs on the regime-C path of the numeric phase (vendor or KA dense calls per front): regime C, or a regime-B bin wider than REGIME_B_MAX_WIDTH (no fused kernel).

tree_height(parent) -> Vector{Int}

Height of every node of a forest given by parent (0 = root) with parent[s] > s: leaves 1, otherwise one more than the highest child.

tree_levels(parent) -> (height::Vector{Int}, nlevels::Int)

Height of every node of the forest parent (leaves have height 1, a node is one above its highest child) and the number of levels, the maximum height (0 for an empty forest). Nodes of equal height can be eliminated concurrently, so nlevels is the length of the critical path in the level schedule.

trsm!(side, uplo, trans, diag, α, A, B; impl = :auto) -> B

Triangular solve with multiple right-hand sides, BLAS trsm: B ← α op(A)⁻¹ B for side = :L, B ← α B op(A)⁻¹ for side = :R; A is lower (uplo = :L) or upper triangular, with unit diagonal if diag = :U; trans ∈ (:N, :T, :C).

trsm_strided_batched!(side, uplo, trans, diag, α, A, B; impl = :auto) -> B

Uniform strided batch of trsm! on 3-D arrays A (n × n × count) and B. Implementations: :vendor, :ka.

unpermute_solution!(X, Y, perm; transposed = false, conjugate = false, bm = nothing, scale = nothing) -> X

X[perm[k], r] = Y[k, r] (conjugated when conjugate; times scale[perm[k]] when a scale vector is given): the inverse of permute_rhs!, from the first nrhs columns of Y into the user array X (same layouts; with a BatchMap bm, only the columns of the active members of X are written). One launch; asynchronous.

vendor_gemm!(transA::Char, transB::Char, α, A, B, β, C) -> C
vendor_syrk!(uplo::Char, α, A, β, C) -> C                  # uplo triangle of α A Aᵀ + β C
vendor_herk!(uplo::Char, α, A, β, C) -> C                  # uplo triangle of α A Aᴴ + β C (complex)
vendor_trsm!(side, uplo, trans, diag, α, A, B) -> B
vendor_potrf!(uplo::Char, A) -> info::Int
vendor_potrf_info!(uplo::Char, A, info, idx) -> info       # LAPACK info into the device Int32 vector info[idx]
vendor_getrf!(A, ipiv) -> info::Int
vendor_sytrf!(uplo::Char, A, ipiv) -> info::Int            # Bunch–Kaufman (capability probe only)
vendor_gemm_strided_batched!(transA, transB, α, A, B, β, C) -> C     # 3-D arrays
vendor_trsm_batched!(side, uplo, trans, diag, α, A, B) -> B          # 3-D arrays
vendor_potrf_batched!(uplo, A, info) -> info                         # 3-D A, Int32 info
vendor_getrf_batched!(A, ipiv, info) -> info                         # 3-D A, Int32 ipiv matrix
vendor_batch_pointers(A) -> ptrs                                     # device vector of the member pointers of 3-D A
vendor_trsm_batched_ptrs!(side, uplo, trans, diag, α, m, n, Aptrs, lda, Bptrs, ldb, count)
vendor_potrf_batched_ptrs!(uplo, n, Aptrs, lda, info, count) -> info # Int32 info, count entries

The _ptrs variants take device vectors of member pointers built once by vendor_batch_pointers (the uniform-batch regime C, T17: no allocation in the numeric phase); vendor_trsm_batched!/vendor_potrf_batched! build them per call.

Vendor dense routines behind the impl = :vendor path of the dense interface. The generic methods throw NotSupportedError; backends add methods for their array types (host BLAS/LAPACK for Arrays, cuBLAS/cuSOLVER for CuArrays).

zero_fronts!(numeric, symbolic, first, count) -> numeric

Zero the panels and the update-stack contribution blocks of the count fronts symbolic.group_nodes[first:(first + count - 1)] (one workgroup per front). Asynchronous.

zero_pivot_pairs!(pairs, cands::PairCandidates, P::SymmetricPattern, perm, u) -> (added, unmatched)

One round of the pivot_pairs = "default" selection (issue #66): pair the 2×2 candidates whose pivot is structurally zero under the ordering perm (structural_zero_pivots), appending (partner, candidate) to pairs (kept sorted by candidate).

Why only those: every other candidate receives a structurally nonzero update of its diagonal from earlier columns before it is eliminated, and in-front pivoting deals with it. A structurally zero pivot stays zero (or tiny) whatever the values, so it needs a partner in its fully-summed block. On MadNLP K2 systems this pairs about a fifth of the candidates and costs 1.3–1.4× nnz(L) (instead of ~2× for every candidate, pivot_pairs), and leaves a few dozen perturbed pivots instead of about a thousand.

Partners: a neighbour j that is in no pair yet with |aᵢⱼ| ≥ u · maxₖ≠ᵢ |aᵢₖ| (and aᵢⱼ ≠ 0), another candidate allowed (a [0 a; a 0] block is a good 2×2 pivot, the slack–dual case of K2 systems). The candidates are taken by increasing number of such free partners at the start of the round, then by position in perm; each takes the free partner with the largest |aᵢⱼ| (ties: earliest in perm), as a weak coupling makes a poor 2×2 pivot. Returns the number of pairs added and the number of structurally zero pivots left without a partner.

Pairing moves the partner, so the analysis re-orders on the compressed graph and calls it again until no pair is added (compute_ordering).