Skip to content

Commit

Permalink
Provide a more thorough explanation of precs, including reinit! with
Browse files Browse the repository at this point in the history
reuse_precs=true.
  • Loading branch information
j-fu committed Oct 22, 2024
1 parent 89ed45b commit 8b8b171
Show file tree
Hide file tree
Showing 2 changed files with 22 additions and 18 deletions.
1 change: 1 addition & 0 deletions docs/Project.toml
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
[deps]
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
IncompleteLU = "40713840-3770-5561-ab4c-a76e7d0d7895"
LinearSolve = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae"

[compat]
Expand Down
39 changes: 21 additions & 18 deletions docs/src/basics/Preconditioners.md
Original file line number Diff line number Diff line change
Expand Up @@ -36,34 +36,30 @@ and `Pr` for right preconditioner, respectively. By default, if no preconditione
the identity ``I``.


In the following, we will use the `DiagonalPreconditioner` to define a two-sided
preconditioned system which first divides by some random numbers and then
multiplies by the same values. This is commonly used in the case where if, instead
of random, `s` is an approximation to the eigenvalues of a system.
In the following, we will use a left sided diagonal (Jacobi) preconditioner.

```@example precon
```@example precon1
using LinearSolve, LinearAlgebra
n = 4
s = rand(n)
Pl = Diagonal(s)
A = rand(n, n)
b = rand(n)
Pl=Diagonal(A)
prob = LinearProblem(A, b)
sol = solve(prob, KrylovJL_GMRES(), Pl = Pl)
sol.u
```

Alternatively, preconditioners can be specified via the `precs` argument to the constructor of
an iterative solver specification. This argument shall deliver a function mapping `A` and a
an iterative solver specification. This argument shall deliver a factory method mapping `A` and a
parameter `p` to a tuple `(Pl,Pr)` consisting a left and a right preconditioner.


```@example precon2
using LinearSolve, LinearAlgebra
n = 4
s = rand(n)
A = rand(n, n)
b = rand(n)
Expand All @@ -73,26 +69,33 @@ sol = solve(prob, KrylovJL_GMRES(precs = (A,p)->(Diagonal(A),I)) )
sol.u
```
This approach has the advantage that the specification of the preconditioner is possible without
the knowledge of a concrete matrix `A`. It also allows to specify the preconditioner via a callable object:
the knowledge of a concrete matrix `A`. It also allows to specify the preconditioner via a callable object
and to pass parameters to the constructor of the preconditioner instances. The example below also shows how
to reuse the preconditioner once constructed for the subsequent solution of a modified problem.

```@example precon2
```@example precon3
using LinearSolve, LinearAlgebra
struct DiagonalPrecs end
Base.@kwdef struct WeightedDiagonalBuilder
w::Float64
end
(::DiagonalPrecs)(A,p) = (Diagonal(A),I)
(builder::WeightedDiagonalBuilder)(A,p) = (builder.w*Diagonal(A),I)
n = 4
s = rand(n)
A = rand(n, n)
A = n*I-rand(n, n)
b = rand(n)
prob = LinearProblem(A, b)
sol = solve(prob, KrylovJL_GMRES(precs = DiagonalPrecs()) )
sol = solve(prob, KrylovJL_GMRES(precs = WeightedDiagonalBuilder(w=0.9)) )
sol.u
```
B=A.+0.1
cache=sol.cache
reinit!(cache,A=B, reuse_precs=true)
sol = solve!(cache, KrylovJL_GMRES(precs = WeightedDiagonalBuilder(w=0.9)) )
sol.u
```
## Preconditioner Interface

To define a new preconditioner you define a Julia type which satisfies the
Expand Down

0 comments on commit 8b8b171

Please sign in to comment.