From 21495f21c9f495f41842e004fb5b1b784012ef0c Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Thu, 27 Jun 2024 18:08:31 -0400 Subject: [PATCH 01/14] make preconditioners part of the solver rather than a random extra --- Project.toml | 3 ++- src/common.jl | 47 ++++++++++++++++++++++++++++++++++++--- src/iterative_wrappers.jl | 23 +++++++++---------- 3 files changed, 56 insertions(+), 17 deletions(-) diff --git a/Project.toml b/Project.toml index f083c4a79..0e2e53ac0 100644 --- a/Project.toml +++ b/Project.toml @@ -135,6 +135,7 @@ IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153" JET = "c3a54625-cd67-489e-a8e7-0a5a0ff4e31b" KernelAbstractions = "63c18a36-062a-441e-b654-da1e3ab1ce7c" KrylovKit = "0b1a1467-8014-51b9-945f-bf0ae24f4b77" +KrylovPreconditioners = "45d422c2-293f-44ce-8315-2cb988662dec" MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195" Metal = "dde4c033-4e86-420c-a63e-0dd931031962" MultiFloats = "bdf0d083-296b-4888-a5b6-7498122e68a5" @@ -148,4 +149,4 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" [targets] -test = ["Aqua", "Test", "IterativeSolvers", "InteractiveUtils", "JET", "KrylovKit", "Pkg", "Random", "SafeTestsets", "MultiFloats", "ForwardDiff", "HYPRE", "MPI", "BlockDiagonals", "Enzyme", "FiniteDiff", "BandedMatrices", "FastAlmostBandedMatrices", "StaticArrays", "AllocCheck", "StableRNGs", "Zygote"] +test = ["Aqua", "Test", "IterativeSolvers", "InteractiveUtils", "JET", "KrylovKit", "KrylovPreconditioners", "Pkg", "Random", "SafeTestsets", "MultiFloats", "ForwardDiff", "HYPRE", "MPI", "BlockDiagonals", "Enzyme", "FiniteDiff", "BandedMatrices", "FastAlmostBandedMatrices", "StaticArrays", "AllocCheck", "StableRNGs", "Zygote"] diff --git a/src/common.jl b/src/common.jl index a53741411..49708df83 100644 --- a/src/common.jl +++ b/src/common.jl @@ -121,6 +121,8 @@ default_alias_b(::Any, ::Any, ::Any) = false default_alias_A(::AbstractKrylovSubspaceMethod, ::Any, ::Any) = true default_alias_b(::AbstractKrylovSubspaceMethod, ::Any, ::Any) = true +DEFAULT_PRECS(A, p) = IdentityOperator(size(A)[1]), IdentityOperator(size(A)[2]) + function __init_u0_from_Ab(A, b) u0 = similar(b, size(A, 2)) fill!(u0, false) @@ -136,12 +138,12 @@ function SciMLBase.init(prob::LinearProblem, alg::SciMLLinearSolveAlgorithm, reltol = default_tol(real(eltype(prob.b))), maxiters::Int = length(prob.b), verbose::Bool = false, - Pl = IdentityOperator(size(prob.A)[1]), - Pr = IdentityOperator(size(prob.A)[2]), + Pl = nothing, + Pr = nothing, assumptions = OperatorAssumptions(issquare(prob.A)), sensealg = LinearSolveAdjoint(), kwargs...) - @unpack A, b, u0, p = prob + (;A, b, u0, p) = prob A = if alias_A || A isa SMatrix A @@ -167,6 +169,18 @@ function SciMLBase.init(prob::LinearProblem, alg::SciMLLinearSolveAlgorithm, reltol = real(eltype(prob.b))(reltol) abstol = real(eltype(prob.b))(abstol) + precs = hasproperty(alg, :precs) ? alg.precs : DEFAULT_PRECS + _Pl, _Pr = precs(A, p) + if isnothing(Pl) + Pl = _Pl + else + @warn "passing Preconditioners at `init`/`solve` time is deprecated. Instead add a `precs` function to your algorithm." + end + if isnothing(Pr) + Pr = _Pr + else + @warn "passing Preconditioners at `init`/`solve` time is deprecated. Instead add a `precs` function to your algorithm." + end cacheval = init_cacheval(alg, A, b, u0_, Pl, Pr, maxiters, abstol, reltol, verbose, assumptions) isfresh = true @@ -179,6 +193,33 @@ function SciMLBase.init(prob::LinearProblem, alg::SciMLLinearSolveAlgorithm, return cache end + +function SciMLBase.reinit!(cache::LinearCache; + A = nothing, + b = cache.b, + u = cache.u, + p = nothing,) + (; alg, cacheval, isfresh, abstol, reltol, maxiters, verbose, assumptions, sensealg) = cache + + precs = hasproperty(alg, :precs) ? alg.precs : DEFAULT_PRECS + Pl, Pr = if isnothing(A) || isnothing(p) + (cache.Pl, cache.Pr) + else + if isnothing(A) + A = cache.A + end + if isnothing(p) + p = cache.p + end + precs(A, p) + end + + return LinearCache{typeof(A), typeof(b), typeof(u), typeof(p), typeof(alg), typeof(cacheval), + typeof(Pl), typeof(Pr), typeof(reltol), typeof(assumptions.issq), + typeof(sensealg)}(A, b, u, p, alg, cacheval, isfresh, Pl, Pr, abstol, reltol, + maxiters, verbose, assumptions, sensealg) +end + function SciMLBase.solve(prob::LinearProblem, args...; kwargs...) return solve(prob, nothing, args...; kwargs...) end diff --git a/src/iterative_wrappers.jl b/src/iterative_wrappers.jl index bb93ba632..ba2016975 100644 --- a/src/iterative_wrappers.jl +++ b/src/iterative_wrappers.jl @@ -10,19 +10,21 @@ KrylovJL(args...; KrylovAlg = Krylov.gmres!, A generic wrapper over the Krylov.jl krylov-subspace iterative solvers. """ -struct KrylovJL{F, I, A, K} <: AbstractKrylovSubspaceMethod +struct KrylovJL{F, I, A, P, K} <: AbstractKrylovSubspaceMethod KrylovAlg::F gmres_restart::I window::I + precs::P args::A kwargs::K end function KrylovJL(args...; KrylovAlg = Krylov.gmres!, gmres_restart = 0, window = 0, + precs = DEFAULT_PRECS, kwargs...) return KrylovJL(KrylovAlg, gmres_restart, window, - args, kwargs) + precs, args, kwargs) end default_alias_A(::KrylovJL, ::Any, ::Any) = true @@ -231,8 +233,7 @@ function SciMLBase.solve!(cache::LinearCache, alg::KrylovJL; kwargs...) cache.isfresh = false end - M = cache.Pl - N = cache.Pr + M, N = cache.Pl, cache.Pr # use no-op preconditioner for Krylov.jl (LinearAlgebra.I) when M/N is identity M = _isidentity_struct(M) ? I : M @@ -258,25 +259,21 @@ function SciMLBase.solve!(cache::LinearCache, alg::KrylovJL; kwargs...) end args = (cacheval, cache.A, cache.b) - kwargs = (atol = atol, rtol = rtol, itmax = itmax, verbose = verbose, + kwargs = (atol = atol, rtol, itmax, verbose, ldiv = true, history = true, alg.kwargs...) if cache.cacheval isa Krylov.CgSolver N !== I && @warn "$(alg.KrylovAlg) doesn't support right preconditioning." - Krylov.solve!(args...; M = M, - kwargs...) + Krylov.solve!(args...; M, kwargs...) elseif cache.cacheval isa Krylov.GmresSolver - Krylov.solve!(args...; M = M, N = N, restart = alg.gmres_restart > 0, - kwargs...) + Krylov.solve!(args...; M, N, restart = alg.gmres_restart > 0, kwargs...) elseif cache.cacheval isa Krylov.BicgstabSolver - Krylov.solve!(args...; M = M, N = N, - kwargs...) + Krylov.solve!(args...; M, N, kwargs...) elseif cache.cacheval isa Krylov.MinresSolver N !== I && @warn "$(alg.KrylovAlg) doesn't support right preconditioning." - Krylov.solve!(args...; M = M, - kwargs...) + Krylov.solve!(args...; M, kwargs...) else Krylov.solve!(args...; kwargs...) end From d863dd0aece120643713523d657fdc82a2612069 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Fri, 28 Jun 2024 08:52:26 -0400 Subject: [PATCH 02/14] address code review --- src/common.jl | 24 +++++++++++++++++------- 1 file changed, 17 insertions(+), 7 deletions(-) diff --git a/src/common.jl b/src/common.jl index 49708df83..a7deb5b47 100644 --- a/src/common.jl +++ b/src/common.jl @@ -198,13 +198,12 @@ function SciMLBase.reinit!(cache::LinearCache; A = nothing, b = cache.b, u = cache.u, - p = nothing,) + p = nothing, + reinit_cache = false,) (; alg, cacheval, isfresh, abstol, reltol, maxiters, verbose, assumptions, sensealg) = cache precs = hasproperty(alg, :precs) ? alg.precs : DEFAULT_PRECS Pl, Pr = if isnothing(A) || isnothing(p) - (cache.Pl, cache.Pr) - else if isnothing(A) A = cache.A end @@ -212,12 +211,23 @@ function SciMLBase.reinit!(cache::LinearCache; p = cache.p end precs(A, p) + else + (cache.Pl, cache.Pr) end - return LinearCache{typeof(A), typeof(b), typeof(u), typeof(p), typeof(alg), typeof(cacheval), - typeof(Pl), typeof(Pr), typeof(reltol), typeof(assumptions.issq), - typeof(sensealg)}(A, b, u, p, alg, cacheval, isfresh, Pl, Pr, abstol, reltol, - maxiters, verbose, assumptions, sensealg) + if reinit_cache + return LinearCache{typeof(A), typeof(b), typeof(u), typeof(p), typeof(alg), typeof(cacheval), + typeof(Pl), typeof(Pr), typeof(reltol), typeof(assumptions.issq), + typeof(sensealg)}(A, b, u, p, alg, cacheval, isfresh, Pl, Pr, abstol, reltol, + maxiters, verbose, assumptions, sensealg) + else + cache.A = A + cache.b = b + cache.u = u + cache.p = p + cache.Pl = Pl + cache.Pr = Pr + end end function SciMLBase.solve(prob::LinearProblem, args...; kwargs...) From dea2e35ce662e267de62c0ed3ff60e53b7728019 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Mon, 1 Jul 2024 11:47:00 -0400 Subject: [PATCH 03/14] modify extension algs --- Project.toml | 1 + docs/src/advanced/custom.md | 2 +- docs/src/basics/FAQ.md | 11 +++++------ ext/LinearSolveIterativeSolversExt.jl | 4 ++-- ext/LinearSolveKrylovKitExt.jl | 5 +++-- src/extension_algs.jl | 6 ++++-- src/iterative_wrappers.jl | 2 +- 7 files changed, 17 insertions(+), 14 deletions(-) diff --git a/Project.toml b/Project.toml index 0e2e53ac0..d9edebbd8 100644 --- a/Project.toml +++ b/Project.toml @@ -90,6 +90,7 @@ JET = "0.8.28" KLU = "0.6" KernelAbstractions = "0.9.16" Krylov = "0.9" +KrylovPreconditioners = "0.2" KrylovKit = "0.8" LazyArrays = "1.8, 2" Libdl = "1.10" diff --git a/docs/src/advanced/custom.md b/docs/src/advanced/custom.md index 1b6dda391..2ac9acbda 100644 --- a/docs/src/advanced/custom.md +++ b/docs/src/advanced/custom.md @@ -33,7 +33,7 @@ The inputs to the function are as follows: - `p`, a set of parameters - `newA`, a `Bool` which is `true` if `A` has been modified since last solve - `Pl`, left-preconditioner - - `Pr`, right-preconditioner + - `Pr`, right-preconditionerz - `solverdata`, solver cache set to `nothing` if solver hasn't been initialized - `kwargs`, standard SciML keyword arguments such as `verbose`, `maxiters`, `abstol`, `reltol` diff --git a/docs/src/basics/FAQ.md b/docs/src/basics/FAQ.md index 108f748b6..10bec1e2d 100644 --- a/docs/src/basics/FAQ.md +++ b/docs/src/basics/FAQ.md @@ -45,8 +45,8 @@ a few ways: ## How do I use IterativeSolvers solvers with a weighted tolerance vector? -IterativeSolvers.jl computes the norm after the application of the left preconditioner -`Pl`. Thus, in order to use a vector tolerance `weights`, one can mathematically +IterativeSolvers.jl computes the norm after the application of the left preconditioner. +Thus, in order to use a vector tolerance `weights`, one can mathematically hack the system via the following formulation: ```@example FAQPrec @@ -57,11 +57,10 @@ A = rand(n, n) b = rand(n) weights = [1e-1, 1] -Pl = LinearSolve.InvPreconditioner(Diagonal(weights)) -Pr = Diagonal(weights) +precs = Returns(LinearSolve.InvPreconditioner(Diagonal(weights)), Diagonal(weights)) prob = LinearProblem(A, b) -sol = solve(prob, KrylovJL_GMRES(), Pl = Pl, Pr = Pr) +sol = solve(prob, KrylovJL_GMRES(precs)) sol.u ``` @@ -84,5 +83,5 @@ Pl = LinearSolve.ComposePreconditioner(LinearSolve.InvPreconditioner(Diagonal(we Pr = Diagonal(weights) prob = LinearProblem(A, b) -sol = solve(prob, KrylovJL_GMRES(), Pl = Pl, Pr = Pr) +sol = solve(prob, KrylovJL_GMRES(precs=Returns((Pl,Pr)))) ``` diff --git a/ext/LinearSolveIterativeSolversExt.jl b/ext/LinearSolveIterativeSolversExt.jl index 507e75ad2..23ce663c8 100644 --- a/ext/LinearSolveIterativeSolversExt.jl +++ b/ext/LinearSolveIterativeSolversExt.jl @@ -12,9 +12,9 @@ end function LinearSolve.IterativeSolversJL(args...; generate_iterator = IterativeSolvers.gmres_iterable!, - gmres_restart = 0, kwargs...) + gmres_restart = 0, precs = DEFAULT_PRECS, kwargs...) return IterativeSolversJL(generate_iterator, gmres_restart, - args, kwargs) + precs, args, kwargs) end function LinearSolve.IterativeSolversJL_CG(args...; kwargs...) diff --git a/ext/LinearSolveKrylovKitExt.jl b/ext/LinearSolveKrylovKitExt.jl index a26e26688..1aa1e5d52 100644 --- a/ext/LinearSolveKrylovKitExt.jl +++ b/ext/LinearSolveKrylovKitExt.jl @@ -1,12 +1,13 @@ module LinearSolveKrylovKitExt using LinearSolve, KrylovKit, LinearAlgebra -using LinearSolve: LinearCache +using LinearSolve: LinearCache, DEFAULT_PRECS function LinearSolve.KrylovKitJL(args...; KrylovAlg = KrylovKit.GMRES, gmres_restart = 0, + precs = DEFAULT_PRECS, kwargs...) - return KrylovKitJL(KrylovAlg, gmres_restart, args, kwargs) + return KrylovKitJL(KrylovAlg, gmres_restart, precs, args, kwargs) end function LinearSolve.KrylovKitJL_CG(args...; kwargs...) diff --git a/src/extension_algs.jl b/src/extension_algs.jl index c5e2db955..7534d2fa1 100644 --- a/src/extension_algs.jl +++ b/src/extension_algs.jl @@ -259,9 +259,10 @@ solvers. Using this solver requires adding the package KrylovKit.jl, i.e. `using KrylovKit` """ -struct KrylovKitJL{F, A, I, K} <: LinearSolve.AbstractKrylovSubspaceMethod +struct KrylovKitJL{F, I, P, A, K} <: LinearSolve.AbstractKrylovSubspaceMethod KrylovAlg::F gmres_restart::I + precs::P args::A kwargs::K end @@ -306,9 +307,10 @@ A generic wrapper over the IterativeSolvers.jl solvers. Using this solver requires adding the package IterativeSolvers.jl, i.e. `using IterativeSolvers` """ -struct IterativeSolversJL{F, I, A, K} <: LinearSolve.AbstractKrylovSubspaceMethod +struct IterativeSolversJL{F, I, P, A, K} <: LinearSolve.AbstractKrylovSubspaceMethod generate_iterator::F gmres_restart::I + precs::P args::A kwargs::K end diff --git a/src/iterative_wrappers.jl b/src/iterative_wrappers.jl index ba2016975..1ba9067fd 100644 --- a/src/iterative_wrappers.jl +++ b/src/iterative_wrappers.jl @@ -10,7 +10,7 @@ KrylovJL(args...; KrylovAlg = Krylov.gmres!, A generic wrapper over the Krylov.jl krylov-subspace iterative solvers. """ -struct KrylovJL{F, I, A, P, K} <: AbstractKrylovSubspaceMethod +struct KrylovJL{F, I, P, A, K} <: AbstractKrylovSubspaceMethod KrylovAlg::F gmres_restart::I window::I From d60e407b29446ea6b27704d35f53396b68231642 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Mon, 1 Jul 2024 16:59:14 -0400 Subject: [PATCH 04/14] fixes --- ext/LinearSolveIterativeSolversExt.jl | 2 +- src/common.jl | 4 +++- 2 files changed, 4 insertions(+), 2 deletions(-) diff --git a/ext/LinearSolveIterativeSolversExt.jl b/ext/LinearSolveIterativeSolversExt.jl index 23ce663c8..cb4589642 100644 --- a/ext/LinearSolveIterativeSolversExt.jl +++ b/ext/LinearSolveIterativeSolversExt.jl @@ -1,7 +1,7 @@ module LinearSolveIterativeSolversExt using LinearSolve, LinearAlgebra -using LinearSolve: LinearCache +using LinearSolve: LinearCache, DEFAULT_PRECS import LinearSolve: IterativeSolversJL if isdefined(Base, :get_extension) diff --git a/src/common.jl b/src/common.jl index a7deb5b47..f93fa3f77 100644 --- a/src/common.jl +++ b/src/common.jl @@ -200,7 +200,7 @@ function SciMLBase.reinit!(cache::LinearCache; u = cache.u, p = nothing, reinit_cache = false,) - (; alg, cacheval, isfresh, abstol, reltol, maxiters, verbose, assumptions, sensealg) = cache + (; alg, cacheval, abstol, reltol, maxiters, verbose, assumptions, sensealg) = cache precs = hasproperty(alg, :precs) ? alg.precs : DEFAULT_PRECS Pl, Pr = if isnothing(A) || isnothing(p) @@ -214,6 +214,7 @@ function SciMLBase.reinit!(cache::LinearCache; else (cache.Pl, cache.Pr) end + isfresh = true if reinit_cache return LinearCache{typeof(A), typeof(b), typeof(u), typeof(p), typeof(alg), typeof(cacheval), @@ -227,6 +228,7 @@ function SciMLBase.reinit!(cache::LinearCache; cache.p = p cache.Pl = Pl cache.Pr = Pr + cache.isfresh = true end end From 38ac7638dede4df023df305820802dd2c642fd76 Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Mon, 29 Jul 2024 08:40:09 -0400 Subject: [PATCH 05/14] Update docs/src/advanced/custom.md Co-authored-by: Dennis Ogiermann --- docs/src/advanced/custom.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/advanced/custom.md b/docs/src/advanced/custom.md index 2ac9acbda..1b6dda391 100644 --- a/docs/src/advanced/custom.md +++ b/docs/src/advanced/custom.md @@ -33,7 +33,7 @@ The inputs to the function are as follows: - `p`, a set of parameters - `newA`, a `Bool` which is `true` if `A` has been modified since last solve - `Pl`, left-preconditioner - - `Pr`, right-preconditionerz + - `Pr`, right-preconditioner - `solverdata`, solver cache set to `nothing` if solver hasn't been initialized - `kwargs`, standard SciML keyword arguments such as `verbose`, `maxiters`, `abstol`, `reltol` From b80c048b63df59702a35f94848db584a72eeb16e Mon Sep 17 00:00:00 2001 From: Christopher Rackauckas Date: Tue, 30 Jul 2024 01:38:40 -0400 Subject: [PATCH 06/14] Don't default to MKL on EPYC (#518) * Don't default to MKL on EPYC * import correctly * Update Project.toml * Update Project.toml Co-authored-by: Chris Elrod * Update qa.jl * Update static_arrays.jl * Update test/static_arrays.jl * Update test/static_arrays.jl * Update adjoint.jl * Update test/adjoint.jl * Update static_arrays.jl * Update static_arrays.jl * Update test/static_arrays.jl * Update static_arrays.jl * Update test/static_arrays.jl * Update test/static_arrays.jl Co-authored-by: Chris Elrod --------- Co-authored-by: Chris Elrod --- src/LinearSolve.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/LinearSolve.jl b/src/LinearSolve.jl index d8490e86a..db5409049 100644 --- a/src/LinearSolve.jl +++ b/src/LinearSolve.jl @@ -41,6 +41,8 @@ import PrecompileTools import Krylov using SciMLBase import Preferences + import CpuId + const CRC = ChainRulesCore @static if Sys.ARCH === :x86_64 || Sys.ARCH === :i686 From cfac4e7a6684b0134e9a36732a9a1503a7c70f2d Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Mon, 1 Jul 2024 11:47:00 -0400 Subject: [PATCH 07/14] modify extension algs --- docs/src/advanced/custom.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/advanced/custom.md b/docs/src/advanced/custom.md index 1b6dda391..2ac9acbda 100644 --- a/docs/src/advanced/custom.md +++ b/docs/src/advanced/custom.md @@ -33,7 +33,7 @@ The inputs to the function are as follows: - `p`, a set of parameters - `newA`, a `Bool` which is `true` if `A` has been modified since last solve - `Pl`, left-preconditioner - - `Pr`, right-preconditioner + - `Pr`, right-preconditionerz - `solverdata`, solver cache set to `nothing` if solver hasn't been initialized - `kwargs`, standard SciML keyword arguments such as `verbose`, `maxiters`, `abstol`, `reltol` From 06a88a711b3517612e3890aea9712001ed2acb23 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Mon, 5 Aug 2024 17:13:48 -0400 Subject: [PATCH 08/14] rebase fix --- src/LinearSolve.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/src/LinearSolve.jl b/src/LinearSolve.jl index db5409049..c2fb6067c 100644 --- a/src/LinearSolve.jl +++ b/src/LinearSolve.jl @@ -41,7 +41,6 @@ import PrecompileTools import Krylov using SciMLBase import Preferences - import CpuId const CRC = ChainRulesCore From 8c6e383211125f3a5c427f8a8d5212b0a84c1491 Mon Sep 17 00:00:00 2001 From: Oscar Smith Date: Mon, 5 Aug 2024 17:49:10 -0400 Subject: [PATCH 09/14] Update docs/src/advanced/custom.md Co-authored-by: Dennis Ogiermann --- docs/src/advanced/custom.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/advanced/custom.md b/docs/src/advanced/custom.md index 2ac9acbda..1b6dda391 100644 --- a/docs/src/advanced/custom.md +++ b/docs/src/advanced/custom.md @@ -33,7 +33,7 @@ The inputs to the function are as follows: - `p`, a set of parameters - `newA`, a `Bool` which is `true` if `A` has been modified since last solve - `Pl`, left-preconditioner - - `Pr`, right-preconditionerz + - `Pr`, right-preconditioner - `solverdata`, solver cache set to `nothing` if solver hasn't been initialized - `kwargs`, standard SciML keyword arguments such as `verbose`, `maxiters`, `abstol`, `reltol` From 80d4993b574fcdcbe54bd6fe298530585ee0ff86 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Wed, 7 Aug 2024 11:46:01 -0400 Subject: [PATCH 10/14] fix typo --- docs/src/basics/FAQ.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/basics/FAQ.md b/docs/src/basics/FAQ.md index 10bec1e2d..8e07a9957 100644 --- a/docs/src/basics/FAQ.md +++ b/docs/src/basics/FAQ.md @@ -57,7 +57,7 @@ A = rand(n, n) b = rand(n) weights = [1e-1, 1] -precs = Returns(LinearSolve.InvPreconditioner(Diagonal(weights)), Diagonal(weights)) +precs = Returns((LinearSolve.InvPreconditioner(Diagonal(weights)), Diagonal(weights))) prob = LinearProblem(A, b) sol = solve(prob, KrylovJL_GMRES(precs)) From 92a6bd346a636017bb8bcf1a6d1a1759ab08721c Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Wed, 7 Aug 2024 13:26:19 -0400 Subject: [PATCH 11/14] add test --- test/basictests.jl | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/test/basictests.jl b/test/basictests.jl index cb64a1246..d27c0a8dc 100644 --- a/test/basictests.jl +++ b/test/basictests.jl @@ -1,6 +1,6 @@ using LinearSolve, LinearAlgebra, SparseArrays, MultiFloats, ForwardDiff using SciMLOperators -using IterativeSolvers, KrylovKit, MKL_jll +using IterativeSolvers, KrylovKit, MKL_jll, KrylovPreconditioners using Test import Random @@ -267,10 +267,12 @@ end @testset "KrylovJL" begin kwargs = (; gmres_restart = 5) + precs = (A,p=nothing) -> (BlockJacobiPreconditioner(A, 2), I) algorithms = ( ("Default", KrylovJL(kwargs...)), ("CG", KrylovJL_CG(kwargs...)), ("GMRES", KrylovJL_GMRES(kwargs...)), + ("GMRES_prec", KrylovJL_GMRES(;precs, ldiv=false, kwargs...)), # ("BICGSTAB",KrylovJL_BICGSTAB(kwargs...)), ("MINRES", KrylovJL_MINRES(kwargs...)) ) From 665350c43ad4bad83c7d5dea40cad4b609be8ceb Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Wed, 7 Aug 2024 14:22:39 -0400 Subject: [PATCH 12/14] fix setproperty, remove deprecation for now --- src/common.jl | 17 +++++++++++++++-- 1 file changed, 15 insertions(+), 2 deletions(-) diff --git a/src/common.jl b/src/common.jl index f93fa3f77..84c98227c 100644 --- a/src/common.jl +++ b/src/common.jl @@ -85,7 +85,18 @@ end function Base.setproperty!(cache::LinearCache, name::Symbol, x) if name === :A + if hasproperty(cache.alg, :precs) + Pl, Pr = cache.alg.precs(A, cache.p) + setfield!(cache, :Pl, Pl) + setfield!(cache, :Pr, Pr) + end setfield!(cache, :isfresh, true) + elseif name === :p + if hasproperty(cache.alg, :precs) + Pl, Pr = cache.alg.precs(cache.A, p) + setfield!(cache, :Pl, Pl) + setfield!(cache, :Pr, Pr) + end elseif name === :b # In case there is something that needs to be done when b is updated update_cacheval!(cache, :b, x) @@ -174,12 +185,14 @@ function SciMLBase.init(prob::LinearProblem, alg::SciMLLinearSolveAlgorithm, if isnothing(Pl) Pl = _Pl else - @warn "passing Preconditioners at `init`/`solve` time is deprecated. Instead add a `precs` function to your algorithm." + # TODO: deprecate once all docs are updated to the new form + #@warn "passing Preconditioners at `init`/`solve` time is deprecated. Instead add a `precs` function to your algorithm." end if isnothing(Pr) Pr = _Pr else - @warn "passing Preconditioners at `init`/`solve` time is deprecated. Instead add a `precs` function to your algorithm." + # TODO: deprecate once all docs are updated to the new form + #@warn "passing Preconditioners at `init`/`solve` time is deprecated. Instead add a `precs` function to your algorithm." end cacheval = init_cacheval(alg, A, b, u0_, Pl, Pr, maxiters, abstol, reltol, verbose, assumptions) From 3563eada0032a05ddc81e34fdce35f4f0a72bb8e Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Wed, 7 Aug 2024 14:57:55 -0400 Subject: [PATCH 13/14] bugfix --- src/common.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/common.jl b/src/common.jl index 84c98227c..f843727c9 100644 --- a/src/common.jl +++ b/src/common.jl @@ -86,14 +86,14 @@ end function Base.setproperty!(cache::LinearCache, name::Symbol, x) if name === :A if hasproperty(cache.alg, :precs) - Pl, Pr = cache.alg.precs(A, cache.p) + Pl, Pr = cache.alg.precs(x, cache.p) setfield!(cache, :Pl, Pl) setfield!(cache, :Pr, Pr) end setfield!(cache, :isfresh, true) elseif name === :p if hasproperty(cache.alg, :precs) - Pl, Pr = cache.alg.precs(cache.A, p) + Pl, Pr = cache.alg.precs(cache.A, x) setfield!(cache, :Pl, Pl) setfield!(cache, :Pr, Pr) end From a456104866cd943b42cf18906b83e4a1573dbdb7 Mon Sep 17 00:00:00 2001 From: oscarddssmith Date: Wed, 7 Aug 2024 15:58:31 -0400 Subject: [PATCH 14/14] fix ODEs --- src/common.jl | 12 ++++++++---- src/iterative_wrappers.jl | 2 +- 2 files changed, 9 insertions(+), 5 deletions(-) diff --git a/src/common.jl b/src/common.jl index f843727c9..3c222eca5 100644 --- a/src/common.jl +++ b/src/common.jl @@ -85,14 +85,14 @@ end function Base.setproperty!(cache::LinearCache, name::Symbol, x) if name === :A - if hasproperty(cache.alg, :precs) + if hasproperty(cache.alg, :precs) && !isnothing(cache.alg.precs) Pl, Pr = cache.alg.precs(x, cache.p) setfield!(cache, :Pl, Pl) setfield!(cache, :Pr, Pr) end setfield!(cache, :isfresh, true) elseif name === :p - if hasproperty(cache.alg, :precs) + if hasproperty(cache.alg, :precs) && !isnothing(cache.alg.precs) Pl, Pr = cache.alg.precs(cache.A, x) setfield!(cache, :Pl, Pl) setfield!(cache, :Pr, Pr) @@ -180,7 +180,11 @@ function SciMLBase.init(prob::LinearProblem, alg::SciMLLinearSolveAlgorithm, reltol = real(eltype(prob.b))(reltol) abstol = real(eltype(prob.b))(abstol) - precs = hasproperty(alg, :precs) ? alg.precs : DEFAULT_PRECS + precs = if hasproperty(alg, :precs) + isnothing(alg.precs) ? DEFAULT_PRECS : alg.precs + else + DEFAULT_PRECS + end _Pl, _Pr = precs(A, p) if isnothing(Pl) Pl = _Pl @@ -215,7 +219,7 @@ function SciMLBase.reinit!(cache::LinearCache; reinit_cache = false,) (; alg, cacheval, abstol, reltol, maxiters, verbose, assumptions, sensealg) = cache - precs = hasproperty(alg, :precs) ? alg.precs : DEFAULT_PRECS + precs = (hasproperty(alg, :precs) && !isnothing(alg.precs)) ? alg.precs : DEFAULT_PRECS Pl, Pr = if isnothing(A) || isnothing(p) if isnothing(A) A = cache.A diff --git a/src/iterative_wrappers.jl b/src/iterative_wrappers.jl index 1ba9067fd..16a50a27c 100644 --- a/src/iterative_wrappers.jl +++ b/src/iterative_wrappers.jl @@ -21,7 +21,7 @@ end function KrylovJL(args...; KrylovAlg = Krylov.gmres!, gmres_restart = 0, window = 0, - precs = DEFAULT_PRECS, + precs = nothing, kwargs...) return KrylovJL(KrylovAlg, gmres_restart, window, precs, args, kwargs)