# Cg! per-iteration trace differs between single call and multiple calls

**URL:** <https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408>\
**Category:** Numerics\
**Created:** [July 22, 2026, 5:24am UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408 "2026-07-22T05:24:34Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![pitsianis](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pitsianis/32/26588_2.png) [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Post date:** [July 22, 2026, 5:24am UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/1 "2026-07-22T05:24:34Z")

</div>

Calling `IterativeSolvers.cg!` for N iterations in one call produces a different per-iteration residual sequence than performing the same total number of iterations via multiple `cg!` calls.  
A single call with 6 iterations produces a different result from two calls of thee iterations, or six calls of one iteration. Continuing with the same `statevars` or reseting them, does not seem to have an effect.

Is this normal behavior in the `cg!` implementation?

See [pcg\_minrepro.jl](https://discourse.julialang.org/uploads/short-url/xC0WfMJoQvmJzHytnRtysUsVKuQ.jl) (2.0 KB)

Output sample:  
Example (n=50 run):

```julia-repl
continued 1x6 = [34.64101615137756, 33.22649545167231, 31.811947441173714, 30.397368307141303, 28.982753492378862, 27.568097504180436]
continued 3+3 = [34.64101615137756, 33.22649545167231, 31.811947441173714, 22.449999274606412, 22.31701560306767, 15.91631301611577]
restarted 3+3 = [34.64101615137756, 33.22649545167231, 31.811947441173714, 22.449999274606412, 22.31701560306767, 15.91631301611577]
continued 6x1 = [34.64101615137756, 17.620484724263306, 13.3463058496797, 11.283385733387556, 10.331019501894122, 9.610386186347737]
restarted 6x1 = [34.64101615137756, 17.620484724263306, 13.3463058496797, 11.283385733387556, 10.331019501894122, 9.610386186347737]

```

## Environment

- Julia: 1.12.6
- OS: macOS (arm64-apple-darwin24.0.0) on Apple M2 Max
- IterativeSolvers package: v0.9.4

---

<div class="post-metadata">

**Author:** ![pitsianis](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pitsianis/32/26588_2.png) [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Post date:** [July 22, 2026, 4:35pm UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/2 "2026-07-22T16:35:13Z")

</div>

It looks like this issue has been raised five years ago 🙂

> <https://github.com/JuliaLinearAlgebra/IterativeSolvers.jl/issues/291>
>
> There are cases where iterative solvers finish prematurely by hitting \`maxiter\`,… before the convergence criterion \`‖rₙ‖ = ‖b - A xₙ‖ ≤ tol\` is met. In that case, it is useful to have the capability to resume the solvers from the point they finished.
> 
> Resuming iterative solvers requires not only the last approximate solution \`xₙ\`, but also other state variables. For example in CG, these state variables are the last residual vector \`rₙ\` and the last search direction vector \`pₙ\`. Luckily, these state variables are already stored in \`CGStateVariables\`, but currently the code is not written to resume the CG algorithm from the state stored in the passed \`statevars::CGStateVariables\`.
> 
> I have implemented the resumability in \`cg.jl\` and it works as expected. Before posting a PR for this new capability, I would like to have feedbacks from other users as to whether the resumability would be a useful feature that the community would welcome. Also, I wonder if it would make sense to implement the resumability in as many iterative solvers as possible.

and the solution outlined:

> Resuming iterative solvers requires not only the last approximate solution xₙ, but also other state variables. For example in CG, these state variables are he last residual vector rₙ and the last search direction vector pₙ. Luckily, these state variables are already stored in `CGStateVariables`, but currently the code is not written to resume the CG algorithm from the state stored in the passed `statevars::CGStateVariables`.

So is there a good reason not to have it?

---

<div class="post-metadata">

**Author:** ![j-fu](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j-fu/32/11373_2.png) [@j-fu](https://discourse.julialang.org/u/j-fu)\
**Post date:** [July 22, 2026, 7:19pm UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/3 "2026-07-22T19:19:56Z")

</div>

CG is a Krylov subspace method which tries to update the solution in orthogonal directions, which means that behind the scenes, a orthogonalization process runs which delivers orthogonal residuals and A-ortohogonal search directions. A restart without the the context accumulated so far loses the orthogonal basis built so far.

CG for symmetric positive definite A has the advantage that the new orthogonal direction only depends on the residual and the previous direction which occurs as “state vector”, and does not need to carry the residual of all the previous iterations (unlike GMRES). If the state vector is lost, the orthogonality information is discarded.

With restart(s), you thus you should observe a slowdown of the convergence.

---

<div class="post-metadata">

**Author:** ![pitsianis](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pitsianis/32/26588_2.png) [@pitsianis](https://discourse.julialang.org/u/pitsianis)\
**Post date:** [July 23, 2026, 7:00pm UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/4 "2026-07-23T19:00:46Z")

</div>

You are absolutely right. In CG, the new orthogonal search direction only depends on the residual rₙ and the previous search direction pₙ. In the `cg!` implementation this information is maintained in the `statevars::CGStateVariables`.

The genius in our group figured out how to use the iterator constructed and used by `cg!` to achieve resumption of iterations from the latest xₙ, rₙ and pₙ. But the question still remains. Why not expose this functionality directly?

---

<div class="post-metadata">

**Author:** ![ziolai](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ziolai/32/23422_2.png) [@ziolai](https://discourse.julialang.org/u/ziolai)\
**Post date:** [July 24, 2026, 11:13am UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/5 "2026-07-24T11:13:38Z")

</div>

If the CG algorithm uses short recurrences to compute iterands as @j-fu and the theory suggests, it should suffice to provide the lastest iterand x\_n at restart. Thus only limited advantage to keep track of r\_n and p\_n.

The situation does change in case of long term recurrence methods like GMRES and friends or in case of solving linear systems with same matrix and various right-hand hands.

Does this address your concerns?

---

<div class="post-metadata">

**Author:** ![j-fu](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/j-fu/32/11373_2.png) [@j-fu](https://discourse.julialang.org/u/j-fu)\
**Post date:** [July 24, 2026, 12:04pm UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/6 "2026-07-24T12:04:20Z")

</div>

BTW there are two more iterative solver packages: Krylov.jl and KrylovKit.jl. I is possible that one of the provides an API which allows to keep the context with restarts.

---

<div class="post-metadata">

**Author:** ![amontoison](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/amontoison/32/218741_2.png) [@amontoison](https://discourse.julialang.org/u/amontoison)\
**Post date:** [August 4, 2026, 7:29am UTC](https://discourse.julialang.org/t/cg-per-iteration-trace-differs-between-single-call-and-multiple-calls/138408/7 "2026-08-04T07:29:13Z")

</div>

For some solvers like CG, it should not be too hard to have an option to resume the solver from the current workspace but it may be way harder for some solvers where we update the solution with a lag like MINRES-QLP or only at the end like GMRES.

I like to compare what you request as a reverse-communication API, where you have the full control on the solve and can do one by one the iteration and check the status / information from the solver.

To be able to do what you request in `Krylov.jl` or `IterativeSolvers.jl`, we need to store a few additional scalar variables in the workspace, they are the key missing part to resume the solver from a current state.  
For example, `\alpha` and `\beta` in CG.  
In practice we only store vectors in the `workspace` and keep the local scalar variable on the stack.

Please open an issue / discussion an `Krylov.jl` if you are interested by this feature.