Merge pull request #57 from EarthSciML/mtk9

Update to MTK v9
EarthSciML · Aug 21, 2024 · 21d7edb · 21d7edb · ctessum · Aug 21, 2024
2 parents 16f48de + a71359b
commit 21d7edb
Show file tree

Hide file tree

Showing 26 changed files with 222 additions and 264 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,14 +1,15 @@
 name = "EarthSciMLBase"
 uuid = "e53f1632-a13c-4728-9402-0c66d48804b0"
 authors = ["EarthSciML Authors and Contributors"]
-version = "0.15.3"
+version = "0.16.0"
 
 [deps]
 BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0"
 Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83"
 DiffEqCallbacks = "459566f4-90b8-5000-8ac3-15dfb0a30def"
 DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
 DomainSets = "5b8099bc-c8ec-5219-889f-1d9e522a28bf"
+DynamicQuantities = "06fc5a27-2a28-4c7c-a15d-362465fb6821"
 Graphs = "86223c79-3864-5bf0-83f7-82e725a168b6"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 MetaGraphsNext = "fa8bd995-216d-47f1-8a91-f3b68fbeb377"
@@ -18,25 +19,24 @@ ProgressLogging = "33c8b6b6-d38a-422a-b730-caa89a2f386c"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SciMLOperators = "c0aeaf25-5076-4817-a8d5-81caf7dfa961"
 Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
-Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [compat]
 BlockBandedMatrices = "0.13"
-Catalyst = "13"
+Catalyst = "14"
 DiffEqCallbacks = "2"
 DifferentialEquations = "7"
 DocStringExtensions = "0.9"
-DomainSets = "0.5, 0.6"
+DomainSets = "0.7"
+DynamicQuantities = "0.13"
 Graphs = "1"
 MetaGraphsNext = "0.5, 0.6, 0.7"
-MethodOfLines = "0.10"
-ModelingToolkit = "8"
+MethodOfLines = "0.11"
+ModelingToolkit = "9"
 OrdinaryDiffEq = "6"
 ProgressLogging = "0.1"
 SciMLBase = "2"
 SciMLOperators = "0.3"
 Symbolics = "5"
-Unitful = "1"
 julia = "1.6"
 
 [extras]

diff --git a/docs/Project.toml b/docs/Project.toml
@@ -4,6 +4,7 @@ Catalyst = "479239e8-5488-4da2-87a7-35f2df7eef83"
 DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 DomainSets = "5b8099bc-c8ec-5219-889f-1d9e522a28bf"
+DynamicQuantities = "06fc5a27-2a28-4c7c-a15d-362465fb6821"
 EarthSciMLBase = "e53f1632-a13c-4728-9402-0c66d48804b0"
 GraphMakie = "1ecd5474-83a3-4783-bb4f-06765db800d2"
 MetaGraphsNext = "fa8bd995-216d-47f1-8a91-f3b68fbeb377"
@@ -12,19 +13,18 @@ ModelingToolkit = "961ee093-0014-501f-94e3-6117800e7a78"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 SciMLOperators = "c0aeaf25-5076-4817-a8d5-81caf7dfa961"
 Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
-Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [compat]
 CairoMakie = "0.12"
-Catalyst = "13"
+Catalyst = "14"
 DifferentialEquations = "7"
 Documenter = "1"
 DomainSets = "0.6,0.7"
+DynamicQuantities = "0.13"
 GraphMakie = "0.5"
 MetaGraphsNext = "0.7"
-MethodOfLines = "0.8,0.9,0.10"
-ModelingToolkit = "8"
+MethodOfLines = "0.11"
+ModelingToolkit = "9"
 Plots = "1"
 SciMLOperators = "0.3"
 Symbolics = "5"
-Unitful = "1"
diff --git a/docs/src/advection.md b/docs/src/advection.md
@@ -9,18 +9,18 @@ Advection is implemented with the [`Advection`](@ref) type.
 To demonstrate how this can work, we will start with a simple system of equations:
 
 ```@example advection
-using EarthSciMLBase, ModelingToolkit, Unitful
+using EarthSciMLBase, ModelingToolkit
+using ModelingToolkit: t_nounits, D_nounits
+t = t_nounits
+D = D_nounits
 
-@parameters t
-
-function ExampleSys(t)
+function ExampleSys()
     @variables y(t)
     @parameters p=2.0
-    D = Differential(t)
     ODESystem([D(y) ~ p], t; name=:ExampleSys)
 end
 
-ExampleSys(t)
+ExampleSys()
 ```
 
 We also need to create our initial and boundary conditions.
