|
| 1 | +# Passing in a Custom Linear Solver |
| 2 | +Julia users are constantly a wide variety of applications in the SciML ecosystem, |
| 3 | +often requiring custom handling. As existing solvers in `LinearSolve.jl` may not |
| 4 | +be optimally suited for novel applications, it is essential for the linear solve |
| 5 | +interface to be easily extendable by users. To that end, the linear solve algorithm |
| 6 | +`LinearSolveFunction()` accepts a user-defined function for handling the solve. A |
| 7 | +user can pass in their custom linear solve function, say `my_linsolve`, to |
| 8 | +`LinearSolveFunction()`. A contrived example of solving a linear system with a custom solver is below. |
| 9 | +```julia |
| 10 | +using LinearSolve, LinearAlgebra |
| 11 | + |
| 12 | +function my_linsolve(A,b,u,p,newA,Pl,Pr,solverdata;verbose=true, kwargs...) |
| 13 | + if verbose == true |
| 14 | + println("solving Ax=b") |
| 15 | + end |
| 16 | + u = A \ b |
| 17 | + return u |
| 18 | +end |
| 19 | + |
| 20 | +prob = LinearProblem(Diagonal(rand(4)), rand(4)) |
| 21 | +alg = LinearSolveFunction(my_linsolve), |
| 22 | +sol = solve(prob, alg) |
| 23 | +``` |
| 24 | +The inputs to the function are as follows: |
| 25 | +- `A`, the linear operator |
| 26 | +- `b`, the right-hand-side |
| 27 | +- `u`, the solution initialized as `zero(b)`, |
| 28 | +- `p`, a set of parameters |
| 29 | +- `newA`, a `Bool` which is `true` if `A` has been modified since last solve |
| 30 | +- `Pl`, left-preconditioner |
| 31 | +- `Pr`, right-preconditioner |
| 32 | +- `solverdata`, solver cache set to `nothing` if solver hasn't been initialized |
| 33 | +- `kwargs`, standard SciML keyword arguments such as `verbose`, `maxiters`, |
| 34 | +`abstol`, `reltol` |
| 35 | +The function `my_linsolve` must accept the above specified arguments, and return |
| 36 | +the solution, `u`. As memory for `u` is already allocated, the user may choose |
| 37 | +to modify `u` in place as follows: |
| 38 | +```julia |
| 39 | +function my_linsolve!(A,b,u,p,newA,Pl,Pr,solverdata;verbose=true, kwargs...) |
| 40 | + if verbose == true |
| 41 | + println("solving Ax=b") |
| 42 | + end |
| 43 | + u .= A \ b # in place |
| 44 | + return u |
| 45 | +end |
| 46 | + |
| 47 | +alg = LinearSolveFunction(my_linsolve!) |
| 48 | +sol = solve(prob, alg) |
| 49 | +``` |
| 50 | +Finally, note that `LinearSolveFunction()` dispatches to the default linear solve |
| 51 | +algorithm handling if no arguments are passed in. |
| 52 | +```julia |
| 53 | +alg = LinearSolveFunction() |
| 54 | +sol = solve(prob, alg) # same as solve(prob, nothing) |
| 55 | +``` |
0 commit comments