refactor: migrate Curvefit from Optim.jl to CurveFit.jl

Project.toml

@@ -12,18 +12,18 @@ RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 
 [weakdeps]
+CurveFit = "5a033b19-8c74-5913-a970-47c3779ef25c"
 ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"
 FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
 Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
-Optim = "429524aa-4258-5aef-a3af-852621145aeb"
 RegularizationTools = "29dad682-9a27-4bc3-9c72-016788665182"
 SparseConnectivityTracer = "9f842d2f-2579-4b1d-911e-f412cf18a3f5"
 Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f"
 
 [extensions]
 DataInterpolationsChainRulesCoreExt = "ChainRulesCore"
-DataInterpolationsOptimExt = "Optim"
+DataInterpolationsCurveFitExt = "CurveFit"
 DataInterpolationsRegularizationToolsExt = "RegularizationTools"
 DataInterpolationsSparseConnectivityTracerExt = ["SparseConnectivityTracer", "FillArrays"]
 DataInterpolationsSymbolicsExt = "Symbolics"
@@ -34,14 +34,15 @@ AllocCheck = "0.2"
 Aqua = "0.8"
 BenchmarkTools = "1"
 ChainRulesCore = "1.25"
+CurveFit = "1.5"
 EnumX = "1.0.4"
 FillArrays = "1.13"
 FindFirstFunctions = "1.3"
 FiniteDifferences = "0.12.31"
 ForwardDiff = "0.10.37, 1"
 LinearAlgebra = "1.10"
 Makie = "0.22, 0.23, 0.24"
-Optim = "1.6, 2"
+NonlinearSolve = "4.16.0"
 PrettyTables = "2.4, 3"
 QuadGK = "2.9.1"
 RecipesBase = "1.3"
@@ -62,10 +63,11 @@ AllocCheck = "9b6a8646-10ed-4001-bbdc-1d2f46dfbb1a"
 Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
 BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
 ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4"
+CurveFit = "5a033b19-8c74-5913-a970-47c3779ef25c"
 FiniteDifferences = "26cc04aa-876d-5657-8c51-4c34ba976000"
 ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210"
 Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a"
-Optim = "429524aa-4258-5aef-a3af-852621145aeb"
+NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
 QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
 RegularizationTools = "29dad682-9a27-4bc3-9c72-016788665182"
 SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f"
@@ -78,4 +80,4 @@ Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f"
 
 [targets]
-test = ["Aqua", "AllocCheck", "BenchmarkTools", "SafeTestsets", "ChainRulesCore", "Optim", "RegularizationTools", "Test", "StableRNGs", "FiniteDifferences", "QuadGK", "ForwardDiff", "StaticArrays", "Symbolics", "Unitful", "Zygote", "SparseConnectivityTracer"]
+test = ["Aqua", "AllocCheck", "BenchmarkTools", "SafeTestsets", "ChainRulesCore", "CurveFit", "NonlinearSolve", "RegularizationTools", "Test", "StableRNGs", "FiniteDifferences", "QuadGK", "ForwardDiff", "StaticArrays", "Symbolics", "Unitful", "Zygote", "SparseConnectivityTracer"]

docs/Project.toml

@@ -1,19 +1,20 @@
 [deps]
+CurveFit = "5a033b19-8c74-5913-a970-47c3779ef25c"
+DataInterpolations = "82cc6244-b520-54b8-b5a6-8a565e85f1d0"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 ModelingToolkitStandardLibrary = "16a59e39-deab-5bd0-87e4-056b12336739"
-Optim = "429524aa-4258-5aef-a3af-852621145aeb"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 RegularizationTools = "29dad682-9a27-4bc3-9c72-016788665182"
 StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3"
 Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 
 [compat]
+CurveFit = "1.5"
 Documenter = "1"
 ModelingToolkit = "9, 10.1, 11.1"
 ModelingToolkitStandardLibrary = "2"
-Optim = "1, 2.0"
 OrdinaryDiffEq = "6"
 Plots = "1"
 RegularizationTools = "0.6"

docs/src/methods.md

@@ -263,8 +263,8 @@ match our data. Let's start with the guess of every `p` being zero, that is
 `p=ones(4)`. Then we would fit this curve using:
 
 ```@example tutorial
-using Optim
-A = Curvefit(u, t, m, ones(4), LBFGS())
+using CurveFit
+A = Curvefit(u, t, m, ones(4))
 plot(A)
 ```
 
@@ -280,7 +280,7 @@ parameters. For example, with `p=zeros(4)` as the initial parameters, the fit
 is not good:
 
 ```@example tutorial
-A = Curvefit(u, t, m, zeros(4), LBFGS())
+A = Curvefit(u, t, m, zeros(4))
 plot(A)
 ```
 

ext/DataInterpolationsOptimExt.jl → ext/DataInterpolationsCurveFitExt.jl

@@ -1,37 +1,36 @@
-module DataInterpolationsOptimExt
+module DataInterpolationsCurveFitExt
 
 using DataInterpolations
 import DataInterpolations: munge_data,
     Curvefit, CurvefitCache, _interpolate, get_show, derivative,
     ExtrapolationError,
     integral, IntegralNotFoundError, DerivativeNotFoundError
 
-using Optim, ForwardDiff
+using CurveFit, ForwardDiff
 
 ### Curvefit
 function Curvefit(
         u,
         t,
         model,
         p0,
-        alg,
+        alg = nothing,
         box = false,
         lb = nothing,
         ub = nothing;
         extrapolate = false
     )
     u, t = munge_data(u, t)
