simplify implementation of named_sensitivity_function

Author: baggepinnen

docs/src/batch_linearization.md

@@ -175,7 +175,7 @@ We will explore these options below, starting with the first option, breaking th
 ```@example BATCHLIN
 kwargs = (; adaptive=false, legend=false)
 plants, _ = trajectory_ss(closed_loop, closed_loop.u, closed_loop.y, sol; t=timepoints, verbose=true, loop_openings=[closed_loop.u]);
-controllersy, ssy = trajectory_ss(closed_loop, closed_loop.r, closed_loop.u, sol; t=timepoints, verbose=true, loop_openings=[closed_loop.y]);
+controllersy, ssy, ops3, resolved_ops3 = trajectory_ss(closed_loop, closed_loop.r, closed_loop.u, sol; t=timepoints, verbose=true, loop_openings=[closed_loop.y]);
 
 closed_loopsy = feedback.(plants .* controllersy)
 bodeplot(closed_loopsy, w; title="Loop open at y", kwargs...)

src/ode_system.jl

@@ -197,12 +197,13 @@ function named_sensitivity_function(
     fun,
     sys::ModelingToolkit.AbstractSystem,
     inputs, args...;
+    descriptor = true,
     kwargs...,
 )
 
     inputs = vcat(inputs)
     inputs = map(inputs) do inp
-        if inp isa System
+    if inp isa System
             @variables u(t)
             if u ∈ Set(unknowns(inp))
                 inp.u
@@ -219,23 +220,10 @@ function named_sensitivity_function(
     unames = symstr.(inputs)
     fm(x) = convert(Matrix{Float64}, x)
     if nu > 0 && size(matrices.B, 2) == 2nu
-        nx = size(matrices.A, 1)
-         # This indicates that input derivatives are present
-        duinds = findall(any(!iszero, eachcol(matrices.B[:, nu+1:end])))
-        B̄ = matrices.B[:, duinds .+ nu]
-        ndu = length(duinds)
-        B = matrices.B[:, 1:nu]
-        Iu = duinds .== (1:nu)'
-        E = [I(nx) -B̄; zeros(ndu, nx+ndu)]
-
-        Ae = cat(matrices.A, -I(ndu), dims=(1,2))
-        Be = [B; Iu]
-        Ce = [fm(matrices.C) zeros(ny, ndu)]
-        De = fm(matrices.D[:, 1:nu])
-        dsys = dss(Ae, E, Be, Ce, De)
-        lsys = ss(RobustAndOptimalControl.DescriptorSystems.dss2ss(dsys)[1])
-        # unames = [unames; Symbol.("der_" .* string.(unames))]
-        # sys = ss(matrices...)
+        # This indicates that input derivatives are present
+        duinds = findall(any(!iszero, eachcol(matrices.B[:, nu+1:end]))) .+ nu
+        u2du = (1:nu) .=> duinds # This maps inputs to their derivatives
+        lsys = causal_simplification(matrices, u2du; descriptor)
     else
         lsys = ss(matrices...)
     end