Finite-budget stable block reordering for Schur and QZ forms

Inputs are an order $d\ge 0$ real or complex Schur matrix $T$, or a regular real or complex generalized Schur pair $(S,T)$, with the declared diagonal block partition. Complex Schur and QZ forms have scalar blocks; real Schur forms have scalar blocks or real $2\times 2$ blocks; real generalized Schur forms have real $1\times 1$ or $2\times 2$ blocks. Each generalized diagonal block pair must be regular. Optional accumulated factors are exact unitary/orthogonal products; if omitted they are identity. Other inputs are a CARE selector (strict left half-plane) or DARE selector (strict unit disk), a required integer selected dimension $0\le r\le d$, nonnegative integer budgets $K$ and $W$, a positive local separation threshold $\eta$, and either an exact arithmetic controller or a bounded outward evaluator.

Use the CARE selector with an ordinary Schur form and the DARE selector with a generalized Schur form; an unmatched form-selector combination returns invalid-form.

Before any charged operation, initialize attempts, accepted swaps and primitive/evaluator counts to zero, diagnostic lists to empty, and every uncomputed diagnostic to unavailable. In particular, the actual selected dimension, block labels, boundary and separation predicates, and initial inversion count remain unavailable until their respective scans or tests finish; the requested dimension $r$ is always retained separately. Report only completed tests and certificates, together with the last committed arrays and factors, on an early exit. Unavailable predicates assert neither truth nor falsity.

The entries of the input arrays are treated as their exact represented real or complex values. All local transformed arrays and accumulated factors are retained exactly as finite products; a rounded dense copy is only a separately evaluated diagnostic. Before each scalar primitive or outward evaluation, return work-budget-exhausted if no units of $W$ remain; otherwise consume one unit before requesting the operation. Charge validation, label and inversion-count scans, local construction, and the final ordering scan to this same budget. For rational or dyadic arrays, the exact core uses rational arithmetic and positive square roots, with exact real-algebraic sign comparisons where needed. Exact decisions on arbitrary real inputs are an exact-real arithmetic model, not a bit-computability claim. The bounded-evaluator branch returns arithmetic-inconclusive whenever a needed zero, sign, division, or normalization remains unresolved.

1. Validate dimensions, the declared block pattern, and regularity of each generalized diagonal pair in one finite scan. For ordinary scalar blocks over either field classify $t_{ii}$ by $\mathrm{Re}{t}_{ii}<0$ or $>0$ for CARE; for real $2\times 2$ blocks first verify negative discriminant and then classify by trace less than or greater than zero. For generalized scalar blocks over either field $[s_{ii}:t_{ii}]$, verify the pair is nonzero and classify by $\mathrm{\mid }{s}_{ii}\mathrm{\mid }<\mathrm{\mid }{t}_{ii}\mathrm{\mid }$ or $\mathrm{\mid }{s}_{ii}\mathrm{\mid }>\mathrm{\mid }{t}_{ii}\mathrm{\mid }$ for DARE; thus $t_{ii}=0$, $s_{ii}\mathrm{\ne }0$ is complementary and represents infinity without division. For a real generalized $2\times 2$ block, verify $\det (T_b)\mathrm{\ne }0$ and that $\det (S_b-\lambda T_b)$ has negative discriminant; its conjugate roots have common squared modulus $\det (S_b)\mathrm{/}\det (T_b)$, so classify by strict comparison of this ratio with one, cross-multiplying with the certified sign of $\det (T_b)$. Exact equality with the selected boundary returns boundary-spectrum. A bounded enclosure touching the boundary returns separation-inconclusive. Require that the selected block sizes sum to $r$, or return dimension-mismatch. Invalid patterns or singular local homogeneous pairs return invalid-form. Compute the initial inversion count once by scanning all ordered block pairs after classification.

2. Give every whole block its selected or complementary label. Locate the leftmost adjacent complementary-selected inversion. If none exists, return ordered with the current exact arrays and accumulated factors. If an inversion exists but $K$ swap attempts have already been made, return swap-budget-exhausted. Otherwise consume one attempt before starting its local construction. Never split a real conjugate pair.

