Merge pull request #301 from SciML/mtk_update

Update to MTK9

Project.toml

@@ -40,7 +40,7 @@ IterTools = "1"
 LinearAlgebra = "1.6, 1.7"
 Logging = "1.6, 1.7"
 MacroTools = "0.5"
-ModelingToolkit = "8.75"
+ModelingToolkit = "9"
 Nemo = "0.43"
 ParamPunPam = "0.4"
 Pkg = "1.6, 1.7"
@@ -49,7 +49,7 @@ Primes = "0.5"
 Random = "1.6, 1.7"
 SpecialFunctions = "2"
 SymbolicUtils = "1.4, 1.5"
-Symbolics = "5.17"
+Symbolics = "5.20"
 Test = "1.6, 1.7"
 TestSetExtensions = "2"
 TimerOutputs = "0.5"

docs/src/tutorials/discrete_time.md

@@ -87,39 +87,6 @@ assess_local_identifiability(dds; funcs_to_check = [β * S])
 As other main functions in the package, `assess_local_identifiability` accepts an optional parameter `loglevel` (default: `Logging.Info`)
 to adjust the verbosity of logging.
 
-If one loads `ModelingToolkit` (and thus the `ModelingToolkitSIExt` extension), one can use `DiscreteSystem` from `ModelingToolkit` to
-describe the input model (now in terms of difference!):
-
-```@example discrete_mtk
-using ModelingToolkit
-using StructuralIdentifiability
-
-@parameters α β
-@variables t S(t) I(t) R(t) y(t)
-D = Difference(t; dt = 1.0)
-
-eqs = [D(S) ~ S - β * S * I, D(I) ~ I + β * S * I - α * I, D(R) ~ R + α * I]
-@named sir = DiscreteSystem(eqs)
-```
-
-Then the same computation can be carried out with the models defined this way:
-
-```@example discrete_mtk
-assess_local_identifiability(sir; measured_quantities = [y ~ I])
-```
-
-In principle, it is not required to give a name to the observable, so one can write this shorter
-
-```@example discrete_mtk
-assess_local_identifiability(sir; measured_quantities = [I])
-```
-
-The same example but with specified functions to check
-
-```@example discrete_mtk
-assess_local_identifiability(sir; measured_quantities = [I], funcs_to_check = [β * S])
-```
-
 The implementation is based on a version of the observability rank criterion and will be described in a forthcoming paper.
 
 [^1]: > S. Nõmm, C. Moog, [*Identifiability of discrete-time nonlinear systems*](https://doi.org/10.1016/S1474-6670(17)31245-4), IFAC Proceedings Volumes, 2004.

ext/ModelingToolkitSIExt.jl

@@ -180,10 +180,10 @@ function __mtk_to_si(
     end
 
     y_functions = [each[2] for each in measured_quantities]
-    inputs = filter(v -> ModelingToolkit.isinput(v), ModelingToolkit.states(de))
+    inputs = filter(v -> ModelingToolkit.isinput(v), ModelingToolkit.unknowns(de))
     state_vars = filter(
         s -> !(ModelingToolkit.isinput(s) || ModelingToolkit.isoutput(s)),
-        ModelingToolkit.states(de),
+        ModelingToolkit.unknowns(de),
     )
     params = ModelingToolkit.parameters(de)
     t = ModelingToolkit.arguments(diff_eqs[1].lhs)[1]
@@ -299,7 +299,7 @@ end
     end
     if length(funcs_to_check) == 0
         funcs_to_check = vcat(
-            [e for e in ModelingToolkit.states(ode) if !ModelingToolkit.isoutput(e)],
+            [e for e in ModelingToolkit.unknowns(ode) if !ModelingToolkit.isoutput(e)],
             ModelingToolkit.parameters(ode),
         )
     end
@@ -413,6 +413,7 @@ Output:
 
 The result is correct with probability at least `prob_threshold`.
 """
+#=
 function StructuralIdentifiability.assess_local_identifiability(
     dds::ModelingToolkit.DiscreteSystem;
     measured_quantities = Array{ModelingToolkit.Equation}[],
@@ -498,6 +499,8 @@ function _assess_local_identifiability(
     end
     return out_dict
 end
+=#
+
 # ------------------------------------------------------------------------------
 
 """

src/global_identifiability.jl

@@ -106,8 +106,8 @@ The function returns a tuple containing the following:
     rational_interpolator = :VanDerHoevenLecerf,
 ) where {T}
     @info "Computing IO-equations"
-    ioeq_time =
-        @elapsed io_equations = find_ioequations(ode; var_change_policy = var_change_policy)
+    ioeq_time = @elapsed io_equations =
+        _find_ioequations(ode; var_change_policy = var_change_policy)
     @debug "Sizes: $(map(length, values(io_equations)))"
     @info "Computed IO-equations in $ioeq_time seconds"
     _runtime_logger[:ioeq_time] = ioeq_time

src/io_equation.jl

@@ -346,6 +346,17 @@ Output:
     ode::ODE{P};
     var_change_policy = :default,
     loglevel = Logging.Info,
+) where {P <: MPolyRingElem{<:FieldElem}}
+    restart_logging(loglevel = loglevel)
+    reset_timings()
+    with_logger(_si_logger[]) do
+        return _find_ioequations(ode, var_change_policy = var_change_policy)
+    end
+end
+
+@timeit _to function _find_ioequations(
+    ode::ODE{P};
+    var_change_policy = :default,
 ) where {P <: MPolyRingElem{<:FieldElem}}
     # Setting the var_change policy
     if (var_change_policy == :yes) ||

test/extensions/modelingtoolkit.jl

@@ -358,7 +358,7 @@ if GROUP == "All" || GROUP == "ModelingToolkitSIExt"
         function getbyname(sys, name)
             println(name)
             return first([
-                v for v in vcat(states(sys), parameters(sys)) if
+                v for v in vcat(unknowns(sys), parameters(sys)) if
                 replace(string(v), "(t)" => "") == name
             ])
         end
@@ -412,6 +412,7 @@ if GROUP == "All" || GROUP == "ModelingToolkitSIExt"
               correct
     end
 
+    #=
     @testset "Discrete local identifiability, ModelingToolkit interface" begin
         cases = []
 
@@ -656,7 +657,7 @@ if GROUP == "All" || GROUP == "ModelingToolkitSIExt"
             ) == c[:res]
         end
     end
-
+    =#
     @testset "Exporting ModelingToolkit Model to SI Model" begin
 
         # Creates MTK model and assesses its identifiability.