Move to MTK v11

CHANGELOG.md

@@ -1,3 +1,7 @@
+# v1.4
+
+- Moved to MTKv11 (ProcessBasedModelling.jl v1.9)
+
 # v1.3
 
 ## `CloudToppedMixedLayerModel`

Project.toml

@@ -1,19 +1,21 @@
 name = "ConceptualClimateModels"
 uuid = "3e640dcb-3446-4ed7-aadb-0366b77312ec"
 authors = ["Datseris <datseris.george@gmail.com>"]
-version = "1.3.0"
+version = "1.4.0"
 
 [deps]
 DynamicalSystemsBase = "6e36e845-645a-534a-86f2-f5d4aa5a06b4"
 NaNMath = "77ba4419-2d1f-58cd-9bb1-8ffee604a2e3"
+OrdinaryDiffEqNonlinearSolve = "127b3ac7-2247-4354-8eb6-78cf4e7c58e8"
 ProcessBasedModelling = "ca969041-2cf3-4b10-bc21-86f4417093eb"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665"
 
 [compat]
 DynamicalSystemsBase = "3.5"
 NaNMath = "1"
-ProcessBasedModelling = "1.8"
+OrdinaryDiffEqNonlinearSolve = "2"
+ProcessBasedModelling = "1.9"
 Reexport = "1.2"
 Roots = "2"
 julia = "1.9"

docs/make.jl

