From 12c0186c7358395ff84c0459d3c9160bfb31bd94 Mon Sep 17 00:00:00 2001 From: Chris Rackauckas Date: Sun, 4 Dec 2022 12:32:02 +0100 Subject: [PATCH 1/4] Add steady state docs and MINPACK --- docs/pages.jl | 4 +- docs/src/basics/FAQ.md | 2 +- docs/src/basics/NonlinearProblem.md | 28 ++++++-- docs/src/basics/NonlinearSolution.md | 5 ++ docs/src/solvers/NonlinearSystemSolvers.md | 75 +++++++++++++++++++++- docs/src/solvers/SteadyStateSolvers.md | 62 ++++++++++++++++++ 6 files changed, 164 insertions(+), 12 deletions(-) create mode 100644 docs/src/basics/NonlinearSolution.md create mode 100644 docs/src/solvers/SteadyStateSolvers.md diff --git a/docs/pages.jl b/docs/pages.jl index 784411c0f..414f6fc15 100644 --- a/docs/pages.jl +++ b/docs/pages.jl @@ -5,7 +5,9 @@ pages = ["index.md", "tutorials/iterator_interface.md"], "Basics" => Any["basics/NonlinearProblem.md", "basics/NonlinearFunctions.md", + "basics/NonlinearSolution.md", "basics/FAQ.md"], "Solvers" => Any["solvers/NonlinearSystemSolvers.md", - "solvers/BracketingSolvers.md"], + "solvers/BracketingSolvers.md", + "solvers/SteadyStateSolvers.md"], ] diff --git a/docs/src/basics/FAQ.md b/docs/src/basics/FAQ.md index fb106259b..2ead2aedd 100644 --- a/docs/src/basics/FAQ.md +++ b/docs/src/basics/FAQ.md @@ -35,4 +35,4 @@ MATLAB 2022a achieves 1.66s. Try this code yourself: we receive 0.06 seconds, or This example is still not optimized in the Julia code and we expect an improvement in a near future version. -For more information on performance of SciML, see the [SciMLBenchmarks](https://docs.sciml.ai/SciMLBenchmarksOutput/stable/) +For more information on performance of SciML, see the [SciMLBenchmarks](https://docs.sciml.ai/SciMLBenchmarksOutput/stable/). diff --git a/docs/src/basics/NonlinearProblem.md b/docs/src/basics/NonlinearProblem.md index 33cfd7b80..80852039c 100644 --- a/docs/src/basics/NonlinearProblem.md +++ b/docs/src/basics/NonlinearProblem.md @@ -1,32 +1,46 @@ # Nonlinear Problems -## The Two Types of Nonlinear Problems +## The Three Types of Nonlinear Problems NonlinearSolve.jl tackles two related types of nonlinear systems: 1. Interval rootfinding problems. I.e., find the ``t in [t_0, t_f]`` such that `f(t) = 0`. 2. Systems of nonlinear equations, i.e. find the `u` such that `f(u) = 0`. +3. Steady state problems, i.e. find the `u` such that `u' = f(u,t)` has reached steady state, + i.e. `0 = f(u, ∞)`. -The former is for solving scalar rootfinding problems, i.e. finding a single number, and +The first is for solving scalar rootfinding problems, i.e. finding a single number, and requires that a bracketing interval is known. For a bracketing interval, one must have that the sign of `f(t_0)` is opposite the sign of `f(t_f)`, thus guaranteeing a root in the interval. -The latter type of nonlinear system can be multidimensional and thus no ordering nor -boundaries are assumed to be known. For a system of nonlinear equations, `f` can return -an array and the solver seeks to find the value of `u` for which all outputs of `f` are -simultaniously zero. - !!! note Interval rootfinding problems allow for `f` to return an array, in which case the interval rootfinding problem is interpreted as finding the first `t` such that any of the components of the array hit zero. +The second type of nonlinear system can be multidimensional and thus no ordering nor +boundaries are assumed to be known. For a system of nonlinear equations, `f` can return +an array and the solver seeks to find the value of `u` for which all outputs of `f` are +simultaniously zero. + +The last type if equivalent to a nonlinear system but with the extra interpretion of +having a potentially preferred unique root. That is, when there are multiple `u` such +that `f(u) = 0`, the `NonlinearProblem` does not have a preferred solution, while for the +`SteadyStateProblem` the preferred solution is the `u(∞)` that would arise from solving the +ODE `u' = f(u,t)`. + +!!! warn + + Most solvers for `SteadyStateProblem` do not guarentee the preferred solution and + instead will solve for some `u` in the set of solutions. The documentation of the + nonlinear solvers will note if they return the preferred solution. ## Problem Construction Details ```@docs SciMLBase.IntervalNonlinearProblem SciMLBase.NonlinearProblem +SciMLBase.SteadyStateProblem ``` diff --git a/docs/src/basics/NonlinearSolution.md b/docs/src/basics/NonlinearSolution.md new file mode 100644 index 000000000..65df15e86 --- /dev/null +++ b/docs/src/basics/NonlinearSolution.md @@ -0,0 +1,5 @@ +# Nonlinear Solutions + +```@docs +SciMLBase.NonlinearSolution +``` diff --git a/docs/src/solvers/NonlinearSystemSolvers.md b/docs/src/solvers/NonlinearSystemSolvers.md index e16a956cf..5285d470a 100644 --- a/docs/src/solvers/NonlinearSystemSolvers.md +++ b/docs/src/solvers/NonlinearSystemSolvers.md @@ -17,7 +17,7 @@ is an implementation which is specialized for small equations. It is non-allocat static arrays and thus really well-optimized for small systems, thus usually outperforming the other methods when such types are used for `u0`. `NLSolveJL`'s `:trust_region` method can be a good choice for high stability, along with -`CMINPACK`.s +`DynamicSS`. For a system which is very non-stiff (i.e., the condition number of the Jacobian is small, or the eigenvalues of the Jacobian are within a few orders of magnitude), @@ -51,9 +51,48 @@ can be used directly to reduce dependencies and improve load times. very efficient and non-allocating. Thus this implmentation is well-suited for small systems of equations. +### SteadyStateDiffEq.jl + +Note that these solvers do not come by default, and thus one needs to install +the package before using these solvers: + +```julia +]add SteadyStateDiffEq +using SteadyStateDiffEq +``` + +SteadyStateDiffEq.jl uses ODE solvers to iteratively approach the steady state. It is a +very stable method for solving nonlinear systems, though in many cases can be more +computationally expensive than direct methods. + +- `DynamicSS` : Uses an ODE solver to find the steady state. Automatically + terminates when close to the steady state. + `DynamicSS(alg;abstol=1e-8,reltol=1e-6,tspan=Inf)` requires that an + ODE algorithm is given as the first argument. The absolute and + relative tolerances specify the termination conditions on the + derivative's closeness to zero. This internally uses the + `TerminateSteadyState` callback from the Callback Library. The + simulated time for which given ODE is solved can be limited by + `tspan`. If `tspan` is a number, it is equivalent to passing + `(zero(tspan), tspan)`. + +Example usage: + +```julia +using NonlinearSolve, SteadyStateDiffEq, OrdinaryDiffEq +sol = solve(prob,DynamicSS(Tsit5())) + +using Sundials +sol = solve(prob,DynamicSS(CVODE_BDF()),dt=1.0) +``` + +!!! note + + If you use `CVODE_BDF` you may need to give a starting `dt` via `dt=....`.* + ### SciMLNLSolve.jl -This is a wrapper package for importing solvers from other packages into this interface. +This is a wrapper package for importing solvers from NLsolve.jl into this interface. Note that these solvers do not come by default, and thus one needs to install the package before using these solvers: @@ -62,7 +101,6 @@ the package before using these solvers: using SciMLNLSolve ``` -- `CMINPACK()`: A wrapper for using the classic MINPACK method through [MINPACK.jl](https://github.com/sglyon/MINPACK.jl) - `NLSolveJL()`: A wrapper for [NLsolve.jl](https://github.com/JuliaNLSolvers/NLsolve.jl) ```julia @@ -134,3 +172,34 @@ The choices for the linear solver are: - `:KLU`: A sparse factorization method. Requires that the user specify a Jacobian. The Jacobian must be set as a sparse matrix in the `ODEProblem` type. + +### NonlinearSolveMINPACK + +This is a wrapper package for importing solvers from MINPACK.jl into this interface. +Note that these solvers do not come by default, and thus one needs to install +the package before using these solvers: + +```julia +]add NonlinearSolveMINPACK +using NonlinearSolveMINPACK +``` + +- `CMINPACK()`: A wrapper for using the classic MINPACK method through [MINPACK.jl](https://github.com/sglyon/MINPACK.jl) + +```julia +CMINPACK(;show_trace::Bool=false, tracing::Bool=false, method::Symbol=:hybr, + io::IO=stdout) +``` + +The keyword argument `method` can take on different value depending on which method of `fsolve` you are calling. The standard choices of `method` are: + +- `:hybr`: Modified version of Powell's algorithm. Uses MINPACK routine [`hybrd1`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/hybrd1.c) +- `:lm`: Levenberg-Marquardt. Uses MINPACK routine [`lmdif1`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/lmdif1.c) +- `:lmdif`: Advanced Levenberg-Marquardt (more options available with `;kwargs...`). See MINPACK routine [`lmdif`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/lmdif.c) for more information +- `:hybrd`: Advacned modified version of Powell's algorithm (more options available with `;kwargs...`). See MINPACK routine [`hybrd`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/hybrd.c) for more information + +If a Jacobian is supplied as part of the [`NonlinearFunction`](@ref nonlinearfunctions), +then the following methods are allowed: + +- `:hybr`: Advacned modified version of Powell's algorithm with user supplied Jacobian. Additional arguments are available via `;kwargs...`. See MINPACK routine [`hybrj`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/hybrj.c) for more information +- `:lm`: Advanced Levenberg-Marquardt with user supplied Jacobian. Additional arguments are available via `;kwargs...`. See MINPACK routine [`lmder`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/lmder.c) for more information diff --git a/docs/src/solvers/SteadyStateSolvers.md b/docs/src/solvers/SteadyStateSolvers.md new file mode 100644 index 000000000..0baafa15e --- /dev/null +++ b/docs/src/solvers/SteadyStateSolvers.md @@ -0,0 +1,62 @@ +# Steady State Solvers + +`solve(prob::SteadyStateProblem,alg;kwargs)` + +Solves for the steady states in the problem defined by `prob` using the algorithm +`alg`. If no algorithm is given, a default algorithm will be chosen. + +## Recommended Methods + +Conversion to a NonlinearProblem is generally the fastest method. However, this will not +guarantee the preferred root, and thus if the preferred root is required, then it's +recommended that one uses `DynamicSS`. For `DynamicSS`, in many cases an adaptive stiff +solver, like a Rosenbrock or BDF method (`Rodas5` or `QNDF`), is a good way to allow for +very large time steps as the steady state approaches. + +!!! note + + The SteadyStateDiffEq.jl methods on a `SteadyStateProblem` respect the time definition + in the nonlinear definition, i.e. `u' = f(u,t)` uses the correct values for `t` as the + solution evolves. A conversion of a `SteadyStateProblem` to a `NonlinearProblem` + replaces this with the nonlinear system `u' = f(u,∞)`, and thus the direct + `SteadyStateProblem` approach can give different answers (i.e. the correct unique + fixed point) on ODEs with non-autonomous dynamics. + +## Full List of Methods + +### Conversion to NonlinearProblem + +Any `SteadyStateProblem` can be trivially converted to a `NonlinearProblem` via +`NonlinearProblem(prob::SteadyStateProblem)`. Using this appraoch, any of the solvers from +the [Nonlinear System Solvers page](@ref nonlinearsystemsolvers) can be used. + +### SteadyStateDiffEq.jl + +SteadyStateDiffEq.jl uses ODE solvers to iteratively approach the steady state. It is a +very stable method for solving nonlinear systems, though in many cases can be more +computationally expensive than direct methods. + +- `DynamicSS` : Uses an ODE solver to find the steady state. Automatically + terminates when close to the steady state. + `DynamicSS(alg;abstol=1e-8,reltol=1e-6,tspan=Inf)` requires that an + ODE algorithm is given as the first argument. The absolute and + relative tolerances specify the termination conditions on the + derivative's closeness to zero. This internally uses the + `TerminateSteadyState` callback from the Callback Library. The + simulated time for which given ODE is solved can be limited by + `tspan`. If `tspan` is a number, it is equivalent to passing + `(zero(tspan), tspan)`. + +Example usage: + +```julia +using NonlinearSolve, SteadyStateDiffEq, OrdinaryDiffEq +sol = solve(prob,DynamicSS(Tsit5())) + +using Sundials +sol = solve(prob,DynamicSS(CVODE_BDF()),dt=1.0) +``` + +!!! note + + If you use `CVODE_BDF` you may need to give a starting `dt` via `dt=....`.* From de7c3037b5887ac6c085b0663fa28661e750e803 Mon Sep 17 00:00:00 2001 From: Chris Rackauckas Date: Sun, 4 Dec 2022 13:56:42 +0100 Subject: [PATCH 2/4] Setup extended API pages --- docs/Project.toml | 10 ++ docs/make.jl | 6 +- docs/pages.jl | 11 +- docs/src/api/minpack.md | 18 +++ docs/src/api/nlsolve.md | 16 ++ docs/src/api/nonlinearsolve.md | 9 ++ docs/src/api/simplenonlinearsolve.md | 11 ++ docs/src/api/steadystatediffeq.md | 19 +++ docs/src/api/sundials.md | 18 +++ docs/src/index.md | 12 +- docs/src/solvers/BracketingSolvers.md | 2 - docs/src/solvers/NonlinearSystemSolvers.md | 173 ++++----------------- src/raphson.jl | 50 ++++++ 13 files changed, 200 insertions(+), 155 deletions(-) create mode 100644 docs/src/api/minpack.md create mode 100644 docs/src/api/nlsolve.md create mode 100644 docs/src/api/nonlinearsolve.md create mode 100644 docs/src/api/simplenonlinearsolve.md create mode 100644 docs/src/api/steadystatediffeq.md create mode 100644 docs/src/api/sundials.md diff --git a/docs/Project.toml b/docs/Project.toml index 4a1c3650b..c39b612b1 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -2,10 +2,20 @@ BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf" Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec" +NonlinearSolveMINPACK = "c100e077-885d-495a-a2ea-599e143bf69d" +SciMLNLSolve = "e9a6253c-8580-4d32-9898-8661bb511710" +SimpleNonlinearSolve = "727e6d20-b764-4bd8-a329-72de5adea6c7" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" +SteadyStateDiffEq = "9672c7b4-1e72-59bd-8a11-6ac3964bc41f" +Sundials = "c3572dad-4567-51f8-b174-8c6c989267f4" [compat] BenchmarkTools = "1" Documenter = "0.27" NonlinearSolve = "1" +NonlinearSolveMINPACK = "0.1" +SciMLNLSolve = "0.1" +SimpleNonlinearSolve = "0.1" StaticArrays = "1" +SteadyStateDiffEq = "1.10" +Sundials = "4.11" diff --git a/docs/make.jl b/docs/make.jl index 631a2b5c9..a5ece1773 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -1,4 +1,5 @@ -using Documenter, NonlinearSolve +using Documenter, NonlinearSolve, SimpleNonlinearSolve, Sundials, SciMLNLSolve, + NonlinearSolveMINPACK, SteadyStateDiffEq cp("./docs/Manifest.toml", "./docs/src/assets/Manifest.toml", force = true) cp("./docs/Project.toml", "./docs/src/assets/Project.toml", force = true) @@ -7,7 +8,8 @@ include("pages.jl") makedocs(sitename = "NonlinearSolve.jl", authors = "Chris Rackauckas", - modules = [NonlinearSolve, NonlinearSolve.SciMLBase], + modules = [NonlinearSolve, NonlinearSolve.SciMLBase, SimpleNonlinearSolve, + Sundials, SciMLNLSolve, NonlinearSolveMINPACK, SteadyStateDiffEq], clean = true, doctest = false, strict = [ :doctest, diff --git a/docs/pages.jl b/docs/pages.jl index 414f6fc15..58a50bfb2 100644 --- a/docs/pages.jl +++ b/docs/pages.jl @@ -7,7 +7,16 @@ pages = ["index.md", "basics/NonlinearFunctions.md", "basics/NonlinearSolution.md", "basics/FAQ.md"], - "Solvers" => Any["solvers/NonlinearSystemSolvers.md", + "Solver Summaries and Recommendations" => Any[ + "solvers/NonlinearSystemSolvers.md", "solvers/BracketingSolvers.md", "solvers/SteadyStateSolvers.md"], + "Detailed Solver APIs" => Any[ + "nonlinearsolve.md", + "simplenonlinearsolve.md", + "minpack.md", + "nlsolve.md", + "sundials.md", + "steadystatediffeq.md" + ] ] diff --git a/docs/src/api/minpack.md b/docs/src/api/minpack.md new file mode 100644 index 000000000..60ec7b293 --- /dev/null +++ b/docs/src/api/minpack.md @@ -0,0 +1,18 @@ +# MINPACK.jl + +This is a wrapper package for importing solvers from Sundials into the SciML interface. +Note that these solvers do not come by default, and thus one needs to install +the package before using these solvers: + +```julia +]add NonlinearSolveMINPACK +using NonlinearSolveMINPACK +``` + +These methods can be used independently of the rest of NonlinearSolve.jl + +## Solver API + +```@docs +CMINPACK +``` diff --git a/docs/src/api/nlsolve.md b/docs/src/api/nlsolve.md new file mode 100644 index 000000000..2a2077e49 --- /dev/null +++ b/docs/src/api/nlsolve.md @@ -0,0 +1,16 @@ +# NLsolve.jl + +This is a wrapper package for importing solvers from NLsolve.jl into the SciML interface. +Note that these solvers do not come by default, and thus one needs to install +the package before using these solvers: + +```julia +]add SciMLNLSolve +using SciMLNLSolve +``` + +## Solver API + +```@docs +NLsolveJL +``` diff --git a/docs/src/api/nonlinearsolve.md b/docs/src/api/nonlinearsolve.md new file mode 100644 index 000000000..3da6bbc51 --- /dev/null +++ b/docs/src/api/nonlinearsolve.md @@ -0,0 +1,9 @@ +# NonlinearSolve.jl Native Solvers + +These are the native solvers of NonlinearSolve.jl. + +## Solver API + +```@docs +NewtonRaphson +``` diff --git a/docs/src/api/simplenonlinearsolve.md b/docs/src/api/simplenonlinearsolve.md new file mode 100644 index 000000000..a101bb98e --- /dev/null +++ b/docs/src/api/simplenonlinearsolve.md @@ -0,0 +1,11 @@ +# SimpleNonlinearSolve.jl + +These methods can be used independently of the rest of NonlinearSolve.jl + +## Solver API + +```@docs +Bisection +Falsi +SimpleNewtonRaphson +``` diff --git a/docs/src/api/steadystatediffeq.md b/docs/src/api/steadystatediffeq.md new file mode 100644 index 000000000..f190c37fb --- /dev/null +++ b/docs/src/api/steadystatediffeq.md @@ -0,0 +1,19 @@ +# SteadyStateDiffEq.jl + +This is a wrapper package for using ODE solvers from +[DifferentialEquations.jl](https://docs.sciml.ai/DiffEqDocs/stable/) into the SciML interface. +Note that these solvers do not come by default, and thus one needs to install +the package before using these solvers: + +```julia +]add SteadyStateDiffEq +using SteadyStateDiffEq +``` + +These methods can be used independently of the rest of NonlinearSolve.jl + +## Solver API + +```@docs +DynamicSS +``` diff --git a/docs/src/api/sundials.md b/docs/src/api/sundials.md new file mode 100644 index 000000000..e3b6b7840 --- /dev/null +++ b/docs/src/api/sundials.md @@ -0,0 +1,18 @@ +## Sundials.jl + +This is a wrapper package for importing solvers from Sundials into the SciML interface. +Note that these solvers do not come by default, and thus one needs to install +the package before using these solvers: + +```julia +]add Sundials +using Sundials +``` + +These methods can be used independently of the rest of NonlinearSolve.jl + +## Solver API + +```@docs +KINSOL +``` diff --git a/docs/src/index.md b/docs/src/index.md index 35419a3c2..c10c1dde2 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -7,8 +7,8 @@ ability to use sparse automatic differentiation for Jacobian construction and Jacobian-vector products. It interfaces with other packages of the Julia ecosystem to make it easy to test alternative solver packages and pass small types to control algorithm swapping. It also interfaces with the -[ModelingToolkit.jl](https://docs.sciml.ai/ModelingToolkit/stable/) world of symbolic modeling to -allow for automatically generating high-performance code. +[ModelingToolkit.jl](https://docs.sciml.ai/ModelingToolkit/stable/) world of symbolic +modeling to allow for automatically generating high-performance code. Performance is key: the current methods are made to be highly performant on scalar and statically sized small problems, with options for large-scale systems. @@ -39,14 +39,8 @@ Pkg.add("NonlinearSolve") - On the [Julia Discourse forums](https://discourse.julialang.org) - See also [SciML Community page](https://sciml.ai/community/) -## Roadmap - -The current algorithms should support automatic differentiation, though improved -adjoint overloads are planned to be added in the current update (which will make -use of the `f(u,p)` form). Future updates will include standard methods for -larger scale nonlinear solving like Newton-Krylov methods. - ## Reproducibility + ```@raw html
The documentation of this SciML package was built using these direct dependencies, ``` diff --git a/docs/src/solvers/BracketingSolvers.md b/docs/src/solvers/BracketingSolvers.md index 0dd6c4375..d043bdd8e 100644 --- a/docs/src/solvers/BracketingSolvers.md +++ b/docs/src/solvers/BracketingSolvers.md @@ -5,8 +5,6 @@ Solves for ``f(t)=0`` in the problem defined by `prob` using the algorithm `alg`. If no algorithm is given, a default algorithm will be chosen. -This page is solely focused on the bracketing methods for scalar nonlinear equations. - ## Recommended Methods `Falsi()` can have a faster convergence and is discretely differentiable, but is diff --git a/docs/src/solvers/NonlinearSystemSolvers.md b/docs/src/solvers/NonlinearSystemSolvers.md index 5285d470a..8e34549e5 100644 --- a/docs/src/solvers/NonlinearSystemSolvers.md +++ b/docs/src/solvers/NonlinearSystemSolvers.md @@ -5,19 +5,16 @@ Solves for ``f(u)=0`` in the problem defined by `prob` using the algorithm `alg`. If no algorithm is given, a default algorithm will be chosen. -This page is solely focused on the methods for nonlinear systems. - ## Recommended Methods `NewtonRaphson` is a good choice for most problems. For large systems it can make use of sparsity patterns for sparse automatic differentiation and sparse linear solving of very large systems. That said, as a classic Newton -method, its stability region can be smaller than other methods. Meanwhile, `SimpleNewtonRaphson` -is an implementation which is specialized for small equations. It is non-allocating on -static arrays and thus really well-optimized for small systems, thus usually outperforming -the other methods when such types are used for `u0`. `NLSolveJL`'s -`:trust_region` method can be a good choice for high stability, along with -`DynamicSS`. +method, its stability region can be smaller than other methods. Meanwhile, +`SimpleNewtonRaphson` is an implementation which is specialized for small equations. It is +non-allocating on static arrays and thus really well-optimized for small systems, thus +usually outperforming the other methods when such types are used for `u0`. `NLSolveJL`'s +`:trust_region` method can be a good choice for high stability, along with `DynamicSS`. For a system which is very non-stiff (i.e., the condition number of the Jacobian is small, or the eigenvalues of the Jacobian are within a few orders of magnitude), @@ -25,26 +22,26 @@ then `NLSolveJL`'s `:anderson` can be a good choice. ## Full List of Methods -### NonlinearSolve.jl +!!! note -These are the core solvers. + For the full details on the capabilities and constructors of the different solvers + see the Detailed Solver APIs section! -- `NewtonRaphson()`: - A Newton-Raphson method with swappable nonlinear solvers and autodiff methods - for high performance on large and sparse systems. +### NonlinearSolve.jl -#### Details on Controlling NonlinearSolve.jl Solvers +These are the core solvers. These methods excel for large-scale problems that need advanced +linear solver, automatic differentiation, abstract array types, GPU, +sparse/structured matrix support, etc. These methods support the largest set of types and +features, but have a bit of overhead on very small problems. -```julia -NewtonRaphson(; chunk_size = Val{0}(), autodiff = Val{true}(), - standardtag = Val{true}(), concrete_jac = nothing, - diff_type = Val{:forward}, linsolve = nothing, precs = DEFAULT_PRECS) -``` +- `NewtonRaphson()`:A Newton-Raphson method with swappable nonlinear solvers and autodiff + methods for high performance on large and sparse systems. ### SimpleNonlinearSolve.jl These methods are included with NonlinearSolve.jl by default, though SimpleNonlinearSolve.jl -can be used directly to reduce dependencies and improve load times. +can be used directly to reduce dependencies and improve load times. SimpleNonlinearSolve.jl's +methods excell at small problems and problems defined with static arrays. - `SimpleNewtonRaphson()`: A simplified implementation of the Newton-Raphson method. Has the property that when used with states `u` as a `Number` or `StaticArray`, the solver is @@ -53,73 +50,20 @@ can be used directly to reduce dependencies and improve load times. ### SteadyStateDiffEq.jl -Note that these solvers do not come by default, and thus one needs to install -the package before using these solvers: - -```julia -]add SteadyStateDiffEq -using SteadyStateDiffEq -``` - SteadyStateDiffEq.jl uses ODE solvers to iteratively approach the steady state. It is a very stable method for solving nonlinear systems, though in many cases can be more computationally expensive than direct methods. - `DynamicSS` : Uses an ODE solver to find the steady state. Automatically terminates when close to the steady state. - `DynamicSS(alg;abstol=1e-8,reltol=1e-6,tspan=Inf)` requires that an - ODE algorithm is given as the first argument. The absolute and - relative tolerances specify the termination conditions on the - derivative's closeness to zero. This internally uses the - `TerminateSteadyState` callback from the Callback Library. The - simulated time for which given ODE is solved can be limited by - `tspan`. If `tspan` is a number, it is equivalent to passing - `(zero(tspan), tspan)`. - -Example usage: - -```julia -using NonlinearSolve, SteadyStateDiffEq, OrdinaryDiffEq -sol = solve(prob,DynamicSS(Tsit5())) - -using Sundials -sol = solve(prob,DynamicSS(CVODE_BDF()),dt=1.0) -``` - -!!! note - - If you use `CVODE_BDF` you may need to give a starting `dt` via `dt=....`.* ### SciMLNLSolve.jl -This is a wrapper package for importing solvers from NLsolve.jl into this interface. -Note that these solvers do not come by default, and thus one needs to install -the package before using these solvers: - -```julia -]add SciMLNLSolve -using SciMLNLSolve -``` +This is a wrapper package for importing solvers from NLsolve.jl into the SciML interface. - `NLSolveJL()`: A wrapper for [NLsolve.jl](https://github.com/JuliaNLSolvers/NLsolve.jl) -```julia -NLSolveJL(; - method=:trust_region, - autodiff=:central, - store_trace=false, - extended_trace=false, - linesearch=LineSearches.Static(), - linsolve=(x, A, b) -> copyto!(x, A\b), - factor = one(Float64), - autoscale=true, - m=10, - beta=one(Float64), - show_trace=false, - ) -``` - -Choices for methods in `NLSolveJL`: +Submethod choices for this algorithm include: - `:fixedpoint`: Fixed-point iteration - `:anderson`: Anderson-accelerated fixed-point iteration @@ -129,77 +73,24 @@ Choices for methods in `NLSolveJL`: For more information on these arguments, consult the [NLsolve.jl documentation](https://github.com/JuliaNLSolvers/NLsolve.jl). -### Sundials.jl +### MINPACK.jl -This is a wrapper package for the SUNDIALS C library, specifically the KINSOL -nonlinear solver included in that ecosystem. Note that these solvers do not come -by default, and thus one needs to install the package before using these solvers: - -```julia -]add Sundials -using Sundials -``` - -- `KINSOL`: The KINSOL method of the SUNDIALS C library - -```julia -KINSOL(; - linear_solver = :Dense, - jac_upper = 0, - jac_lower = 0, - userdata = nothing, -) -``` - -The choices for the linear solver are: - -- `:Dense`: A dense linear solver -- `:Band`: A solver specialized for banded Jacobians. If used, you must set the - position of the upper and lower non-zero diagonals via `jac_upper` and - `jac_lower`. -- `:LapackDense`: A version of the dense linear solver that uses the Julia-provided - OpenBLAS-linked LAPACK for multithreaded operations. This will be faster than - `:Dense` on larger systems but has noticeable overhead on smaller (<100 ODE) systems. -- `:LapackBand`: A version of the banded linear solver that uses the Julia-provided - OpenBLAS-linked LAPACK for multithreaded operations. This will be faster than - `:Band` on larger systems but has noticeable overhead on smaller (<100 ODE) systems. -- `:Diagonal`: This method is specialized for diagonal Jacobians. -- `:GMRES`: A GMRES method. Recommended first choice Krylov method. -- `:BCG`: A biconjugate gradient method -- `:PCG`: A preconditioned conjugate gradient method. Only for symmetric - linear systems. -- `:TFQMR`: A TFQMR method. -- `:KLU`: A sparse factorization method. Requires that the user specify a - Jacobian. The Jacobian must be set as a sparse matrix in the `ODEProblem` - type. - -### NonlinearSolveMINPACK - -This is a wrapper package for importing solvers from MINPACK.jl into this interface. -Note that these solvers do not come by default, and thus one needs to install -the package before using these solvers: - -```julia -]add NonlinearSolveMINPACK -using NonlinearSolveMINPACK -``` +MINPACK.jl methods are good for medium-sized nonlinear solves. It does not scale due to +the lack of sparse Jacobian support, though the methods are very robust and stable. - `CMINPACK()`: A wrapper for using the classic MINPACK method through [MINPACK.jl](https://github.com/sglyon/MINPACK.jl) -```julia -CMINPACK(;show_trace::Bool=false, tracing::Bool=false, method::Symbol=:hybr, - io::IO=stdout) -``` +Submethod choices for this algorithm include: -The keyword argument `method` can take on different value depending on which method of `fsolve` you are calling. The standard choices of `method` are: +- `:hybr`: Modified version of Powell's algorithm. +- `:lm`: Levenberg-Marquardt. +- `:lmdif`: Advanced Levenberg-Marquardt +- `:hybrd`: Advacned modified version of Powell's algorithm -- `:hybr`: Modified version of Powell's algorithm. Uses MINPACK routine [`hybrd1`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/hybrd1.c) -- `:lm`: Levenberg-Marquardt. Uses MINPACK routine [`lmdif1`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/lmdif1.c) -- `:lmdif`: Advanced Levenberg-Marquardt (more options available with `;kwargs...`). See MINPACK routine [`lmdif`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/lmdif.c) for more information -- `:hybrd`: Advacned modified version of Powell's algorithm (more options available with `;kwargs...`). See MINPACK routine [`hybrd`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/hybrd.c) for more information +### Sundials.jl -If a Jacobian is supplied as part of the [`NonlinearFunction`](@ref nonlinearfunctions), -then the following methods are allowed: +Sundials.jl are a classic set of C/Fortran methods which are known for good scaling of the +Newton-Krylov form. However, KINSOL is known to be less stable than some of the other +implementations as it has no line search or globalizer (trust region). -- `:hybr`: Advacned modified version of Powell's algorithm with user supplied Jacobian. Additional arguments are available via `;kwargs...`. See MINPACK routine [`hybrj`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/hybrj.c) for more information -- `:lm`: Advanced Levenberg-Marquardt with user supplied Jacobian. Additional arguments are available via `;kwargs...`. See MINPACK routine [`lmder`](https://github.com/devernay/cminpack/blob/d1f5f5a273862ca1bbcf58394e4ac060d9e22c76/lmder.c) for more information +- `KINSOL()`: The KINSOL method of the SUNDIALS C library diff --git a/src/raphson.jl b/src/raphson.jl index 9e7ae6c88..365726fcc 100644 --- a/src/raphson.jl +++ b/src/raphson.jl @@ -1,3 +1,53 @@ +""" +```julia +NewtonRaphson(; chunk_size = Val{0}(), autodiff = Val{true}(), + standardtag = Val{true}(), concrete_jac = nothing, + diff_type = Val{:forward}, linsolve = nothing, precs = DEFAULT_PRECS) +``` + +An advanced NewtonRaphson implementation with support for efficient handling of sparse +matrices via colored automatic differentiation and preconditioned linear solvers. Designed +for large-scale and numerically-difficult nonlinear systems. + +### Keyword Arguments + +- `chunk_size`: the chunk size used by the internal ForwardDiff.jl automatic differentiation + system. This allows for multiple derivative columns to be computed simultaniously, + improving performance. Defaults to `0`, which is equivalent to using ForwardDiff.jl's + default chunk size mechanism. For more details, see the documentation for + [ForwardDiff.jl](https://juliadiff.org/ForwardDiff.jl/stable/). +- `autodiff`: whether to use forward-mode automatic differentiation for the Jacobian. + Note that this argument is ignored if an analytical Jacobian is passed as that will be + used instead. Defaults to `Val{true}`, which means ForwardDiff.jl via + SparseDiffTools.jl is used by default. If `Val{false}`, then FiniteDiff.jl is used for + finite differencing. +- `standardtag`: whether to use a standardized tag definition for the purposes of automatic + differentiation. Defaults to true, which thus uses the `NonlinearSolveTag`. If `Val{false}`, + then ForwardDiff's default function naming tag is used, which results in larger stack + traces. +- `concrete_jac`: whether to build a concrete Jacobian. If a Krylov-subspace method is used, + then the Jacobian will not be constructed and instead direct Jacobian-vector products + `J*v` are computed using forward-mode automatic differentiation or finite differencing + tricks (without ever constructing the Jacobian). However, if the Jacobian is still needed, + for example for a preconditioner, `concrete_jac = true` can be passed in order to force + the construction of the Jacobian. +- `diff_type`: the type of finite differencing used if `autodiff = false`. Defaults to + `Val{:forward}` for forward finite differences. For more details on the choices, see the + [FiniteDiff.jl](https://github.com/JuliaDiff/FiniteDiff.jl) documentation. +- `linsolve`: the [LinearSolve.jl](https://github.com/SciML/LinearSolve.jl) used for the + linear solves within the Newton method. Defaults to `nothing`, which menans it uses the + LinearSolve.jl default algorithm choice. For more information on available algorithm + choices, see the [LinearSolve.jl documentation](https://docs.sciml.ai/LinearSolve/stable/) +- `precs`: the choice of preconditioners for the linear solver. Defaults to using no + preconditioners. For more information on specifying preconditioners for LinearSolve + algorithms, consult the + [LinearSolve.jl documentation](https://docs.sciml.ai/LinearSolve/stable/) + +!!! note + + Currently the linear solver and chunk size choice only applies to in-place defined + `NonlinearProblem`s. That is expected to change in the future. +""" struct NewtonRaphson{CS, AD, FDT, L, P, ST, CJ} <: AbstractNewtonAlgorithm{CS, AD, FDT, ST, CJ} linsolve::L From de589103e0df4d915cfc35114470419072e0619a Mon Sep 17 00:00:00 2001 From: Chris Rackauckas Date: Sun, 4 Dec 2022 13:57:25 +0100 Subject: [PATCH 3/4] simplify --- docs/src/solvers/NonlinearSystemSolvers.md | 3 --- 1 file changed, 3 deletions(-) diff --git a/docs/src/solvers/NonlinearSystemSolvers.md b/docs/src/solvers/NonlinearSystemSolvers.md index 8e34549e5..c84ae06cd 100644 --- a/docs/src/solvers/NonlinearSystemSolvers.md +++ b/docs/src/solvers/NonlinearSystemSolvers.md @@ -70,9 +70,6 @@ Submethod choices for this algorithm include: - `:newton`: Classical Newton method with an optional line search - `:trust_region`: Trust region Newton method (the default choice) -For more information on these arguments, consult the -[NLsolve.jl documentation](https://github.com/JuliaNLSolvers/NLsolve.jl). - ### MINPACK.jl MINPACK.jl methods are good for medium-sized nonlinear solves. It does not scale due to From 0000f7d44a0e1fcb6383b8870dda0b2486a172b2 Mon Sep 17 00:00:00 2001 From: Chris Rackauckas Date: Sun, 4 Dec 2022 14:01:57 +0100 Subject: [PATCH 4/4] format and fix links --- docs/make.jl | 2 +- docs/pages.jl | 21 +++++++++------------ 2 files changed, 10 insertions(+), 13 deletions(-) diff --git a/docs/make.jl b/docs/make.jl index a5ece1773..01b0f78a3 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -9,7 +9,7 @@ include("pages.jl") makedocs(sitename = "NonlinearSolve.jl", authors = "Chris Rackauckas", modules = [NonlinearSolve, NonlinearSolve.SciMLBase, SimpleNonlinearSolve, - Sundials, SciMLNLSolve, NonlinearSolveMINPACK, SteadyStateDiffEq], + Sundials, SciMLNLSolve, NonlinearSolveMINPACK, SteadyStateDiffEq], clean = true, doctest = false, strict = [ :doctest, diff --git a/docs/pages.jl b/docs/pages.jl index 58a50bfb2..6f0a7445c 100644 --- a/docs/pages.jl +++ b/docs/pages.jl @@ -7,16 +7,13 @@ pages = ["index.md", "basics/NonlinearFunctions.md", "basics/NonlinearSolution.md", "basics/FAQ.md"], - "Solver Summaries and Recommendations" => Any[ - "solvers/NonlinearSystemSolvers.md", - "solvers/BracketingSolvers.md", - "solvers/SteadyStateSolvers.md"], - "Detailed Solver APIs" => Any[ - "nonlinearsolve.md", - "simplenonlinearsolve.md", - "minpack.md", - "nlsolve.md", - "sundials.md", - "steadystatediffeq.md" - ] + "Solver Summaries and Recommendations" => Any["solvers/NonlinearSystemSolvers.md", + "solvers/BracketingSolvers.md", + "solvers/SteadyStateSolvers.md"], + "Detailed Solver APIs" => Any["api/nonlinearsolve.md", + "api/simplenonlinearsolve.md", + "api/minpack.md", + "api/nlsolve.md", + "api/sundials.md", + "api/steadystatediffeq.md"], ]