@@ -35,8 +35,8 @@ Now we convert add advection to each of the state variables.
 We're also adding a constant wind ([`ConstantWind`](@ref)) in the x-direction, with a speed of 1.0.
 
 ```@example advection
-sys_advection = couple(ExampleSys(t), domain, ConstantWind(t, 1.0), Advection())
-sys_mtk = get_mtk(sys_advection)
+sys_advection = couple(ExampleSys(), domain, ConstantWind(t, 1.0), Advection())
+sys_mtk = convert(PDESystem, sys_advection)
 ```
 
 Finally, we can discretize the system and solve it:

diff --git a/docs/src/composition.md b/docs/src/composition.md
@@ -7,10 +7,11 @@ To demonstrate how this works, below we define three model components, `Photolys
 
 ```@example composition
 using ModelingToolkit, Catalyst, EarthSciMLBase
-@parameters t
+using ModelingToolkit: t_nounits
+t = t_nounits
 
 struct PhotolysisCoupler sys end
-function Photolysis(t; name=:Photolysis)
+function Photolysis(; name=:Photolysis)
     @variables j_NO2(t)
     eqs = [
         j_NO2 ~ max(sin(t/86400),0)
@@ -19,7 +20,7 @@ function Photolysis(t; name=:Photolysis)
         metadata=Dict(:coupletype=>PhotolysisCoupler))
 end
 
-Photolysis(t)
+Photolysis()
 ```
 
 You can see that the system defined above is mostly a standard ModelingToolkit ODE system,
@@ -44,34 +45,34 @@ Let's follow the same process for some additional components:
 
 ```@example composition
 struct ChemistryCoupler sys end
-function Chemistry(t; name=:Chemistry)
+function Chemistry(; name=:Chemistry)
     @parameters jNO2
     @species NO2(t)
     rxs = [
-        Reaction(jNO2, [NO2], [], [1], [1])
+        Reaction(jNO2, [NO2], [], [1], [])
     ]
     rsys = ReactionSystem(rxs, t, [NO2], [jNO2]; 
         combinatoric_ratelaws=false, name=name)
-    convert(ODESystem, rsys, metadata=Dict(:coupletype=>ChemistryCoupler))
+    convert(ODESystem, complete(rsys), metadata=Dict(:coupletype=>ChemistryCoupler))
 end
 
-Chemistry(t)
+Chemistry()
 ```
 
 For our chemistry component above, because it's is originally a ReactionSystem instead of an
 ODESystem, we convert it to an ODESystem before adding the metadata.
 
 ```@example composition
 struct EmissionsCoupler sys end
-function Emissions(t; name=:Emissions)
+function Emissions(; name=:Emissions)
     @parameters emis = 1
     @variables NO2(t)
     eqs = [NO2 ~ emis]
     ODESystem(eqs, t; name=name,
         metadata=Dict(:coupletype=>EmissionsCoupler))
 end
 
-Emissions(t)
+Emissions()
 ```
 
 Now, we need to define ways to couple the model components together.
@@ -104,11 +105,11 @@ end
 nothing # hide
 ```
 