@@ -34,6 +34,6 @@ bib = CitationBibliography(
 
 build_docs_with_style(pages, ConceptualClimateModels, GlobalMeanEBM, CloudToppedMixedLayerModel, ProcessBasedModelling;
     authors = "George Datseris <datseris.george@gmail.com>",
-    bib, warnonly = [:doctest, :missing_docs, :cross_references],
+    bib, warnonly = [:doctest, :missing_docs, :cross_references, :linkcheck],
     repo = Remotes.GitHub("JuliaDynamics", "ConceptualClimateModels.jl")
 )

src/dynamicalsystems.jl

@@ -7,33 +7,42 @@ export all_equations
 export named_current_parameters
 export try_set_parameter!
 
+# This is needed because otherwise we get the error:
+# This ODE requires a DAE initialization and thus a nonlinear solve but no nonlinear solve
+# has been loaded. To solve this problem, do `using OrdinaryDiffEqNonlinearSolve` or pass
+# a custom `nlsolve` choice into the `initializealg`.
+using OrdinaryDiffEqNonlinearSolve
+
 DEFAULT_DIFFEQ = DynamicalSystemsBase.DEFAULT_DIFFEQ
 
 """
     processes_to_coupledodes(processes [, default]; kw...)
 
-Convert a given `Vector` of processes to a `DynamicalSystem`, in particular `CoupledODEs`.
+Convert a given `Vector` of processes to `CoupledODEs`.
 All processes represent symbolic equations managed by ModelingToolkit.jl.
 `default` is a vector for default processes that "process-less" variables
 introduced in `processes` will obtain.
-Use [`processes_to_mtkmodel`](@ref) to obtain the MTK model before it is structurally
-simplified and converted to a `DynamicalSystem`.
+Use [`processes_to_mtkmodel`](@ref) to obtain the MTK model before it is compiled
+and converted to `CoupledODEs`.
 See also [`processes_to_mtkmodel`](@ref) for more details on what `processes` is,
 or see the online [Tutorial](@ref).
 
 ## Keyword arguments
 
-- `diffeq`: options passed to DifferentialEquations.jl ODE solving
+- `diffeq`: options passed to OrdinaryDiffEq.jl ODE solving
   when constructing the `CoupledODEs`.
 - `inplace`: whether the dynamical system will be in place or not.
   Defaults to `true` if the system dimension is ≤ 5.
 - `split = false`: whether to split parameters as per ModelingToolkit.jl.
   Note the default is not ModelingToolkit's default, i.e., no splitting occurs.
-  This accelerates parameter access, assuming all parameters are of the same type.
+  This accelerates parameter access by assuming all parameters are of the same type.
+- `probkw = (warn_initialize_determined = false, )`: keyword arguments expanded into
+  the creation of `ODEProblem`.
 - `kw...`: all other keywords are propagated to `processes_to_mtkmodel`.
 """
 function processes_to_coupledodes(proc, default = [];
-        diffeq = DEFAULT_DIFFEQ, inplace::Bool = false, split::Bool = false, kwargs...
+        diffeq = DEFAULT_DIFFEQ, inplace = nothing, split::Bool = false,
+        probkw = (warn_initialize_determined = false, ), kwargs...
     )
     sys = processes_to_mtkmodel(proc, default; kwargs...)
     ssys = mtkcompile(sys; split)
@@ -45,9 +54,9 @@ function processes_to_coupledodes(proc, default = [];
     # The usage of `nothing` for the initial state assumes all state variables
     # and all parameters have been defined with a default value. This also means
     if inplace
-        prob = ODEProblem(ssys, nothing, (0.0, Inf))
+        prob = ODEProblem(ssys, nothing, (0.0, Inf); probkw...)
     else
-        prob = ODEProblem{false}(ssys, nothing, (0.0, Inf); u0_constructor = x->SVector(x...))
+        prob = ODEProblem{false}(ssys, nothing, (0.0, Inf); u0_constructor = x->SVector(x...), probkw...)
     end
     ds = CoupledODEs(prob, diffeq)
     return ds

test/Project.toml

@@ -1,5 +1,5 @@
 [deps]
-DynamicalSystems = "61744808-ddfa-5f27-97ff-6e42cc95d634"
+Attractors = "f3fd9213-ca85-4dba-9dfd-7fc91308fec7"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"

test/runtests.jl

@@ -3,88 +3,85 @@ using ConceptualClimateModels.GlobalMeanEBM
 using Test
 
 @testset "all processes" begin
-function test_symbolic_var(mtk, var)
-    if has_symbolic_var(mtk, var)
-        @test true
-    else
-        error("var $var doesn't exist")
+    function test_symbolic_var(mtk, var)
+        if has_symbolic_var(mtk, var)
+            @test true
+        else
+            error("var $var doesn't exist")
+        end
     end
-end
-
-@variables begin
-    ε1(t) = 0.1
-    ε2(t) = 0.1
-    ε3(t) = 0.1
-    ε4(t) = 0.1
-end
 
-p1 = [
-    BasicRadiationBalance(),
-    EmissivityStefanBoltzmanOLR(),
-    IceAlbedoFeedback(; min = 0.3, max = 0.7),
-    α ~ α_ice,
-    ParameterProcess(ε), # emissivity is a parameter
-    f ~ 0, # no external forcing
-    # absorbed solar radiation has a default process
-    EmissivityFeedbackTanh(ε = ε1),
-    SoedergrenClearSkyEmissivity(ε = ε2),
-    EmissivitySellers1969(ε = ε3),
-    EmissivityLogConcentration(ε = ε4),
-]
-
-ds = processes_to_coupledodes(p1, GlobalMeanEBM)
-mtk = referrenced_sciml_model(ds)
-eqs = all_equations(mtk)
-
-# This also tests a bunch of default processes
-for var in (T, ε, ε1, ε2, ε3, ε4, q, :ε1_T_sigmoid_ref, mtk.ε_0)
-    test_symbolic_var(eqs, var)
-end
+    @variables begin
+        ε1(t) = 0.1
+        ε2(t) = 0.1
+        ε3(t) = 0.1
+        ε4(t) = 0.1
+    end
 
+    p1 = [
+        BasicRadiationBalance(),
+        EmissivityStefanBoltzmanOLR(),
+        IceAlbedoFeedback(; min = 0.3, max = 0.7),
+        α ~ α_ice,
+        ParameterProcess(ε), # emissivity is a parameter
+        f ~ 0, # no external forcing
+        # absorbed solar radiation has a default process
+        EmissivityFeedbackTanh(ε = ε1),
+        SoedergrenClearSkyEmissivity(ε = ε2),
+        EmissivitySellers1969(ε = ε3),
+        EmissivityLogConcentration(ε = ε4),
+    ]
 
-@variables begin
-    OLR1(t) = 10.0
-    OLR2(t) = 10.0
-    OLR3(t) = 10.0
-    a1(t) = 0.1
-    a2(t) = 0.1
-    a3(t) = 0.1
-    a4(t) = 0.1
-    ΔS1(t) = 0.0
-    ΔT1(t) = 0.0
-end
+    ds = processes_to_coupledodes(p1, GlobalMeanEBM)
+    mtk = referrenced_sciml_model(ds)
+    eqs = all_equations(mtk)
 
-p2 = [
-    LinearOLR(),
-    LinearClearSkyOLR(; OLR = OLR1),
-    BudykoOLR(; OLR = OLR2),
-    CloudAlbedoExponential(),
-    CloudAlbedoLinear(; α_cloud = a1),
-    DirectAlbedoAddition(α = a4),
-    CoAlbedoProduct(α = a2),
-    SeparatedClearAllSkyAlbedo(α = a3),
-    C ~ 0.6,
-    ΔTLinearRelaxation(),
-    ΔTStommelModel(ΔT = ΔT1, ΔS = ΔS1),
-]
-
-ds = processes_to_coupledodes(p2, GlobalMeanEBM)
-mtk = referrenced_sciml_model(ds)
-for var in (T, C, OLR1, OLR2, :α_bg, a1, a2, a3, a4, ΔS1, ΔT)
-    if has_symbolic_var(mtk, var)
-        @test true
-    else
-        error("var $var doesn't exist")
+    # This also tests a bunch of default processes
+    for var in (T, ε, ε1, ε2, ε3, ε4, q, :ε1_T_sigmoid_ref, mtk.ε_0)
+        test_symbolic_var(eqs, var)
     end
-end
 
 
+    @variables begin
+        OLR1(t) = 10.0
+        OLR2(t) = 10.0
+        OLR3(t) = 10.0
+        a1(t) = 0.1
+        a2(t) = 0.1
+        a3(t) = 0.1
+        a4(t) = 0.1
+        ΔS1(t) = 0.0
+        ΔT1(t) = 0.0
+    end
+
+    p2 = [
+        LinearOLR(),
+        LinearClearSkyOLR(; OLR = OLR1),
+        BudykoOLR(; OLR = OLR2),
+        CloudAlbedoExponential(),
+        CloudAlbedoLinear(; α_cloud = a1),
+        DirectAlbedoAddition(α = a4),
+        CoAlbedoProduct(α = a2),
+        SeparatedClearAllSkyAlbedo(α = a3),
+        C ~ 0.6,
+        ΔTLinearRelaxation(),
+        ΔTStommelModel(ΔT = ΔT1, ΔS = ΔS1),
+    ]
 
+    ds = processes_to_coupledodes(p2, GlobalMeanEBM)
+    mtk = referrenced_sciml_model(ds)
+    for var in (T, C, OLR1, OLR2, :α_bg, a1, a2, a3, a4, ΔS1, ΔT)
+        if has_symbolic_var(mtk, var)
+            @test true
+        else
+            error("var $var doesn't exist")
+        end
+    end
 end
 
 
 @testset "dynamical systems integration" begin
-    using DynamicalSystems
+    using Attractors
 
     budyko_processes = [
         BasicRadiationBalance(),
@@ -100,7 +97,6 @@ end
 
     grid = plausible_grid(budyko)
     mapper = AttractorsViaRecurrences(budyko, grid)
-    rfam = RecurrencesFindAndMatch(mapper)
     sampler = plausible_ic_sampler(budyko)
 
     set_parameter!(budyko, :ε_0, 0.5)
@@ -109,5 +105,4 @@ end
     attractors = extract_attractors(mapper)
 
     @test length(attractors) == 2
-
 end