Merge pull request #298 from SciML/fb/ap_updates

Analysis point fixes

Author: baggepinnen

src/Blocks/analysis_points.jl

@@ -18,6 +18,10 @@ function Base.hash(ap::AnalysisPoint, seed::UInt)
     h3 ⊻ (0xd29cdc51aa6562d4 % UInt)
 end
 
+function ModelingToolkit.get_unit(ap::AnalysisPoint)
+    ModelingToolkit.unitless
+end
+
 function ap_var(sys)
     if hasproperty(sys, :u)
         # collect to turn symbolic arrays into arrays of symbols
@@ -28,6 +32,23 @@ function ap_var(sys)
     error("Could not determine the analysis-point variable in system $(nameof(sys)). To use an analysis point, apply it to a connection between two causal blocks containing connectors of type `RealInput/RealOutput` from ModelingToolkitStandardLibrary.Blocks.")
 end
 
+"""
+    find_analysis_points(sys)
+
+Return a list of all analysis points in `sys`. If none are found, the list is empty.
+"""
+function find_analysis_points(sys)
+    sys = ModelingToolkit.flatten(sys)
+    eqs = equations(sys)
+    aps = []
+    for eq in eqs
+        if eq.rhs isa AnalysisPoint
+            push!(aps, eq.rhs)
+        end
+    end
+    aps
+end
+
 """
     AnalysisPoint(in, out, name::Symbol)
     AnalysisPoint(in, out; name::Symbol)
@@ -418,16 +439,26 @@ function ModelingToolkit.linearization_function(sys::ModelingToolkit.AbstractSys
             ui = get_perturbation_var(ap_var(ap.out), "u")
             push!(multiplicities_u, length(ui)) # one ap may yield several new vars
             append!(u, ui)
-            [ap_var(ap.out) .~ ap_var(ap.in) + ui;], ui
+            if loop_openings !== nothing && ap.name ∈ loop_openings
+                # In thise case, we break the existing connection.
+                [ap_var(ap.out) .~ ui;], ui
+            else
+                [ap_var(ap.out) .~ ap_var(ap.in) + ui;], ui
+            end
             #input.in.u ~ 0] # We only need to ground one of the ends, hence not including this equation
         elseif output_name isa SymOrVec && namespaced_ap_match(ap, ns, output_name, nothing)
             push!(namespace_y, ns) # Save the namespace to make it available for renamespace below
             push!(aps_y, ap)
             yi = get_perturbation_var(ap_var(ap.in), "y")
             push!(multiplicities_y, length(yi))
             append!(y, yi)
-            [ap_var(ap.in) .~ yi;
-             ap_var(ap.out) .~ ap_var(ap.in)], yi
+            if loop_openings !== nothing && ap.name ∈ loop_openings
+                [ap_var(ap.in) .~ yi;
+                 ap_var(ap.out) .~ 0], yi # In thise case, we break the existing connection.
+            else
+                [ap_var(ap.in) .~ yi;
+                 ap_var(ap.out) .~ ap_var(ap.in)], yi
+            end
         else # loop opening
             [ap_var(ap.out) .~ 0;], []
         end

test/Blocks/test_analysis_points.jl

@@ -242,6 +242,14 @@ Si = ss(matrices...)
     t,
     systems = [P_inner, feedback, ref])
 
+@test_nowarn Blocks.find_analysis_points(sys_inner)
+P_not_broken, _ = linearize(sys_inner, :u, :y)
+@test P_not_broken.A[] == -2
+P_broken, _ = linearize(sys_inner, :u, :y, loop_openings = [:u])
+@test P_broken.A[] == -1
+P_broken, _ = linearize(sys_inner, :u, :y, loop_openings = [:y])
+@test P_broken.A[] == -1
+
 Sinner = sminreal(ss(get_sensitivity(sys_inner, :u)[1]...))
 
 @named sys_inner = ODESystem(