Merge pull request #580 from ChrisRackauckas-Claude/agent/shooting-explicit-public-api

Remove BoundaryValueDiffEqShooting blanket reexports

Author: ChrisRackauckas

docs/Project.toml

@@ -16,6 +16,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 LinearSolve = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae"
 OptimizationIpopt = "43fad042-7963-4b32-ab19-e2a4f9a67124"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
+OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SimpleBoundaryValueDiffEq = "be0294bd-f90f-4760-ac4e-3421ce2b2da0"
@@ -38,6 +39,7 @@ LinearAlgebra = "1.10"
 LinearSolve = "2.36.2, 3, 4.2, 5.0"
 OptimizationIpopt = "0.2.0, 1.0"
 OrdinaryDiffEq = "7"
+OrdinaryDiffEqTsit5 = "1, 2"
 Plots = "1"
 SciMLBase = "2.60.0, 3"
 SimpleBoundaryValueDiffEq = "1.1.0"

docs/make.jl

@@ -27,12 +27,11 @@ makedocs(;
     modules = [
         BoundaryValueDiffEqCore, BoundaryValueDiffEqMIRK, BoundaryValueDiffEqFIRK,
         BoundaryValueDiffEqMIRKN, BoundaryValueDiffEqShooting, BoundaryValueDiffEqAscher,
-        SciMLBase, DiffEqBase, BoundaryValueDiffEq, SimpleBoundaryValueDiffEq,
+        BoundaryValueDiffEq,
     ],
     clean = true,
-    doctest = false,
+    doctest = true,
     checkdocs = :exports,
-    warnonly = [:missing_docs, :cross_references],
     plugins = [bib, interlinks],
     format = Documenter.HTML(
         assets = ["assets/favicon.ico"],

docs/src/basics/autodiff.md

@@ -12,8 +12,8 @@ BVPJacobianAlgorithm
 
 In BoundaryValueDiffEq.jl, automatic differentiation backend should only be wrapped in `BVPJacobianAlgorithm(diffmode, bc_diffmode, nonbc_diffmode)`. `BVPJacobianAlgorithm(diffmode, bc_diffmode, nonbc_diffmode)` supports user-specified mixed automatic differentiation backends in different part of a boundary value problem, and AD choice should depended on the type of boundary value problem:
 
-  - [`BVProblem`](@ref SciMLBase.BVProblem): Differentiation mode for boundary condition part and non boundary condition part should be specified, for example, `BVPJacobianAlgorithm(; bc_diffmode, nonbc_diffmode)`, default to `BVPJacobianAlgorithm(; bc_diffmode = AutoForwardDiff(), nonbc_diffmode = AutoSparse(AutoForwardDiff()))`.
-  - [`TwoPointBVProblem`](@ref SciMLBase.TwoPointBVProblem): Differentiation mode for overall solving should be specified, for example, `BVPJacobianAlgorithm(; diffmode)`, default to `BVPJacobianAlgorithm(; diffmode = AutoSparse(AutoForwardDiff()))`.
+  - [`BVProblem`](https://docs.sciml.ai/SciMLBase/stable/): Differentiation mode for boundary condition part and non boundary condition part should be specified, for example, `BVPJacobianAlgorithm(; bc_diffmode, nonbc_diffmode)`, default to `BVPJacobianAlgorithm(; bc_diffmode = AutoForwardDiff(), nonbc_diffmode = AutoSparse(AutoForwardDiff()))`.
+  - [`TwoPointBVProblem`](https://docs.sciml.ai/SciMLBase/stable/): Differentiation mode for overall solving should be specified, for example, `BVPJacobianAlgorithm(; diffmode)`, default to `BVPJacobianAlgorithm(; diffmode = AutoSparse(AutoForwardDiff()))`.
 
 In BoundaryValueDiffEq.jl, we require AD to obtain the Jacobian of the loss function which contains the collocation equation and boundary condition equations. For `TwoPointBVProblem`, the Jacobian of the loss function is a sparse banded matrix with known sparsity pattern, but for general multi-points `BVProblem`, the Jacobian of the loss function is an almost banded matrix, which has the first several rows as the boundary conditions and the rest as a sparse banded matrix with known sparsity pattern but without the first several rows. In this case, we can specify mixed AD backend in `BVPJacobianAlgorithm` to make the most of the different sparsity pattern to accelerate BVP solving process.
 

docs/src/basics/bvp_functions.md

@@ -4,7 +4,6 @@ The SciML ecosystem provides an extensive interface for declaring extra function
 
 ## Function Type Definitions
 
-```@docs
-SciMLBase.BVPFunction
-SciMLBase.DynamicalBVPFunction
-```
+`BVPFunction` and `DynamicalBVPFunction` are defined and documented by
+[SciMLBase](https://docs.sciml.ai/SciMLBase/stable/). Use their owner documentation for
+the full problem-function interface.

docs/src/basics/bvp_problem.md

@@ -12,9 +12,7 @@ BoundaryValueDiffEq.jl tackles five related types of boundary value problems:
 
 ## Problem Construction Details
 
-```@docs
-BVProblem
-TwoPointBVProblem
-SecondOrderBVProblem
-TwoPointSecondOrderBVProblem
-```
+`BVProblem`, `TwoPointBVProblem`, `SecondOrderBVProblem`, and
+`TwoPointSecondOrderBVProblem` are defined and documented by
+[SciMLBase](https://docs.sciml.ai/SciMLBase/stable/). Use their owner documentation for
+the complete constructor interfaces.

docs/src/solvers/shooting.md

@@ -7,20 +7,31 @@ using Pkg
 Pkg.add("BoundaryValueDiffEqShooting")
 ```
 
-!!! note "Require OrdinaryDiffEq"
-    
-    Shooting methods require OrdinaryDiffEq.jl loaded to use the ODE solvers
-
-```julia
-solve(prob::BVProblem, alg; kwargs...)
-solve(prob::TwoPointBVProblem, alg; kwargs...)
-```
-
-Shooting methods should be use together with ODE solvers:
-
-```
-BoundaryValueDiffEqShooting.Shooting
-BoundaryValueDiffEqShooting.MultipleShooting
+Shooting algorithms operate on problem definitions and solver functions owned by SciMLBase,
+and require an ODE algorithm from its owning solver package:
+
+```jldoctest
+using BoundaryValueDiffEqShooting: MultipleShooting, Shooting
+using OrdinaryDiffEqTsit5: Tsit5
+using SciMLBase: BVProblem, ReturnCode, solve
+
+function f!(du, u, p, t)
+    du[1] = u[2]
+    du[2] = 0
+end
+
+function bc!(residual, u, p, t)
+    residual[1] = u(0.0)[1] - 1
+    residual[2] = u(1.0)[1]
+end
+
+prob = BVProblem(f!, bc!, [1.0, -1.0], (0.0, 1.0); nlls = Val(false))
+sol = solve(prob, Shooting(Tsit5()); abstol = 1e-8)
+
+@assert sol.retcode == ReturnCode.Success
+@assert isapprox(sol(0.0)[1], 1.0; atol = 1e-6)
+@assert isapprox(sol(1.0)[1], 0.0; atol = 1e-6)
+# output
 ```
 
 ## Full List of Methods

lib/BoundaryValueDiffEqShooting/Project.toml

@@ -17,7 +17,6 @@ OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
 PreallocationTools = "d236fae5-4411-538c-8e31-a6e3d9e00b46"
 PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"
 Preferences = "21216c6a-2e73-6563-6e65-726566657250"
-Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SciMLLogging = "a6db7da4-7206-11f0-1eab-35f2a5dbe1d1"
 Setfield = "efcf1570-3423-57d1-acb7-fd33fddbac46"
@@ -35,6 +34,7 @@ DifferentiationInterface = "0.7.15"
 FastClosures = "0.3.2"
 ForwardDiff = "0.10.38, 1"
 LinearAlgebra = "1.10"
+NonlinearSolveFirstOrder = "1.2, 2"
 OrdinaryDiffEqLowOrderRK = "1, 2"
 OrdinaryDiffEqRosenbrock = "1, 2"
 OrdinaryDiffEqTsit5 = "1, 2"
@@ -44,7 +44,6 @@ Preferences = "1.5"
 PreallocationTools = "1.2"
 Random = "1.10"
 RecursiveArrayTools = "3.27.0, 4"
-Reexport = "1.2"
 SciMLBase = "3.30"
 SciMLLogging = "1.10.1, 2"
 SciMLTesting = "2.4"
@@ -57,6 +56,7 @@ julia = "1.10"
 [extras]
 SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f"
 DiffEqDevTools = "f3b72e0c-5b89-59e1-b016-84e28bfd966d"
+NonlinearSolveFirstOrder = "5959db7a-ea39-4486-b5fe-2dd0bf03d60d"
 OrdinaryDiffEqLowOrderRK = "1344f307-1e59-4825-a18e-ace9aa3fa4c6"
 OrdinaryDiffEqRosenbrock = "43230ef6-c299-4910-a778-202eb28ce4ce"
 OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
@@ -68,4 +68,4 @@ StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
-test = ["DiffEqDevTools", "OrdinaryDiffEqLowOrderRK", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqTsit5", "OrdinaryDiffEqVerner", "Random", "RecursiveArrayTools", "SciMLTesting", "StaticArrays", "Test", "SafeTestsets"]
+test = ["DiffEqDevTools", "NonlinearSolveFirstOrder", "OrdinaryDiffEqLowOrderRK", "OrdinaryDiffEqRosenbrock", "OrdinaryDiffEqTsit5", "OrdinaryDiffEqVerner", "Random", "RecursiveArrayTools", "SciMLTesting", "StaticArrays", "Test", "SafeTestsets"]

lib/BoundaryValueDiffEqShooting/src/BoundaryValueDiffEqShooting.jl

@@ -19,15 +19,14 @@ using BoundaryValueDiffEqCore: BoundaryValueDiffEqCore,
     __extract_u0,
     __initial_guess_on_mesh,
     __get_non_sparse_ad, __build_solution, get_dense_ad,
-    __internal_solve, _process_verbose_param, BVPVerbosity, _unwrap_val
+    __internal_solve, _process_verbose_param, _unwrap_val
 
 using ConcreteStructs: @concrete
 using DifferentiationInterface: DifferentiationInterface,
     overloaded_input_type
 using FastClosures: @closure
 using ForwardDiff: ForwardDiff, pickchunksize
 using LinearAlgebra: LinearAlgebra
-using Reexport: @reexport
 using SciMLBase: SciMLBase, BVProblem, EnsembleSerial, EnsembleThreads,
     NonlinearFunction, ODEProblem, StandardBVProblem, TwoPointBVProblem,
     __solve, isinplace, remake, solve, solve!
@@ -40,8 +39,6 @@ using Preferences: Preferences
 
 const DI = DifferentiationInterface
 
-@reexport using ADTypes, BoundaryValueDiffEqCore, SciMLBase
-
 include("algorithms.jl")
 include("single_shooting.jl")
 include("multiple_shooting.jl")

lib/BoundaryValueDiffEqShooting/src/algorithms.jl

@@ -4,26 +4,52 @@ abstract type AbstractShooting <: AbstractBoundaryValueDiffEqAlgorithm end
 """
     Shooting(ode_alg; kwargs...)
     Shooting(ode_alg, nlsolve; kwargs...)
-    Shooting(; ode_alg = nothing, nlsolve = nothing, optimize = nothing, jac_alg = nothing)
+    Shooting(; ode_alg = nothing, nlsolve = nothing, optimize = nothing, jac_alg = nothing) -> Shooting
 
-Single shooting method, reduces BVP to an initial value problem and solves the IVP.
+Configures the single-shooting algorithm for a boundary value problem. Single shooting
+integrates one initial value problem and solves for the initial condition that satisfies the
+boundary conditions.
 
 ## Arguments
 
-  - `ode_alg`: ODE algorithm to use for solving the IVP. Any solver which conforms to the
-    SciML `ODEProblem` interface can be used! (Defaults to `nothing` which will use
-    poly-algorithm if `DifferentialEquations.jl` is loaded else this must be supplied)
-  - `nlsolve`: Internal Nonlinear solver. Any solver which conforms to the SciML
-    `NonlinearProblem` interface can be used. Note that any autodiff argument for the solver
-    will be ignored and a custom jacobian algorithm will be used.
-  - `optimize`: Internal Optimization solver. Any solver which conforms to the SciML
-    `OptimizationProblem` interface can be used. Note that any autodiff argument for the solver
-    will be ignored and a custom jacobian algorithm will be used.
-  - `jac_alg`: Jacobian Algorithm used for the Nonlinear Solver. If this is not set, we
-    check if `nlsolve.ad` exists and is not nothing. If it is, we use that to construct
-    the jacobian. If not, we try to use the best algorithm based on the input types
-    and problem type. If `BVPJacobianAlgorithm` is provided, only `diffmode` is used
-    (defaults to `AutoForwardDiff` if possible else `AutoFiniteDiff`).
+  - `ode_alg`: algorithm used to solve the internal `SciMLBase.ODEProblem`. Pass this as the
+    first positional argument or keyword argument. `nothing` selects a loaded polyalgorithm;
+    otherwise an ODE algorithm must be supplied.
+  - `nlsolve`: nonlinear-solver algorithm for the shooting residual. Its autodiff setting is
+    superseded by `jac_alg` when a Jacobian algorithm is materialized.
+
+## Keywords
+
+  - `ode_alg = nothing`: ODE algorithm, as described above.
+  - `nlsolve = nothing`: nonlinear-solver algorithm, as described above.
+  - `optimize = nothing`: optimization-solver algorithm used when the selected BVP solve path
+    formulates the residual as an optimization problem.
+  - `jac_alg = nothing`: `BVPJacobianAlgorithm` configuration. When omitted, the constructor
+    derives it from `nlsolve` and the problem during solve initialization. For single shooting,
+    only its `diffmode` setting is used; the default is `AutoForwardDiff` when applicable and
+    otherwise `AutoFiniteDiff`.
+
+## Fields
+
+  - `ode_alg`: configured ODE algorithm or `nothing`.
+  - `nlsolve`: configured nonlinear-solver algorithm or `nothing`.
+  - `optimize`: configured optimization-solver algorithm or `nothing`.
+  - `jac_alg::BVPJacobianAlgorithm`: materialized Jacobian-algorithm configuration.
+
+## Returns
+
+  - `Shooting`: an algorithm object accepted by `SciMLBase.solve` for a boundary value problem.
+
+## Examples
+
+```jldoctest
+using BoundaryValueDiffEqShooting: Shooting
+using OrdinaryDiffEqTsit5: Tsit5
+
+alg = Shooting(Tsit5())
+@assert alg isa Shooting
+# output
+```
 """
 @concrete struct Shooting{J <: BVPJacobianAlgorithm} <: AbstractShooting
     ode_alg
@@ -46,37 +72,40 @@ end
 
 """
     MultipleShooting(; nshoots::Int, ode_alg = nothing, nlsolve = nothing,
-        optimize = nothing, grid_coarsening = true, jac_alg = nothing)
+        optimize = nothing, grid_coarsening = true, jac_alg = nothing) -> MultipleShooting
     MultipleShooting(nshoots::Int; kwargs...)
     MultipleShooting(nshoots::Int, ode_alg; kwargs...)
     MultipleShooting(nshoots::Int, ode_alg, nlsolve; kwargs...)
 
-Multiple Shooting method, reduces BVP to an initial value problem and solves the IVP.
-Significantly more stable than Single Shooting.
+Configures the multiple-shooting algorithm for a boundary value problem. Multiple shooting
+integrates an IVP on `nshoots` subintervals and solves for their matching initial conditions;
+it is generally more stable than [`Shooting`](@ref).
 
 ## Arguments
 
-  - `nshoots`: Number of shooting points.
+  - `nshoots::Int`: number of shooting subintervals.
+  - `ode_alg`: algorithm used to solve each internal `SciMLBase.ODEProblem`. Pass this as the
+    second positional argument or keyword argument. `nothing` selects a loaded polyalgorithm;
+    otherwise an ODE algorithm must be supplied.
+  - `nlsolve`: nonlinear-solver algorithm for the multiple-shooting residual.
 
-  - `ode_alg`: ODE algorithm to use for solving the IVP. Any solver which conforms to the
-    SciML `ODEProblem` interface can be used! (Defaults to `nothing` which will use
-    poly-algorithm if `DifferentialEquations.jl` is loaded else this must be supplied)
-  - `nlsolve`: Internal Nonlinear solver. Any solver which conforms to the SciML
-    `NonlinearProblem` interface can be used.
-  - `optimize`: Internal Optimization solver. Any solver which conforms to the SciML
-    `OptimizationProblem` interface can be used.
-  - `jac_alg`: Jacobian Algorithm used for the nonlinear solver. Defaults to
-    `BVPJacobianAlgorithm()`, which automatically decides the best algorithm to use based
-    on the input types and problem type.
+## Keywords
+
+  - `ode_alg = nothing`: ODE algorithm, as described above.
+  - `nlsolve = nothing`: nonlinear-solver algorithm, as described above.
+  - `optimize = nothing`: optimization-solver algorithm used when the selected BVP solve path
+    formulates the residual as an optimization problem.
+  - `jac_alg = nothing`: `BVPJacobianAlgorithm` configuration. When omitted, the constructor
+    derives it from `nlsolve` and the problem during solve initialization.
 
       + For `TwoPointBVProblem`, only `diffmode` is used (defaults to
         `AutoSparse(AutoForwardDiff())` if possible else `AutoSparse(AutoFiniteDiff())`).
       + For `BVProblem`, `bc_diffmode` and `nonbc_diffmode` are used. For `nonbc_diffmode`
         we default to `AutoSparse(AutoForwardDiff())` if possible else
         `AutoSparse(AutoFiniteDiff())`. For `bc_diffmode`, we default to `AutoForwardDiff`
         if possible else `AutoFiniteDiff`.
-  - `grid_coarsening`: Coarsening the multiple-shooting grid to generate a stable IVP
-    solution. Possible Choices:
+  - `grid_coarsening = true`: coarsens the multiple-shooting grid while generating a stable
+    IVP solution. Supported values are:
 
       + `true`: Halve the grid size, till we reach a grid size of 1.
       + `false`: Do not coarsen the grid. Solve a Multiple Shooting Problem and finally
@@ -87,6 +116,31 @@ Significantly more stable than Single Shooting.
       + `Function`: Takes the current number of shooting points and returns the next number
         of shooting points. For example, if `nshoots = 10` and
         `grid_coarsening = n -> n ÷ 2`, then the grid will be coarsened to `[5, 2]`.
+
+## Fields
+
+  - `ode_alg`: configured ODE algorithm or `nothing`.
+  - `nlsolve`: configured nonlinear-solver algorithm or `nothing`.
+  - `optimize`: configured optimization-solver algorithm or `nothing`.
+  - `jac_alg::BVPJacobianAlgorithm`: materialized Jacobian-algorithm configuration.
+  - `nshoots::Int`: configured number of shooting subintervals.
+  - `grid_coarsening`: configured grid-coarsening strategy.
+
+## Returns
+
+  - `MultipleShooting`: an algorithm object accepted by `SciMLBase.solve` for a boundary value
+    problem.
+
+## Examples
+
+```jldoctest
+using BoundaryValueDiffEqShooting: MultipleShooting
+using OrdinaryDiffEqTsit5: Tsit5
+
+alg = MultipleShooting(8, Tsit5(); grid_coarsening = true)
+@assert alg isa MultipleShooting
+# output
+```
 """
 @concrete struct MultipleShooting{J <: BVPJacobianAlgorithm} <: AbstractShooting
     ode_alg

lib/BoundaryValueDiffEqShooting/test/Core/basic_problems_tests.jl

@@ -1,8 +1,14 @@
 using BoundaryValueDiffEqShooting
+using ADTypes: AutoFiniteDiff, AutoForwardDiff, AutoSparse
+using BoundaryValueDiffEqCore: BVPJacobianAlgorithm
+using NonlinearSolveFirstOrder: NewtonRaphson, TrustRegion
+using SciMLBase: BVPFunction, BVProblem, TwoPointBVProblem, solve
 using Test
 
+import SciMLBase
+
 @testset "Basic Shooting" begin
-    using BoundaryValueDiffEqShooting, LinearAlgebra, OrdinaryDiffEqTsit5
+    using LinearAlgebra, OrdinaryDiffEqTsit5
 
     # JET tests have been moved to the separate QA test group (test/qa/)
 
@@ -111,7 +117,7 @@ using Test
 end
 
 @testset "Shooting with Complex Values" begin
-    using BoundaryValueDiffEqShooting, OrdinaryDiffEqVerner, LinearAlgebra
+    using OrdinaryDiffEqVerner, LinearAlgebra
 
     SOLVERS = [
         Shooting(Vern7(), NewtonRaphson(; autodiff = AutoFiniteDiff())), Shooting(Vern7()),
@@ -148,7 +154,7 @@ end
 end
 
 @testset "Flow in a Channel" begin
-    using BoundaryValueDiffEqShooting, OrdinaryDiffEqTsit5, OrdinaryDiffEqRosenbrock,
+    using OrdinaryDiffEqTsit5, OrdinaryDiffEqRosenbrock,
         LinearAlgebra
 
     function flow_in_a_channel!(du, u, p, t)
@@ -206,7 +212,7 @@ end
 end
 #FIXME: MultipleShooting fails for large out-of-place BVP systems
 @testset "Ray Tracing" begin
-    using BoundaryValueDiffEqShooting, OrdinaryDiffEqVerner, OrdinaryDiffEqRosenbrock,
+    using OrdinaryDiffEqVerner, OrdinaryDiffEqRosenbrock,
         LinearAlgebra
 
     @inline v(x, y, z, p) = 1 / (4 + cos(p[1] * x) + sin(p[2] * y) - cos(p[3] * z))
@@ -380,7 +386,7 @@ end
 end
 
 @testset "Shooting with heterogeneous initial guess" begin
-    using BoundaryValueDiffEqShooting, OrdinaryDiffEqVerner, LinearAlgebra
+    using OrdinaryDiffEqVerner, LinearAlgebra
     g = 9.81
     L = 1.0
     tspan = (0.0, pi / 2.0)