fix: fix linearization tests

Author: vyudu

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 = Any[];
             implicit_dae       = false,
             simplify           = false,
         )
@@ -289,7 +289,7 @@ function inputs_to_parameters!(state::TransformationState, inputsyms)
             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))
@@ -318,14 +318,13 @@ function inputs_to_parameters!(state::TransformationState, inputsyms)
     @set! sys.unknowns = setdiff(unknowns(sys), keys(input_to_parameters))
     ps = parameters(sys)
 
-    if io !== nothing
-        inputs, = io
+    if inputsyms !== nothing
         # 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]
+        permutation = [d[i] for i in inputsyms]
         new_parameters = new_parameters[permutation]
     end
 
@@ -334,8 +333,7 @@ function inputs_to_parameters!(state::TransformationState, inputsyms)
     @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
 
 """
@@ -361,7 +359,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).
 
@@ -410,7 +408,7 @@ 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

src/linearization.jl

@@ -120,7 +120,7 @@ function linearization_function(sys::AbstractSystem, inputs = unbound_inputs(sys
     end
 
     lin_fun = LinearizationFunction(
-        diff_idxs, alge_idxs, length(unknowns(sys)),
+        diff_idxs, alge_idxs, inputs, length(unknowns(sys)),
         prob, h, u0 === nothing ? nothing : similar(u0), uf_jac, h_jac, pf_jac,
         hp_jac, initializealg, initialization_kwargs)
     return lin_fun
@@ -397,7 +397,7 @@ Construct a `LinearizationProblem` for linearizing the system `sys` with the giv
 All other keyword arguments are forwarded to `linearization_function`.
 """
 function LinearizationProblem(sys::AbstractSystem, inputs, outputs; t = 0.0, kwargs...)
-    linfun, _ = linearization_function(sys, inputs, outputs; kwargs...)
+    linfun = linearization_function(sys, inputs, outputs; kwargs...)
     return LinearizationProblem(linfun, t)
 end
 
@@ -764,7 +764,7 @@ Permute the state representation of `sys` obtained from [`linearize`](@ref) so t
 Example:
 
 ```
-lsys, ssys = linearize(pid, [reference.u, measurement.u], [ctr_output.u])
+lsys = linearize(pid, [reference.u, measurement.u], [ctr_output.u])
 desired_order = [int.x, der.x] # Unknowns that are present in unknowns(ssys)
 lsys = ModelingToolkit.reorder_unknowns(lsys, unknowns(ssys), desired_order)
 ```

src/systems/analysis_points.jl

@@ -596,9 +596,13 @@ Add an input without an additional output variable.
 PerturbOutput(ap::AnalysisPoint) = PerturbOutput(ap, false)
 
 function apply_transformation(tf::PerturbOutput, sys::AbstractSystem)
+    @show "ok"
+    @show tf.ap
     modify_nested_subsystem(sys, tf.ap) do ap_sys
         # get analysis point
+        @show tf.ap
         ap_idx = analysis_point_index(ap_sys, tf.ap)
+        @show ap_idx
         ap_idx === nothing &&
             error("Analysis point $(nameof(tf.ap)) not found in system $(nameof(sys)).")
         # modified equations

src/systems/diffeqs/odesystem.jl

@@ -472,7 +472,7 @@ Generates a function that computes the observed value(s) `ts` in the system `sys
 - `eval_expression = false`: If true and `expression = false`, evaluates the returned function in the module `eval_module`
 - `output_type = Array` the type of the array generated by a out-of-place vector-valued function
 - `param_only = false` if true, only allow the generated function to access system parameters
-- `inputs = nothing` additinoal symbolic variables that should be provided to the generated function
+- `inputs = Any[]` additional symbolic variables that should be provided to the generated function
 - `checkbounds = true` checks bounds if true when destructuring parameters
 - `op = Operator` sets the recursion terminator for the walk done by `vars` to identify the variables that appear in `ts`. See the documentation for `vars` for more detail.
 - `throw = true` if true, throw an error when generating a function for `ts` that reference variables that do not exist.
@@ -502,8 +502,8 @@ For example, a function `g(op, unknowns, p..., inputs, t)` will be the in-place
 an array of inputs `inputs` is given, and `param_only` is false for a time-dependent system.
 """
 function build_explicit_observed_function(sys, ts;
-        inputs = nothing,
-        disturbance_inputs = nothing,
+        inputs = Any[],
+        disturbance_inputs = Any[],
         disturbance_argument = false,
         expression = false,
         eval_expression = false,
@@ -576,13 +576,13 @@ function build_explicit_observed_function(sys, ts;
     else
         (unknowns(sys),)
     end
-    if inputs === nothing
+    if isempty(inputs)
         inputs = ()
     else
         ps = setdiff(ps, inputs) # Inputs have been converted to parameters, remove those from the parameter list
         inputs = (inputs,)
     end
-    if disturbance_inputs !== nothing
+    if !isempty(disturbance_inputs)
         # Disturbance inputs may or may not be included as inputs, depending on disturbance_argument
         ps = setdiff(ps, disturbance_inputs)
     end

src/systems/systems.jl

@@ -29,8 +29,8 @@ topological sort of the observed equations in `sys`.
 function mtkbuild(
         sys::AbstractSystem; additional_passes = [], simplify = false, split = true,
         allow_symbolic = false, allow_parameter = true, conservative = false, fully_determined = true,
-        inputs = nothing, outputs = nothing,
-        disturbance_inputs = nothing,
+        inputs = Any[], outputs = Any[],
+        disturbance_inputs = Any[],
         kwargs...)
     isscheduled(sys) && throw(RepeatedStructuralSimplificationError())
     newsys′ = __structural_simplification(sys; simplify,
@@ -74,8 +74,8 @@ function __structural_simplification(sys::SDESystem, args...; kwargs...)
 end
 
 function __structural_simplification(sys::AbstractSystem; simplify = false,
-        inputs = nothing, outputs = nothing,
-        disturbance_inputs = nothing,
+        inputs = Any[], outputs = Any[],
+        disturbance_inputs = Any[],
         kwargs...)
     sys = expand_connections(sys)
     state = TearingState(sys; sort_eqs)

src/systems/systemstructure.jl

@@ -659,8 +659,8 @@ end
 
 function structural_simplification!(state::TearingState; simplify = false,
         check_consistency = true, fully_determined = true, warn_initialize_determined = true,
-        inputs = nothing, outputs = nothing, 
-        disturbance_inputs = nothing,
+        inputs = Any[], outputs = Any[], 
+        disturbance_inputs = Any[],
         kwargs...)
 
     if state.sys isa ODESystem
@@ -672,7 +672,7 @@ function structural_simplification!(state::TearingState; simplify = false,
         cont_inputs = [inputs; clocked_inputs[continuous_id]]
         sys = _structural_simplification!(tss[continuous_id]; simplify,
             check_consistency, fully_determined,
-            cont_inputs, outputs, disturbance_inputs,
+            inputs = cont_inputs, outputs, disturbance_inputs,
             kwargs...)
         if length(tss) > 1
             if continuous_id > 0
@@ -690,7 +690,7 @@ function structural_simplification!(state::TearingState; simplify = false,
                 end
                 disc_inputs = [inputs; clocked_inputs[i]]
                 ss, = _structural_simplification!(state; simplify, check_consistency,
-                    disc_inputs, outputs, disturbance_inputs,
+                    inputs = disc_inputs, outputs, disturbance_inputs,
                     fully_determined, kwargs...)
                 append!(appended_parameters, inputs[i], unknowns(ss))
                 discrete_subsystems[i] = ss
@@ -717,8 +717,8 @@ end
 function _structural_simplification!(state::TearingState; simplify = false,
         check_consistency = true, fully_determined = true, warn_initialize_determined = false,
         dummy_derivative = true,
-        inputs = nothing, outputs = nothing,
-        disturbance_inputs = nothing,
+        inputs = Any[], outputs = Any[],
+        disturbance_inputs = Any[],
         kwargs...)
     if fully_determined isa Bool
         check_consistency &= fully_determined

test/downstream/linearize.jl

@@ -25,14 +25,14 @@ lsys3 = linearize(sys, [r], [y]; autodiff = AutoFiniteDiff())
 @test lsys.C[] == lsys2.C[] == lsys3.C[] == 1
 @test lsys.D[] == lsys2.D[] == lsys3.D[] == 0
 
-lsys, ssys = linearize(sys, [r], [r])
+lsys = linearize(sys, [r], [r])
 
 @test lsys.A[] == -2
 @test lsys.B[] == 1
 @test lsys.C[] == 0
 @test lsys.D[] == 1
 
-lsys, ssys = linearize(sys, r, r) # Test allow scalars
+lsys = linearize(sys, r, r) # Test allow scalars
 
 @test lsys.A[] == -2
 @test lsys.B[] == 1
@@ -89,19 +89,19 @@ connections = [f.y ~ c.r # filtered reference to controller reference
 @named cl = ODESystem(connections, t, systems = [f, c, p])
 cl = mtkbuild(cl, inputs = [f.u], outputs = [p.x])
 
-lsys0, ssys = linearize(cl)
+lsys0 = linearize(cl, [f.u], [p.x])
 desired_order = [f.x, p.x]
-lsys = ModelingToolkit.reorder_unknowns(lsys0, unknowns(ssys), desired_order)
-lsys1, ssys = linearize(cl; autodiff = AutoFiniteDiff())
-lsys2 = ModelingToolkit.reorder_unknowns(lsys1, unknowns(ssys), desired_order)
+lsys = ModelingToolkit.reorder_unknowns(lsys0, unknowns(cl), desired_order)
+lsys1 = linearize(cl, [f.u], [p.x]; autodiff = AutoFiniteDiff())
+lsys2 = ModelingToolkit.reorder_unknowns(lsys1, unknowns(cl), desired_order)
 
 @test lsys.A == lsys2.A == [-2 0; 1 -2]
 @test lsys.B == lsys2.B == reshape([1, 0], 2, 1)
 @test lsys.C == lsys2.C == [0 1]
 @test lsys.D[] == lsys2.D[] == 0
 
 ## Symbolic linearization
-lsyss, _ = ModelingToolkit.linearize_symbolic(cl, [f.u], [p.x])
+lsyss = ModelingToolkit.linearize_symbolic(cl, [f.u], [p.x])
 
 @test ModelingToolkit.fixpoint_sub(lsyss.A, ModelingToolkit.defaults(cl)) == lsys.A
 @test ModelingToolkit.fixpoint_sub(lsyss.B, ModelingToolkit.defaults(cl)) == lsys.B
@@ -116,11 +116,12 @@ Nd = 10
 @named pid = LimPID(; k, Ti, Td, Nd)
 
 @unpack reference, measurement, ctr_output = pid
-lsys0, ssys = linearize(pid, [reference.u, measurement.u], [ctr_output.u];
+pid = mtkbuild(pid, inputs = [reference.u, measurement.u], outputs = [ctr_output.u])
+lsys0 = linearize(pid, [reference.u, measurement.u], [ctr_output.u];
     op = Dict(reference.u => 0.0, measurement.u => 0.0))
 @unpack int, der = pid
 desired_order = [int.x, der.x]
-lsys = ModelingToolkit.reorder_unknowns(lsys0, unknowns(ssys), desired_order)
+lsys = ModelingToolkit.reorder_unknowns(lsys0, unknowns(pid), desired_order)
 
 @test lsys.A == [0 0; 0 -10]
 @test lsys.B == [2 -2; 10 -10]
@@ -150,12 +151,12 @@ lsys = ModelingToolkit.reorder_unknowns(lsys, desired_order, reverse(desired_ord
 
 ## Test that there is a warning when input is misspecified
 if VERSION >= v"1.8"
-    @test_throws "Some specified inputs were not found" linearize(pid,
+    @test_throws "Some parameters are missing from the variable map." linearize(pid,
         [
             pid.reference.u,
             pid.measurement.u
         ], [ctr_output.u])
-    @test_throws "Some specified outputs were not found" linearize(pid,
+    @test_throws "Some parameters are missing from the variable map." linearize(pid,
         [
             reference.u,
             measurement.u
@@ -186,15 +187,16 @@ function saturation(; y_max, y_min = y_max > 0 ? -y_max : -Inf, name)
     ODESystem(eqs, t, name = name)
 end
 @named sat = saturation(; y_max = 1)
+sat = mtkbuild(sat, inputs = [u], outputs = [y])
 # inside the linear region, the function is identity
 @unpack u, y = sat
-lsys, ssys = linearize(sat, [u], [y])
+lsys = linearize(sat, [u], [y])
 @test isempty(lsys.A) # there are no differential variables in this system
 @test isempty(lsys.B)
 @test isempty(lsys.C)
 @test lsys.D[] == 1
 
-@test_skip lsyss, _ = ModelingToolkit.linearize_symbolic(sat, [u], [y]) # Code gen replaces ifelse with if statements causing symbolic evaluation to fail
+@test_skip lsyss = ModelingToolkit.linearize_symbolic(sat, [u], [y]) # Code gen replaces ifelse with if statements causing symbolic evaluation to fail
 # @test substitute(lsyss.A, ModelingToolkit.defaults(sat)) == lsys.A
 # @test substitute(lsyss.B, ModelingToolkit.defaults(sat)) == lsys.B
 # @test substitute(lsyss.C, ModelingToolkit.defaults(sat)) == lsys.C
@@ -267,9 +269,9 @@ closed_loop = ODESystem(connections, t, systems = [model, pid, filt, sensor, r,
         filt.x => 0.0,
         filt.xd => 0.0
     ])
-closed_loop = mtkbuild(closed_loop, inputs = :r, outputs = :y)
+closed_loop = mtkbuild(closed_loop)
 
-@test_nowarn linearize(closed_loop; warn_empty_op = false)
+@test_nowarn linearize(closed_loop, :r, :y; warn_empty_op = false)
 
 # https://discourse.julialang.org/t/mtk-change-in-linearize/115760/3
 @mtkmodel Tank_noi begin