Generalized Residual Cutting (GRC)

• GRC is a Krylov subspace method that generalizes residual cutting by replacing relaxation steps with matrix-vector recurrences to solve large, sparse linear systems.
• It constructs correction vectors through residual minimization over a short recurrence of current and previous directions, ensuring stable convergence even for nonsymmetric matrices.
• When integrated as a stabilization layer for BiCGSTAB, GRC effectively filters unstable inner iterations to overcome stagnation and breakdown in challenging sparse problems.

Generalized Residual Cutting (GRC) is an iterative method for solving large sparse linear systems that was introduced as a generalization of the earlier residual cutting (RC) method. In the linear-solver literature, its defining feature is the construction of each correction by residual minimization over a short recurrence built from the current residual information and a limited set of previous directions, with the result that GRC is both a Krylov subspace method and a low-storage alternative to methods such as GMRES. The method was developed to extend RC beyond elliptic-problem settings governed by convergent relaxation schemes, and later work used GRC as a stabilization layer for BiCGSTAB on large, sparse, nonsymmetric systems (Abe et al., 2017).

1. Origins, motivation, and problem class

GRC emerged from the residual cutting method, which had originally been proposed for linear equations obtained from elliptic partial differential equations. In RC, a relaxation method such as SOR or Gauss–Seidel provides an approximate inner solve, and the outer iteration forms a correction by minimizing the new residual norm over a small span of current and previous directions. This construction is effective when the inner relaxation converges, but it is correspondingly restricted: in its original form, RC is not generally applicable to arbitrary sparse nonsymmetric systems (Abe et al., 2017).

The central motivation for GRC was therefore to retain the residual-minimizing structure of RC while removing its dependence on relaxation-friendly matrices. The generalization replaces the relaxation-based inner step by a matrix-vector recurrence, making the method applicable to general sparse matrix problems. In later work, the same residual-cutting idea was further framed as a stabilization mechanism that can be attached to broader Krylov solvers. In that setting, GRC is used not primarily to accelerate a method that already behaves well, but to improve robustness when convergence is irregular or when stagnation and breakdown occur, as in problematic BiCGSTAB runs on large, sparse, nonsymmetric systems (Abe, 19 Aug 2025).

The target problem is the sparse linear system

${\bf H}{\bf U}={\bf b},$

or equivalently

${\bf H}u=b,$

with residual