-Finally, we can compose the model components together to create our complete model. To do so, we just initialize each model component and add the components together using the [`couple`](@ref) function. We can use the [`get_mtk`](@ref) function to convert the composed model to a [ModelingToolkit](https://mtk.sciml.ai/dev/) model so we can see what the combined equations look like.
+Finally, we can compose the model components together to create our complete model. To do so, we just initialize each model component and add the components together using the [`couple`](@ref) function. We can use the `convert` function to convert the composed model to a [ModelingToolkit](https://mtk.sciml.ai/dev/) `ODESystem` so we can see what the combined equations look like.
 
 ```@example composition
-model = couple(Photolysis(t), Chemistry(t), Emissions(t))
-get_mtk(model)
+model = couple(Photolysis(), Chemistry(), Emissions())
+convert(ODESystem, model)
 ```
 
 Finally, we can use the [`graph`](@ref) function to visualize the model components and their connections.

diff --git a/docs/src/coord_transforms.md b/docs/src/coord_transforms.md
@@ -8,12 +8,12 @@ To handle this, the [`DomainInfo`](@ref) type can be used to define coordinate s
 ```@example trans
 using EarthSciMLBase
 using ModelingToolkit
+using ModelingToolkit: t, D
 using DomainSets
-using Unitful
+using DynamicQuantities
 
 @parameters lon [unit = u"rad"]
 @parameters lat [unit = u"rad"]
-@parameters t [unit = u"s"]
 @parameters lev
 
 partialderivatives_δxyδlonlat([lev, lon, lat])
@@ -57,7 +57,7 @@ This returns a list of functions, one corresponding to each coordinate in our do
 Then we can calculate the symbolic partial derivative of a variable by just calling each function:
 
 ```@example trans
-@variables u
+@variables u(t)
 
 [δs[i](u) for i ∈ eachindex(δs)]
 ```

diff --git a/docs/src/example_all_together.md b/docs/src/example_all_together.md
@@ -11,37 +11,39 @@ Our first example system is a simple reaction system:
 ```@example ex1
 using EarthSciMLBase
 using ModelingToolkit, Catalyst, DomainSets, MethodOfLines, DifferentialEquations
+using ModelingToolkit: t_nounits, D_nounits
+t = t_nounits
+D = D_nounits
 using Plots
 
 # Create our independent variable `t` and our partially-independent variables `x` and `y`.
-@parameters t x y
+@parameters x y
 
 struct ExampleSys1Coupler sys end
-function ExampleSys1(t)
+function ExampleSys1()
     @species c₁(t)=5.0 c₂(t)=5.0
     rs = ReactionSystem(
         [Reaction(2.0, [c₁], [c₂])],
         t; name=:Sys1, combinatoric_ratelaws=false)
-    convert(ODESystem, rs, metadata=Dict(:coupletype=>ExampleSys1Coupler))
+    convert(ODESystem, complete(rs), metadata=Dict(:coupletype=>ExampleSys1Coupler))
 end
 
-ExampleSys1(t)
+ExampleSys1()
 ```
 
 Our second example system is a simple ODE system, with the same two variables.
 
 ```@example ex1
 struct ExampleSys2Coupler sys end
-function ExampleSys2(t)
+function ExampleSys2()
     @variables c₁(t)=5.0 c₂(t)=5.0
     @parameters p₁=1.0 p₂=0.5
-    D = Differential(t)
     ODESystem(
         [D(c₁) ~ -p₁, D(c₂) ~ p₂],
         t; name=:Sys2, metadata=Dict(:coupletype=>ExampleSys2Coupler))
 end
 
-ExampleSys2(t)
+ExampleSys2()
 ```
 
 Now, we specify what should happen when we couple the two systems together.
@@ -62,18 +64,18 @@ nothing # hide
 Once we specify all of the above, it is simple to create our two individual systems and then couple them together. 
 
 ```@example ex1
-sys1 = ExampleSys1(t)
-sys2 = ExampleSys2(t)
+sys1 = ExampleSys1()
+sys2 = ExampleSys2()
 sys = couple(sys1, sys2)
 
-get_mtk(sys)
+convert(ODESystem, sys)
 ```
 
 At this point we have an ODE system that is composed of two other ODE systems.
 We can inspect its observed variables using the `observed` function.
 
 ```@example ex1
-simplified_sys = structural_simplify(get_mtk(sys))
+simplified_sys = structural_simplify(convert(ODESystem, sys))
 ```
 
 ```@example ex1
@@ -105,7 +107,7 @@ domain = DomainInfo(
 
 sys_pde = couple(sys, domain, ConstantWind(t, 1.0, 1.0), Advection())
 
-sys_pde_mtk = get_mtk(sys_pde)
+sys_pde_mtk = convert(PDESystem, sys_pde)
 ```
 
 Now we can inspect this new system that we've created:

diff --git a/docs/src/icbc.md b/docs/src/icbc.md
@@ -8,20 +8,22 @@ To demonstrate how to do this, we will use the following simple system of ordina
 ```@example icbc
 using EarthSciMLBase
 using ModelingToolkit
+using ModelingToolkit: t_nounits, D_nounits
+t = t_nounits
+D = D_nounits
 
-@parameters x y t
+@parameters x y
 
-function ExampleSys(t)
+function ExampleSys()
     @variables u(t) v(t)
-    Dt = Differential(t)
     eqs = [
-        Dt(u) ~ √abs(v),
-        Dt(v) ~ √abs(u),
+        D(u) ~ √abs(v),
+        D(v) ~ √abs(u),
     ]
     ODESystem(eqs, t; name=:Docs₊Example)
 end
 
-ExampleSys(t)
+ExampleSys()
 ```
 Next, we specify our initial and boundary conditions using the [`DomainInfo`](@ref) type.
 We initialize [`DomainInfo`](@ref) with sets of initial and boundary conditions.
@@ -49,9 +51,9 @@ It is also possible to use periodic boundary conditions with [`periodicBC`](@ref
 Finally, we combine our initial and boundary conditions with our system of equations using the [`couple`](@ref) function.
 
 ```@example icbc
-model = couple(ExampleSys(t), icbc)
+model = couple(ExampleSys(), icbc)
 
-eq_sys = get_mtk(model)
+eq_sys = convert(PDESystem, model)
 ```
 
 We can also look at the expanded boundary conditions of the resulting equation system: