[v10] refactor: change inputs/outputs handling in structural_simplify

.github/workflows/Documentation.yml

@@ -4,6 +4,7 @@ on:
   push:
     branches:
       - master
+      - v10
     tags: '*'
   pull_request:
 

.github/workflows/FormatCheck.yml

@@ -4,6 +4,7 @@ on:
   push:
     branches:
       - 'master'
+      - v10
     tags: '*'
   pull_request:
 

.github/workflows/Tests.yml

@@ -5,6 +5,7 @@ on:
     branches:
       - master
       - 'release-'
+      - v10
     paths-ignore:
       - 'docs/**'
   push:

Project.toml

@@ -1,7 +1,7 @@
 name = "ModelingToolkit"
 uuid = "961ee093-0014-501f-94e3-6117800e7a78"
 authors = ["Yingbo Ma <mayingbo5@gmail.com>", "Chris Rackauckas <accounts@chrisrackauckas.com> and contributors"]
-version = "9.75.0"
+version = "10.0.0"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
@@ -174,6 +174,7 @@ Logging = "56ddb016-857b-54e1-b83d-db4d58db5568"
 ModelingToolkitStandardLibrary = "16a59e39-deab-5bd0-87e4-056b12336739"
 NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
 Optimization = "7f7a1694-90dd-40f0-9382-eb1efda571ba"
+OptimizationBase = "bca83a33-5cc9-4baa-983d-23429ab6bcbb"
 OptimizationMOI = "fd9f6733-72f4-499f-8506-86b2bdd0dea1"
 OptimizationOptimJL = "36348300-93cb-4f02-beb5-3c3902f8871e"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
@@ -193,5 +194,10 @@ StochasticDiffEq = "789caeaf-c7a9-5a7d-9973-96adeb23e2a0"
 Sundials = "c3572dad-4567-51f8-b174-8c6c989267f4"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
+[sources]
+ModelingToolkitStandardLibrary = { url = "https://github.com/SciML/ModelingToolkitStandardLibrary.jl/", rev = "mtk-v10" }
+OptimizationBase = { url = "https://github.com/AayushSabharwal/OptimizationBase.jl", rev = "as/mtk-v10" }
+OptimizationMOI = { url = "https://github.com/AayushSabharwal/Optimization.jl", subdir = "lib/OptimizationMOI", rev = "as/mtk-v10" }
+
 [targets]
-test = ["AmplNLWriter", "BenchmarkTools", "BoundaryValueDiffEqMIRK", "BoundaryValueDiffEqAscher", "ControlSystemsBase", "DataInterpolations", "DelayDiffEq", "NonlinearSolve", "ForwardDiff", "Ipopt", "Ipopt_jll", "ModelingToolkitStandardLibrary", "Optimization", "OptimizationOptimJL", "OptimizationMOI", "OrdinaryDiffEq", "OrdinaryDiffEqCore", "OrdinaryDiffEqDefault", "REPL", "Random", "ReferenceTests", "SafeTestsets", "StableRNGs", "Statistics", "SteadyStateDiffEq", "Test", "StochasticDiffEq", "Sundials", "StochasticDelayDiffEq", "Pkg", "JET", "OrdinaryDiffEqNonlinearSolve", "Logging"]
+test = ["AmplNLWriter", "BenchmarkTools", "BoundaryValueDiffEqMIRK", "BoundaryValueDiffEqAscher", "ControlSystemsBase", "DataInterpolations", "DelayDiffEq", "NonlinearSolve", "ForwardDiff", "Ipopt", "Ipopt_jll", "ModelingToolkitStandardLibrary", "Optimization", "OptimizationOptimJL", "OptimizationMOI", "OrdinaryDiffEq", "OrdinaryDiffEqCore", "OrdinaryDiffEqDefault", "REPL", "Random", "ReferenceTests", "SafeTestsets", "StableRNGs", "Statistics", "SteadyStateDiffEq", "Test", "StochasticDiffEq", "Sundials", "StochasticDelayDiffEq", "Pkg", "JET", "OrdinaryDiffEqNonlinearSolve", "Logging", "OptimizationBase"]

docs/Project.toml

@@ -40,7 +40,7 @@ Documenter = "1"
 DynamicQuantities = "^0.11.2, 0.12, 1"
 FMI = "0.14"
 FMIZoo = "1"
-ModelingToolkit = "8.33, 9"
+ModelingToolkit = "10"
 ModelingToolkitStandardLibrary = "2.19"
 NonlinearSolve = "3, 4"
 Optim = "1.7"

docs/src/basics/Composition.md

@@ -135,42 +135,34 @@ sys.y = u * 1.1
 In a hierarchical system, variables of the subsystem get namespaced by the name of the system they are in. This prevents naming clashes, but also enforces that every unknown and parameter is local to the subsystem it is used in. In some cases it might be desirable to have variables and parameters that are shared between subsystems, or even global. This can be accomplished as follows.
 
 ```julia
-@parameters a b c d e f
+@parameters a b c d
 
 # a is a local variable
 b = ParentScope(b) # b is a variable that belongs to one level up in the hierarchy
 c = ParentScope(ParentScope(c)) # ParentScope can be nested
-d = DelayParentScope(d) # skips one level before applying ParentScope
-e = DelayParentScope(e, 2) # second argument allows skipping N levels
-f = GlobalScope(f)
+d = GlobalScope(d)
 
-p = [a, b, c, d, e, f]
+p = [a, b, c, d]
 
 level0 = ODESystem(Equation[], t, [], p; name = :level0)
 level1 = ODESystem(Equation[], t, [], []; name = :level1) ∘ level0
 parameters(level1)
 #level0₊a
 #b
 #c
-#level0₊d
-#level0₊e
-#f
+#d
 level2 = ODESystem(Equation[], t, [], []; name = :level2) ∘ level1
 parameters(level2)
 #level1₊level0₊a
 #level1₊b
 #c
-#level0₊d
-#level1₊level0₊e
-#f
+#d
 level3 = ODESystem(Equation[], t, [], []; name = :level3) ∘ level2
 parameters(level3)
 #level2₊level1₊level0₊a
 #level2₊level1₊b
 #level2₊c
-#level2₊level0₊d
-#level1₊level0₊e
-#f
+#d
 ```
 
 ## Structural Simplify

src/ModelingToolkit.jl

@@ -296,7 +296,7 @@ export PDESystem
 export Differential, expand_derivatives, @derivatives
 export Equation, ConstrainedEquation
 export Term, Sym
-export SymScope, LocalScope, ParentScope, DelayParentScope, GlobalScope
+export SymScope, LocalScope, ParentScope, GlobalScope
 export independent_variable, equations, controls, observed, full_equations
 export initialization_equations, guesses, defaults, parameter_dependencies, hierarchy
 export structural_simplify, expand_connections, linearize, linearization_function,

src/inputoutput.jl

@@ -163,7 +163,7 @@ has_var(ex, x) = x ∈ Set(get_variables(ex))
     (f_oop, f_ip), x_sym, p_sym, io_sys = generate_control_function(
             sys::AbstractODESystem,
             inputs             = unbound_inputs(sys),
-            disturbance_inputs = nothing;
+            disturbance_inputs = disturbances(sys);
             implicit_dae       = false,
             simplify           = false,
         )
@@ -179,9 +179,6 @@ The return values also include the chosen state-realization (the remaining unkno
 
 If `disturbance_inputs` is an array of variables, the generated dynamics function will preserve any state and dynamics associated with disturbance inputs, but the disturbance inputs themselves will (by default) not be included as inputs to the generated function. The use case for this is to generate dynamics for state observers that estimate the influence of unmeasured disturbances, and thus require unknown variables for the disturbance model, but without disturbance inputs since the disturbances are not available for measurement. To add an input argument corresponding to the disturbance inputs, either include the disturbance inputs among the control inputs, or set `disturbance_argument=true`, in which case an additional input argument `w` is added to the generated function `(x,u,p,t,w)->rhs`.
 
-!!! note "Un-simplified system"
-    This function expects `sys` to be un-simplified, i.e., `structural_simplify` or `@mtkbuild` should not be called on the system before passing it into this function. `generate_control_function` calls a special version of `structural_simplify` internally.
-
 # Example
 
 ```
@@ -200,16 +197,18 @@ function generate_control_function(sys::AbstractODESystem, inputs = unbound_inpu
         simplify = false,
         eval_expression = false,
         eval_module = @__MODULE__,
+        check_simplified = true,
         kwargs...)
+    # Remove this when the ControlFunction gets merged.
+    if check_simplified && !iscomplete(sys) 
+        error("A completed `ODESystem` is required. Call `complete` or `structural_simplify` on the system before creating the control function.")
+    end
     isempty(inputs) && @warn("No unbound inputs were found in system.")
-
     if disturbance_inputs !== nothing
         # add to inputs for the purposes of io processing
         inputs = [inputs; disturbance_inputs]
     end
 
-    sys, _ = io_preprocessing(sys, inputs, []; simplify, kwargs...)
-
     dvs = unknowns(sys)
     ps = parameters(sys; initial_parameters = true)
     ps = setdiff(ps, inputs)
@@ -218,7 +217,7 @@ function generate_control_function(sys::AbstractODESystem, inputs = unbound_inpu
         inputs = setdiff(inputs, disturbance_inputs)
         # ps = [ps; disturbance_inputs]
     end
-    inputs = map(x -> time_varying_as_func(value(x), sys), inputs)
+    inputs = map(value, inputs)
     disturbance_inputs = unwrap.(disturbance_inputs)
 
     eqs = [eq for eq in full_equations(sys)]
@@ -257,8 +256,11 @@ function generate_control_function(sys::AbstractODESystem, inputs = unbound_inpu
     (; f, dvs, ps, io_sys = sys)
 end
 
-function inputs_to_parameters!(state::TransformationState, io)
-    check_bound = io === nothing
+"""
+Turn input variables into parameters of the system.
+"""
+function inputs_to_parameters!(state::TransformationState, inputsyms)
+    check_bound = inputsyms === nothing
     @unpack structure, fullvars, sys = state
     @unpack var_to_diff, graph, solvable_graph = structure
     @assert solvable_graph === nothing
@@ -287,7 +289,7 @@ function inputs_to_parameters!(state::TransformationState, io)
             push!(new_fullvars, v)
         end
     end
-    ninputs == 0 && return (state, 1:0)
+    ninputs == 0 && return state
 
     nvars = ndsts(graph) - ninputs
     new_graph = BipartiteGraph(nsrcs(graph), nvars, Val(false))
@@ -316,24 +318,11 @@ function inputs_to_parameters!(state::TransformationState, io)
     @set! sys.unknowns = setdiff(unknowns(sys), keys(input_to_parameters))
     ps = parameters(sys)
 
-    if io !== nothing
-        inputs, = io
-        # Change order of new parameters to correspond to user-provided order in argument `inputs`
-        d = Dict{Any, Int}()
-        for (i, inp) in enumerate(new_parameters)
-            d[inp] = i
-        end
-        permutation = [d[i] for i in inputs]
-        new_parameters = new_parameters[permutation]
-    end
-
     @set! sys.ps = [ps; new_parameters]
-
     @set! state.sys = sys
     @set! state.fullvars = new_fullvars
     @set! state.structure = structure
-    base_params = length(ps)
-    return state, (base_params + 1):(base_params + length(new_parameters)) # (1:length(new_parameters)) .+ base_params
+    return state
 end
 
 """
@@ -359,7 +348,7 @@ function get_disturbance_system(dist::DisturbanceModel{<:ODESystem})
 end
 
 """
-    (f_oop, f_ip), augmented_sys, dvs, p = add_input_disturbance(sys, dist::DisturbanceModel, inputs = nothing)
+    (f_oop, f_ip), augmented_sys, dvs, p = add_input_disturbance(sys, dist::DisturbanceModel, inputs = Any[])
 
 Add a model of an unmeasured disturbance to `sys`. The disturbance model is an instance of [`DisturbanceModel`](@ref).
 
@@ -408,13 +397,13 @@ model_outputs = [model.inertia1.w, model.inertia2.w, model.inertia1.phi, model.i
 
 `f_oop` will have an extra state corresponding to the integrator in the disturbance model. This state will not be affected by any input, but will affect the dynamics from where it enters, in this case it will affect additively from `model.torque.tau.u`.
 """
-function add_input_disturbance(sys, dist::DisturbanceModel, inputs = nothing; kwargs...)
+function add_input_disturbance(sys, dist::DisturbanceModel, inputs = Any[]; kwargs...)
     t = get_iv(sys)
     @variables d(t)=0 [disturbance = true]
     @variables u(t)=0 [input = true] # New system input
     dsys = get_disturbance_system(dist)
 
-    if inputs === nothing
+    if isempty(inputs)
         all_inputs = [u]
     else
         i = findfirst(isequal(dist.input), inputs)
@@ -429,8 +418,9 @@ function add_input_disturbance(sys, dist::DisturbanceModel, inputs = nothing; kw
            dist.input ~ u + dsys.output.u[1]]
     augmented_sys = ODESystem(eqs, t, systems = [dsys], name = gensym(:outer))
     augmented_sys = extend(augmented_sys, sys)
+    ssys = structural_simplify(augmented_sys, inputs = all_inputs, disturbance_inputs = [d])
 
-    (f_oop, f_ip), dvs, p, io_sys = generate_control_function(augmented_sys, all_inputs,
+    (f_oop, f_ip), dvs, p, io_sys = generate_control_function(ssys, all_inputs,
         [d]; kwargs...)
     (f_oop, f_ip), augmented_sys, dvs, p, io_sys
 end