${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$

The method is discussed in connection with discretized PDEs, sparse matrix test problems from the University of Florida sparse matrix collection, and nonsymmetric and sometimes non-positive-definite matrices (Abe, 19 Aug 2025).

2. Core algorithmic structure

The basic residual-cutting update has the form

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$

The correction vector is not taken directly from the current residual. Instead, it is assembled from a temporary inner vector and a short history of previous corrections:

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$

The coefficients are chosen to minimize

$\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$

This residual minimization over a low-dimensional span is the defining operational principle inherited from RC (Abe et al., 2017).

The generalization from RC to GRC is expressed through the choice of the inner vector. In the formulation emphasized in the convergence paper, GRC uses

${\bf \Psi}^m=({\bf 1-H})\phi^{m-1}+{\bf r}^m,$

and more generally

$${\bf \Psi}_{\rm GRC}^m=G(\phi^{m-1},{\bf H},{\bf r}^m),$$

with $G$ linear. The chosen form replaces the inner relaxation solve of RC by a sparse matrix-vector recurrence, which is what makes the method applicable to general sparse linear systems rather than only to systems suited to relaxation methods (Abe et al., 2017).

A later description, used for stabilization, presents a modified Gram–Schmidt variant, GRC-MGS. At outer iteration $m$, one first computes an approximate solution ${\bf H}u=b,$0 to

${\bf H}u=b,$1

using an inner iterative process, and then forms a stabilized direction

${\bf H}u=b,$2

together with

${\bf H}u=b,$3

The outer step length is

${\bf H}u=b,$4

and the update becomes

${\bf H}u=b,$5

Here the stabilized direction is explicitly a linear combination of the current inner solve and several previous directions (Abe, 19 Aug 2025).

3. Krylov characterization and relation to conjugate residual methods

A central mathematical result is that GRC is a Krylov subspace method. The proof proceeds by showing that all correction vectors and residuals can be written as matrix polynomials in ${\bf H}u=b,$6 applied to the initial residual ${\bf H}u=b,$7. The first correction satisfies

${\bf H}u=b,$8

and the next residual is

${\bf H}u=b,$9

Inductively,

${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$0

so that

${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$1

This establishes that GRC is not merely a residual-smoothing heuristic but a genuine Krylov subspace method (Abe et al., 2017).

The method also has a distinct orthogonality structure. The new correction is constructed so that

${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$2

or equivalently

${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$3

In the paper’s description, this is a Gram–Schmidt-like orthogonalization in the ${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$4 inner product. The practical role of that construction is residual minimization within the available truncated subspace (Abe et al., 2017).

The closest formal relation is to the conjugate residual (CR) method. The CR search directions ${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$5 also lie in the same Krylov spaces and satisfy the same type of ${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$6-orthogonality. For symmetric ${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$7, the paper states that the GRC and CR directions coincide up to a scalar factor:

${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$8

Accordingly, for symmetric matrices, GRC and CR generate the same Krylov subspaces, produce the same direction vectors up to scaling, and have matching residual-minimizing behavior. The paper further notes that ${\bf r}^m={\bf b}-{\bf H}{\bf U}^m.$9 appears sufficient for symmetric $${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$0 in the form using

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$1

whereas if one uses the simpler choice

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$2

the procedure becomes completely identical to CR and $${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$3 is sufficient (Abe et al., 2017).

4. GRC as a stabilization framework for BiCGSTAB

A later development treats GRC as an outer stabilization framework wrapped around BiCGSTAB. In this construction, BiCGSTAB serves as the inner solver, but its iterates are not trusted directly as the final outer step. Instead, GRC combines the current BiCGSTAB-produced approximate solution with several recent directions, thereby filtering unstable inner behavior before the outer update is applied (Abe, 19 Aug 2025).

The resulting method, GRC-BiCGSTAB, follows this structure. Starting from $${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$4 and

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$5

the algorithm approximately solves

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$6

with a partial BiCGSTAB run. The inner BiCGSTAB process is stopped early when the residual norm is reduced to half of its initial value,

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$7

The current inner result is then combined with prior directions to form $${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$8, after which the outer solution and residual are updated by

$${\bf U}^{m+1}={\bf U}^m+\phi^m,\qquad {\bf r}^{m+1}={\bf r}^m-{\bf H}\phi^m.$$9

The stabilization paper describes this as convergence that is “further stabilized by GRC” (Abe, 19 Aug 2025).

The motivation is directly tied to BiCGSTAB’s failure modes. BiCGSTAB is smoother than the original BiCG method because it reduces irregular convergence behavior by stabilizing residuals, but it can still fail due to stagnation or breakdown. GRC acts as a stabilization layer because it does not rely on a single raw inner BiCGSTAB step; instead it constructs a new outer step from a linear combination of current and recent BiCGSTAB-generated directions. A plausible implication is that the method reduces sensitivity to unlucky recurrences by replacing a single unstable step with a controlled residual-cutting projection, but the paper itself formulates the point in terms of robustness rather than a new convergence theorem (Abe, 19 Aug 2025).

For reference, the BiCGSTAB inner recurrence used in the paper is the standard one:

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$0

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$1

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$2

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$3

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$4

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$5

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$6

In GRC-BiCGSTAB, these iterates serve as inner-loop approximations rather than as the sole driver of the final update (Abe, 19 Aug 2025).

5. Numerical behavior, robustness, and storage profile

The numerical picture presented in the GRC literature is consistently one of robustness-oriented trade-offs rather than universal speed dominance. In the convergence study, three benchmark situations are especially characteristic. For a Poisson-type problem on a nonuniform grid, RC performs very well in time because the inner SOR relaxation converges effectively; GRC converges more slowly at first and then accelerates; BiCGSTAB shows fluctuations in residual norm; and GMRES residual stops decreasing before convergence in that case. For a PDE matrix with convection term

$$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$7

RC fails because the inner SOR diverges, whereas GRC remains usable. For the matrix raefsky2, BiCGSTAB converges fastest, GMRES is slower, and GRC is slower than GMRES, but the residual does not blow up (Abe et al., 2017).

The large-scale experiments strengthen the memory-and-robustness argument. On a PDE matrix family whose size increases from about $$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$8 to about $$\phi^m=\alpha_1^m\Psi^m+\sum_{k=2}^{L}\alpha_k^m\phi^{m-k+1}.$$9, BiCGSTAB can diverge as size increases, GMRES may terminate with out-of-memory for the largest cases, and GRC remains effective; in the largest case reported, it is the only method said to converge. For these large matrices, GRC typically requires less than half the steps and time of GMRES in the reported tests (Abe et al., 2017).

The BiCGSTAB-stabilization study reaches a similar conclusion from a different angle. Its experiments were implemented in GNU Octave in double precision, with convergence criterion

$\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$0

GRC parameter

$\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$1

and inner BiCGSTAB stopping threshold

$\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$2

On some problems, plain BiCGSTAB is faster, including bfwa62 and pores_1. But when BiCGSTAB becomes unstable or breaks down, GRC-BiCGSTAB succeeds where BiCGSTAB fails: in “Problem 1” BiCGSTAB breaks down at iteration 104; in olm100 at iteration 410; in toeplitz100_2.0 at iteration 45; and in toeplitz100_2.3 at iteration 72. The paper reports that the GRC outer layer smooths the residual history and suppresses the sharp fluctuations seen in raw BiCGSTAB (Abe, 19 Aug 2025).

The storage and per-step operation counts given in the convergence paper summarize GRC’s intended computational niche:

Method	memory (vectors)	MATVEC per step
GRC	$\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$3	1
BiCGSTAB	5	2
GMRES	$\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$4	1

The same comparison also reports DOT per step as $\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$5 for GRC, $\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$6 for BiCGSTAB, and $\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$7 for GMRES, with the experiments using $\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$8 for GRC and restart $\|{\bf r}^{m+1}\|=\|{\bf r}^m-{\bf H}\phi^m\|.$9 for GMRES (Abe et al., 2017). This places GRC as a low-storage residual-minimizing method with a much smaller retained basis than GMRES.

6. Scope, misconceptions, and acronym ambiguity

Within numerical linear algebra, “Generalized Residual Cutting” refers to the residual-cutting Krylov method described above. A common source of confusion is that the acronym “GRC” has been used in unrelated areas of machine learning. In Transformer architecture work, GRC denotes “Gated Residual Connections,” a feature-wise gating mechanism for residual additions,

${\bf \Psi}^m=({\bf 1-H})\phi^{m-1}+{\bf r}^m,$0

introduced as part of a broader transformer-enhancement framework alongside the Evaluator Adjuster Unit; despite the shared acronym, this is not the linear-solver method (Dhayalkar, 2024).

A second unrelated use appears in decoder-only language modeling, where “GRC” stands for “Unifying Reasoning-Driven Generation, Retrieval and Compression.” There the term names a framework built around trainable meta latent tokens, a combined objective

${\bf \Psi}^m=({\bf 1-H})\phi^{m-1}+{\bf r}^m,$1

and one-forward-pass unification of generation, retrieval, and compression. This acronym is likewise independent of Generalized Residual Cutting in sparse linear algebra (Miao et al., 9 May 2026).

A further neighboring but distinct literature concerns generalized closure operators in cutting-plane theory, defined by

${\bf \Psi}^m=({\bf 1-H})\phi^{m-1}+{\bf r}^m,$2

That work studies polyhedrality, rationality, finite generation, and the limit of iterated closures, and it subsumes the Chvátal closure and split closure as special cases. Although the word “cutting” is shared, these geometric closure operators are not the Krylov residual-cutting algorithm for sparse linear systems (Averkov, 2011).

In the linear-solver sense, the settled characterization from the cited literature is narrower and more precise. GRC is a residual-minimizing, truncated Krylov method derived from RC, generalized to sparse matrix problems by a matrix-vector recurrence, closely related to CR, equivalent to CR for symmetric matrices, more memory-efficient than GMRES, and often more robust than BiCGSTAB on difficult nonsymmetric problems, especially when used as an outer stabilization layer (Abe et al., 2017).