# Any conjugate gradient method allowing a custom inner product?

**URL:** <https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982>\
**Category:** Numerics\
**Tags:** question\
**Created:** [July 29, 2022, 3:06pm UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982 "2022-07-29T15:06:39Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![ranocha](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ranocha/32/35588_2.png) [@ranocha](https://discourse.julialang.org/u/ranocha)\
**Post date:** [July 29, 2022, 3:06pm UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982/1 "2022-07-29T15:06:39Z")

</div>

Is there any conjugate gradient (CG) implementation in Julia allowing one to specify the inner product (with a fixed symmetric and positive definite matrix `M`, say)?

In my use case, I know my system matrix `A` is symmetric and positive definite with respect to `M`, i.e., `x' * M * A * x > 0` for all `x != 0` and `y' * M * (A * x) == (A * y)' * M * x` for all `x, y`. Ideally, I would like to avoid scaling the system matrix and vectors myself for using a CG method. Just curious to know whether there’s such a flexible package out there.

---

<div class="post-metadata">

**Author:** ![mstewart](https://avatars.discourse-cdn.com/v4/letter/m/b5a626/32.png) [@mstewart](https://discourse.julialang.org/u/mstewart)\
**Post date:** [July 29, 2022, 6:38pm UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982/2 "2022-07-29T18:38:27Z")

</div>

I don’t know the answer to the question about specific implementations of CG, but since your property of being symmetric and positive definite with respect to M is equivalent to MA being symmetric positive definite, it seems reasonable to just solve MA x = Mb using the [CG routine](https://iterativesolvers.julialinearalgebra.org/v0.8/linear_systems/cg/) from `IterativeSolvers.jl` and form MA as a product of `LinearMaps`:

```
using IterativeSolvers
using LinearMaps
C = LinearMap(M) * LinearMap(A)
x=cg(C, M*b)

```

That’s pretty easy. In fact, since CG includes more inner products than matrix-vector multiplies, I think this actually involves fewer multiplies by M than if you tried to implement CG for a different inner product. If multiplying by M is nontrivial, this is probably faster than what you were hoping to do.

---

<div class="post-metadata">

**Author:** ![SteffenPL](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/steffenpl/32/206270_2.png) [@SteffenPL](https://discourse.julialang.org/u/SteffenPL)\
**Post date:** [July 29, 2022, 6:41pm UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982/3 "2022-07-29T18:41:39Z")

</div>

[Conjugate Gradients · IterativeSolvers.jl](https://iterativesolvers.julialinearalgebra.org/dev/linear_systems/cg/) has the optional argument `Pl` to add a left precondition. (That should correspond to a change to the inner product as @mstewart explained)  
As the documentation say it uses a different algorithm which is specialised for this setting.

---

<div class="post-metadata">

**Author:** ![mstewart](https://avatars.discourse-cdn.com/v4/letter/m/b5a626/32.png) [@mstewart](https://discourse.julialang.org/u/mstewart)\
**Post date:** [July 29, 2022, 7:11pm UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982/4 "2022-07-29T19:11:06Z")

</div>

It’s a tempting idea, and there may be some trick to make use of a preconditioner, but doing it in a direct way for this problem doesn’t seem like an option, given that PCG assumes that both the preconditioner and A are symmetric positive definite. That’s not the case here.

---

<div class="post-metadata">

**Author:** ![antoine-levitt](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/antoine-levitt/32/4008_2.png) [@antoine-levitt](https://discourse.julialang.org/u/antoine-levitt)\
**Post date:** [July 29, 2022, 7:42pm UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982/5 "2022-07-29T19:42:21Z")

</div>

KrylovKit supports any vector type that have custom dot and norm methods, without assuming they are the standard ones. See eg [Introduction · KrylovKit.jl](https://jutho.github.io/KrylovKit.jl/stable/man/intro/#KrylovKit.InnerProductVec)

---

<div class="post-metadata">

**Author:** ![ranocha](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/ranocha/32/35588_2.png) [@ranocha](https://discourse.julialang.org/u/ranocha)\
**Post date:** [July 30, 2022, 5:57am UTC](https://discourse.julialang.org/t/any-conjugate-gradient-method-allowing-a-custom-inner-product/84982/6 "2022-07-30T05:57:01Z")

</div>

Thank you all! I can only select one answer and have a hard time selecting one 😅