-    errfun(t, u, p) = sum(abs2.(u .- model(t, p)))
-    if box == false
-        mfit = optimize(p -> errfun(t, u, p), p0, alg)
+    if box && (lb === nothing || ub === nothing)
+        error("lower or upper bound should not be nothing")
+    end
+    prob = CurveFit.NonlinearCurveFitProblem((p, x) -> model(x, p), p0, t, u; lb = lb, ub = ub)
+    sol = if isnothing(alg)
+        CurveFit.solve(prob)
     else
-        if lb === nothing || ub === nothing
-            error("lower or upper bound should not be nothing")
-        end
-        od = OnceDifferentiable(p -> errfun(t, u, p), p0, autodiff = :finite)
-        mfit = optimize(od, lb, ub, p0, Fminbox(alg))
+        CurveFit.solve(prob, alg)
     end
-    pmin = Optim.minimizer(mfit)
+    pmin = sol.u
     return CurvefitCache(u, t, model, p0, ub, lb, alg, pmin, extrapolate)
 end
 
@@ -65,7 +64,7 @@ end
 
 function get_show(A::CurvefitCache)
     return "Curvefit" *
-        " with $(length(A.t)) points, using $(nameof(typeof(A.alg)))\n"
+        " with $(length(A.t)) points.\n"
 end
 
 function integral(A::CurvefitCache{<:AbstractVector{<:Number}}, t::Number)

src/DataInterpolations.jl

@@ -215,9 +215,9 @@ struct CurvefitCache{
     end
 end
 
-# Define an empty function, so that it can be extended via `DataInterpolationsOptimExt`
+# Define an empty function, so that it can be extended via `DataInterpolationsCurveFitExt`
 function Curvefit()
-    error("CurveFit requires loading Optim and ForwardDiff, e.g. `using Optim, ForwardDiff`")
+    error("CurveFit requires loading CurveFit, e.g. `using CurveFit`")
 end
 
 export Curvefit

test/derivative_tests.jl

@@ -5,7 +5,7 @@ using DataInterpolations: derivative, get_transition_ts
 using Symbolics
 using StableRNGs
 using RegularizationTools
-using Optim
+using CurveFit
 import ForwardDiff
 using LinearAlgebra
 
@@ -397,7 +397,7 @@ end
     t = range(-10, stop = 10, length = 40)
     u = model(t, [1.0, 2.0]) + 0.01 * randn(rng, length(t))
     p0 = [0.5, 0.5]
-    test_derivatives(Curvefit; args = [u, t, model, p0, LBFGS()], name = "Curvefit")
+    test_derivatives(Curvefit; args = [u, t, model, p0], name = "Curvefit")
 end
 
 @testset "Symbolic derivatives" begin

test/integral_tests.jl

@@ -1,7 +1,7 @@
 using DataInterpolations, Test
 using QuadGK
 using DataInterpolations: integral
-using Optim, ForwardDiff
+using CurveFit, ForwardDiff
 using RegularizationTools
 using StableRNGs
 using Unitful
@@ -247,7 +247,7 @@ end
     t = range(-10, stop = 10, length = 40)
     u = model(t, [1.0, 2.0]) + 0.01 * randn(rng, length(t))
     p0 = [0.5, 0.5]
-    A = Curvefit(u, t, model, p0, LBFGS())
+    A = Curvefit(u, t, model, p0)
     @test_throws DataInterpolations.IntegralNotFoundError integral(A, 0.0, 1.0)
     @test_throws DataInterpolations.IntegralNotFoundError integral(A, 5.0)
 end

test/interpolation_tests.jl

@@ -1,7 +1,7 @@
 using DataInterpolations
 using FindFirstFunctions: searchsortedfirstcorrelated
 using StableRNGs
-using Optim, ForwardDiff
+using CurveFit, NonlinearSolve, ForwardDiff
 using BenchmarkTools
 using Unitful
 using LinearAlgebra
@@ -1225,7 +1225,7 @@ end
     u = model(t, [1.0, 2.0]) + 0.01 * randn(rng, length(t))
     p0 = [0.5, 0.5]
 
-    A = Curvefit(u, t, model, p0, LBFGS())
+    A = Curvefit(u, t, model, p0, LevenbergMarquardt())
 
     ts = [-7.0, -2.0, 0.0, 2.5, 5.0]
     vs = [
@@ -1241,10 +1241,14 @@ end
     @test @inferred(output_size(A)) == ()
 
     # Test extrapolation
-    A = Curvefit(u, t, model, p0, LBFGS(); extrapolate = true)
+    A = Curvefit(u, t, model, p0, LevenbergMarquardt(); extrapolate = true)
     @test A(15.0) == model(15.0, A.pmin)
-    A = Curvefit(u, t, model, p0, LBFGS())
+    A = Curvefit(u, t, model, p0, LevenbergMarquardt())
     @test_throws DataInterpolations.ExtrapolationError A(15.0)
+
+    # With lb, ub
+    A = Curvefit(u, t, model, p0, LevenbergMarquardt(), false, [0.0, 0.0], [1.0, 1.0])
+    @test all(0.0 .<= A.pmin .<= 1.0)
 end
 
 @testset "Type of vector returned" begin

test/show.jl

@@ -1,5 +1,5 @@
 using DataInterpolations
-using Optim, StableRNGs
+using CurveFit, StableRNGs
 using RegularizationTools
 
 t = [1.0, 2.0, 3.0, 4.0, 5.0]
@@ -67,10 +67,10 @@ end
     t = range(-10, stop = 10, length = 40)
     u = model(t, [1.0, 2.0]) + 0.01 * randn(rng, length(t))
     p0 = [0.5, 0.5]
-    A = Curvefit(u, t, model, p0, LBFGS())
+    A = Curvefit(u, t, model, p0)
     @test startswith(
         sprint(io -> show(io, MIME"text/plain"(), A)),
-        "Curvefit with 40 points, using LBFGS\n"
+        "Curvefit with 40 points.\n"
     )
 end