3. Let $p,q\in \{1,2\}$ be the two adjacent block sizes (over $\mathbb{C}$, $p=q=1$), with local blocks $T_{JJ}=\left[\begin{array}{cc}A& C\\ 0& B\end{array}\right]$ for ordinary Schur form, or $S_{JJ}=\left[\begin{array}{cc}A& C\\ 0& B\end{array}\right],T_{JJ}=\left[\begin{array}{cc}D& F\\ 0& E\end{array}\right]$ for QZ form. With column-vectorization, form the local coefficient system $H=I_q\otimes A-B^T\otimes I_p,h=-\mathrm{vec}(C)$ for ordinary Schur, with unknown $\mathrm{vec}(R)$, and $H=\left[\begin{array}{cc}{I}_q\otimes A& -B^T\otimes I_p\\ I_q\otimes D& -E^T\otimes I_p\end{array}\right],h=-\left[\begin{array}{c}\mathrm{vec}(C)\\ \mathrm{vec}(F)\end{array}\right]$ for QZ, with unknown $[\mathrm{vec}(R);\mathrm{vec}(L)]$. The transpose in $B^T$ and $E^T$ is intentional, also over $\mathbb{C}$. Solve $Hx=h$ and $HY=I$ by finite Gaussian elimination with complete pivoting: at each stage scan the finite remaining submatrix, choose an entry of largest certified modulus with lexicographic ties, swap rows and columns, eliminate, and back-substitute while undoing column permutations. An exact all-zero remaining submatrix returns swap-singular. An unresolved pivot or range failure returns separation-inconclusive or range-failure. Verify $Hx=h$ and $HY=I$ exactly, or verify a usable approximate-inverse certificate $\mathrm{\parallel }I-YH{\mathrm{\parallel }}_2\le \theta <1$ and still solve $Hx=h$ exactly. In the latter case it gives $\mathrm{\parallel }{H}^{-1}{\mathrm{\parallel }}_2\le \mathrm{\parallel }Y{\mathrm{\parallel }}_2\mathrm{/}(1-\theta )$. In the exact-inverse branch take $\beta =\mathrm{\parallel }Y{\mathrm{\parallel }}_F$ in exact arithmetic, or an outward upper bound for that norm; in the approximate-inverse branch take an outward bound for $\mathrm{\parallel }Y{\mathrm{\parallel }}_2\mathrm{/}(1-\theta )$. Require $\eta \beta \le 1$; otherwise return swap-separation-inconclusive. Record the threshold, solves, pivot decisions and certificate.

4. For ordinary Schur form set $V=[R;I_q]$. For QZ set $V=[R;I_q]$ and $W_L=[L;I_q]$. Complete the columns of $V$, and separately $W_L$, to exact orthonormal bases by Gram–Schmidt: process the graph columns first, then $e_1,\dots ,e_{p+q}$, skipping exactly zero residuals and normalizing each nonzero residual by the positive square root of its squared norm. The graph columns are independent because their last $q$ rows form $I_q$; the coordinate scan completes a basis. An unresolved zero or normalization decision in bounded mode returns arithmetic-inconclusive.

5. Let $Z_{\mathrm{l}\mathrm{o}\mathrm{c}}$ be the ordinary orthogonal/unitary completion; for QZ let $Z_{\mathrm{l}\mathrm{o}\mathrm{c}}$ be the right completion and $Q_{\mathrm{l}\mathrm{o}\mathrm{c}}$ the left completion. Propose the local transform $Z_{\mathrm{l}\mathrm{o}\mathrm{c}}^{\ast }T{Z}_{\mathrm{l}\mathrm{o}\mathrm{c}}$ in the ordinary case, or $Q_{\mathrm{l}\mathrm{o}\mathrm{c}}^{\ast }S{Z}_{\mathrm{l}\mathrm{o}\mathrm{c}}$ and $Q_{\mathrm{l}\mathrm{o}\mathrm{c}}^{\ast }T{Z}_{\mathrm{l}\mathrm{o}\mathrm{c}}$ in the QZ case. Apply it to the entire affected row and column panels, including upper off-diagonal blocks. The lower-left block is set to exact zero only by the local equations solved in step 3. For a real QZ form, factor each newly formed real $2\times 2$ diagonal block as $T_b=L_b{R}_b$ by Gram–Schmidt on its two columns, apply the embedded $L_b^{\ast }$ to the corresponding rows of both arrays, and update $Q_{\mathrm{l}\mathrm{o}\mathrm{c}}\leftarrow Q_{\mathrm{l}\mathrm{o}\mathrm{c}}{L}_b$ with the same embedded factor; such a block is nonsingular because its roots are finite and nonreal. A size-one infinite pair requires no division by its zero second coordinate. Commit the provisional panels and factor updates only after every required local operation succeeds, incrementing accepted swaps by one as part of that same commit; otherwise return the previous committed arrays with a failed-attempt diagnostic and leave accepted swaps unchanged.

6. After a committed exchange, swap the two labels and sizes and return to step 2. Ordinary accumulated factors are updated by right multiplication with $Z_{\mathrm{l}\mathrm{o}\mathrm{c}}$; generalized accumulated factors are updated as $Q\leftarrow Q{Q}_{\mathrm{l}\mathrm{o}\mathrm{c}}$ and $Z\leftarrow Z{Z}_{\mathrm{l}\mathrm{o}\mathrm{c}}$. The first $r$ columns of the updated factor are reported as the selected basis candidate; in QZ the right factor $Z$, not the left factor $Q$, supplies the right deflating basis.

On every exit return the current partial form, exact product factors, status, selected dimension, boundary and separation predicates, initial inversion count, attempts, accepted swaps, primitive/evaluator counts, original residual certificates, local inverse-norm certificates, and any residual transport data. Other statuses are work-budget-exhausted, arithmetic-inconclusive, range-failure, invalid-form, dimension-mismatch, boundary-spectrum, separation-inconclusive, swap-singular, swap-separation-inconclusive, and swap-budget-exhausted. ordered is distinct from any Riccati acceptance status. Reordering a computed form does not assert a structured Hamiltonian or symplectic backward-error theorem.