From 79e3f39fec8be459e983bd2d8bcf71723114e716 Mon Sep 17 00:00:00 2001 From: Qingyu Qu <2283984853@qq.com> Date: Sun, 23 Mar 2025 15:23:27 +0800 Subject: [PATCH 1/4] More generic safe similar --- lib/NonlinearSolveBase/src/jacobian.jl | 6 +++--- lib/NonlinearSolveBase/src/utils.jl | 8 ++++---- 2 files changed, 7 insertions(+), 7 deletions(-) diff --git a/lib/NonlinearSolveBase/src/jacobian.jl b/lib/NonlinearSolveBase/src/jacobian.jl index 363b7b77a..8e23b6632 100644 --- a/lib/NonlinearSolveBase/src/jacobian.jl +++ b/lib/NonlinearSolveBase/src/jacobian.jl @@ -41,7 +41,7 @@ function construct_jacobian_cache( (linsolve === nothing || needs_concrete_A(linsolve))) needs_jac = linsolve_needs_jac || concrete_jac(alg) - @bb fu_cache = similar(fu) + fu_cache = Utils.safe_similar(fu) if !has_analytic_jac && needs_jac if autodiff === nothing @@ -80,9 +80,9 @@ function construct_jacobian_cache( end else if eltype(f.jac_prototype) <: Bool - similar(f.jac_prototype, promote_type(eltype(fu), eltype(u))) + Utils.safe_similar(f.jac_prototype, promote_type(eltype(fu), eltype(u))) else - similar(f.jac_prototype) + Utils.safe_similar(f.jac_prototype) end end end diff --git a/lib/NonlinearSolveBase/src/utils.jl b/lib/NonlinearSolveBase/src/utils.jl index 99e8e50c5..12db8d961 100644 --- a/lib/NonlinearSolveBase/src/utils.jl +++ b/lib/NonlinearSolveBase/src/utils.jl @@ -101,13 +101,13 @@ restructure(y, x) = ArrayInterface.restructure(y, x) function safe_similar(x, args...; kwargs...) y = similar(x, args...; kwargs...) - return init_bigfloat_array!!(y) + return init_similar_array!!(y) end -init_bigfloat_array!!(x) = x +init_similar_array!!(x) = x -function init_bigfloat_array!!(x::AbstractArray{<:BigFloat}) - ArrayInterface.can_setindex(x) && fill!(x, BigFloat(0)) +function init_similar_array!!(x::AbstractArray{<:T}) where {T} + ArrayInterface.can_setindex(x) && fill!(x, T(0)) return x end From 54f981c13c7ee5b55815c4b7547e1ba967e2c79c Mon Sep 17 00:00:00 2001 From: Qingyu Qu <2283984853@qq.com> Date: Sun, 23 Mar 2025 21:50:43 +0800 Subject: [PATCH 2/4] Fix docs failure --- docs/src/basics/faq.md | 46 +++++++++---------- docs/src/tutorials/iterator_interface.md | 2 +- .../src/NonlinearSolveBase.jl | 3 +- lib/NonlinearSolveBase/src/linear_solve.jl | 3 +- lib/SimpleNonlinearSolve/src/utils.jl | 2 +- test/qa_tests.jl | 3 +- 6 files changed, 31 insertions(+), 28 deletions(-) diff --git a/docs/src/basics/faq.md b/docs/src/basics/faq.md index ca35771c6..8411477fe 100644 --- a/docs/src/basics/faq.md +++ b/docs/src/basics/faq.md @@ -71,29 +71,29 @@ differentiate the function based on the input types. However, this function has `xx = [1.0, 2.0, 3.0, 4.0]` followed by a `xx[1] = var[1] - v_true[1]` where `var` might be a Dual number. This causes the error. To fix it: -1. Specify the `autodiff` to be `AutoFiniteDiff` - - ```@example dual_error_faq - sol = solve(prob_oop, LevenbergMarquardt(; autodiff = AutoFiniteDiff()); - maxiters = 10000, abstol = 1e-8) - ``` - - This worked but, Finite Differencing is not the recommended approach in any scenario. - -2. Rewrite the function to use - [PreallocationTools.jl](https://github.com/SciML/PreallocationTools.jl) or write it as - - ```@example dual_error_faq - function fff_correct(var, p) - v_true = [1.0, 0.1, 2.0, 0.5] - xx = eltype(var)[1.0, 2.0, 3.0, 4.0] - xx[1] = var[1] - v_true[1] - return xx - v_true - end - - prob_oop = NonlinearLeastSquaresProblem{false}(fff_correct, v_init) - sol = solve(prob_oop, LevenbergMarquardt(); maxiters = 10000, abstol = 1e-8) - ``` + 1. Specify the `autodiff` to be `AutoFiniteDiff` + + ```@example dual_error_faq + sol = solve(prob_oop, LevenbergMarquardt(; autodiff = AutoFiniteDiff()); + maxiters = 10000, abstol = 1e-8) + ``` + + This worked but, Finite Differencing is not the recommended approach in any scenario. + + 2. Rewrite the function to use + [PreallocationTools.jl](https://github.com/SciML/PreallocationTools.jl) or write it as + + ```@example dual_error_faq + function fff_correct(var, p) + v_true = [1.0, 0.1, 2.0, 0.5] + xx = eltype(var)[1.0, 2.0, 3.0, 4.0] + xx[1] = var[1] - v_true[1] + return xx - v_true + end + + prob_oop = NonlinearLeastSquaresProblem{false}(fff_correct, v_init) + sol = solve(prob_oop, LevenbergMarquardt(); maxiters = 10000, abstol = 1e-8) + ``` ## I thought NonlinearSolve.jl was type-stable and fast. But it isn't, why? diff --git a/docs/src/tutorials/iterator_interface.md b/docs/src/tutorials/iterator_interface.md index 1b6aee101..26873fe35 100644 --- a/docs/src/tutorials/iterator_interface.md +++ b/docs/src/tutorials/iterator_interface.md @@ -20,7 +20,7 @@ to iterate the solver. The iterator interface supports: ```@docs -step!(nlcache::NonlinearSolve.AbstractNonlinearSolveCache, args...; kwargs...) +step!(nlcache::NonlinearSolveBase.AbstractNonlinearSolveCache, args...; kwargs...) ``` We can perform 10 steps of the Newton-Raphson solver with the following: diff --git a/lib/NonlinearSolveBase/src/NonlinearSolveBase.jl b/lib/NonlinearSolveBase/src/NonlinearSolveBase.jl index dd522829f..d86718329 100644 --- a/lib/NonlinearSolveBase/src/NonlinearSolveBase.jl +++ b/lib/NonlinearSolveBase/src/NonlinearSolveBase.jl @@ -19,7 +19,8 @@ using RecursiveArrayTools: AbstractVectorOfArray, ArrayPartition using SciMLBase: SciMLBase, ReturnCode, AbstractODEIntegrator, AbstractNonlinearProblem, AbstractNonlinearAlgorithm, AbstractNonlinearFunction, NonlinearProblem, NonlinearLeastSquaresProblem, StandardNonlinearProblem, - NonlinearFunction, NullParameters, NLStats, LinearProblem, LinearAliasSpecifier + NonlinearFunction, NullParameters, NLStats, LinearProblem, + LinearAliasSpecifier using SciMLJacobianOperators: JacobianOperator, StatefulJacobianOperator using SciMLOperators: AbstractSciMLOperator, IdentityOperator using SymbolicIndexingInterface: SymbolicIndexingInterface diff --git a/lib/NonlinearSolveBase/src/linear_solve.jl b/lib/NonlinearSolveBase/src/linear_solve.jl index d97e047e2..a6f0f34d9 100644 --- a/lib/NonlinearSolveBase/src/linear_solve.jl +++ b/lib/NonlinearSolveBase/src/linear_solve.jl @@ -74,7 +74,8 @@ function construct_linear_solver(alg, linsolve, A, b, u; stats, kwargs...) linprob = LinearProblem(A, b; u0 = u_cache, kwargs...) # unlias here, we will later use these as caches - lincache = init(linprob, linsolve; alias = LinearAliasSpecifier(alias_A = false, alias_b = false)) + lincache = init( + linprob, linsolve; alias = LinearAliasSpecifier(alias_A = false, alias_b = false)) return LinearSolveJLCache(lincache, linsolve, nothing, stats) end diff --git a/lib/SimpleNonlinearSolve/src/utils.jl b/lib/SimpleNonlinearSolve/src/utils.jl index c19e7538c..1090ac5f1 100644 --- a/lib/SimpleNonlinearSolve/src/utils.jl +++ b/lib/SimpleNonlinearSolve/src/utils.jl @@ -160,7 +160,7 @@ end function compute_hvvp(prob, autodiff, _, x::Number, dir::Number) H = DI.second_derivative(prob.f, autodiff, x, Constant(prob.p)) - return H*dir + return H * dir end function compute_hvvp(prob, autodiff, fx, x, dir) jvp_fn = if SciMLBase.isinplace(prob) diff --git a/test/qa_tests.jl b/test/qa_tests.jl index db7dcc073..7dce47324 100644 --- a/test/qa_tests.jl +++ b/test/qa_tests.jl @@ -1,7 +1,8 @@ @testitem "Aqua" tags=[:misc] begin using NonlinearSolve, SimpleNonlinearSolve, Aqua - Aqua.test_all(NonlinearSolve; ambiguities = false, piracies = false, stale_deps = false, deps_compat = false) + Aqua.test_all(NonlinearSolve; ambiguities = false, piracies = false, + stale_deps = false, deps_compat = false) Aqua.test_ambiguities(NonlinearSolve; recursive = false) Aqua.test_stale_deps(SimpleNonlinearSolve; ignore = [:SciMLJacobianOperators]) Aqua.test_deps_compat(SimpleNonlinearSolve; ignore = [:SciMLJacobianOperators]) From 8d0d8d91c09ff3e83854afede76bc71840b60d5d Mon Sep 17 00:00:00 2001 From: Qingyu Qu <2283984853@qq.com> Date: Mon, 24 Mar 2025 02:14:47 +0800 Subject: [PATCH 3/4] Fix typo in docs --- docs/src/basics/diagnostics_api.md | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/src/basics/diagnostics_api.md b/docs/src/basics/diagnostics_api.md index 795348bd6..b3f8dea4c 100644 --- a/docs/src/basics/diagnostics_api.md +++ b/docs/src/basics/diagnostics_api.md @@ -20,7 +20,7 @@ solve process. This is controlled by 3 keyword arguments to `solve`: All the native NonlinearSolve.jl algorithms come with in-built [TimerOutputs.jl](https://github.com/KristofferC/TimerOutputs.jl) support. However, this -is disabled by default and can be enabled via [`NonlinearSolve.enable_timer_outputs`](@ref). +is disabled by default and can be enabled via [`NonlinearSolveBase.enable_timer_outputs`](@ref). Note that you will have to restart Julia to disable the timer outputs once enabled. From 4ad1c97e65fa6f4f4850cc683d0aed8ef97c6791 Mon Sep 17 00:00:00 2001 From: Qingyu Qu <2283984853@qq.com> Date: Mon, 24 Mar 2025 10:15:16 +0800 Subject: [PATCH 4/4] Restrict type of T --- lib/NonlinearSolveBase/src/utils.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/lib/NonlinearSolveBase/src/utils.jl b/lib/NonlinearSolveBase/src/utils.jl index 12db8d961..933dbc773 100644 --- a/lib/NonlinearSolveBase/src/utils.jl +++ b/lib/NonlinearSolveBase/src/utils.jl @@ -106,7 +106,7 @@ end init_similar_array!!(x) = x -function init_similar_array!!(x::AbstractArray{<:T}) where {T} +function init_similar_array!!(x::AbstractArray{<:T}) where {T <: Number} ArrayInterface.can_setindex(x) && fill!(x, T(0)) return x end