Clamp default frequency grid to Nyquist for discrete diskmargin and mu

src/mimo_diskmargin.jl

@@ -247,7 +247,9 @@ function structured_singular_value(M0::LTISystem, w::AbstractVector; kwargs...)
 end
 
 function structured_singular_value(M0::LTISystem; kwargs...)
-    w = exp10.(LinRange(-3, 3, 1500))
+    # Clamp the upper end of the default grid to the Nyquist frequency for discrete `M0`.
+    u = isdiscrete(M0) ? min(3.0, log10(π/M0.Ts)) : 3.0
+    w = exp10.(LinRange(-3, u, 1500))
     μ = structured_singular_value(M0, w; kwargs...)
     maximum(μ)
 end
@@ -290,6 +292,10 @@ Return the smallest simultaneous diskmargin over the grid l:u
 See also [`ncfmargin`](@ref).
 """
 function sim_diskmargin(L, σ::Real=0, l::Real=1e-3, u::Real=1e3; kwargs...)
+    # For discrete L the default upper bound 1e3 typically exceeds the Nyquist frequency π/Ts.
+    # Above Nyquist `freqresp` aliases (cis(ω·Ts) is 2π/Ts-periodic) so margins computed there
+    # reflect the response at an aliased frequency rather than the intended one.
+    u = isdiscrete(L) ? min(u, π/L.Ts) : u
     m = sim_diskmargin(L, σ, exp10.(LinRange(log10(l), log10(u), 500)); kwargs...)
     m = argmin(d->d.α, m)
 end

test/test_diskmargin.jl

@@ -124,6 +124,25 @@ let
     @test maximum(abs, αs_c .- αs_d) < 0.02
 end
 
+# Discrete-time sim_diskmargin: default frequency grid must respect Nyquist.
+# Regression: with slow sampling (here Ts=1.0, Nyquist ≈ 3.14 rad/s) the default
+# upper bound u=1e3 used to alias under freqresp. The default-grid result must
+# agree with an explicit Nyquist-respecting grid.
+let
+    a = 10
+    Pc = ss([0 a; -a 0], I(2), [1 a; -a 1], 0)
+    Ts = 1.0
+    Pd = c2d(Pc, Ts)
+    Kd = ss(1.0I(2), Ts)
+    dm_default = sim_diskmargin(Pd, Kd)
+    @test dm_default isa Diskmargin
+    @test dm_default.ω0 ≤ π/Ts + sqrt(eps())
+    wgrid = exp10.(LinRange(log10(1e-3), log10(π/Ts), 500))
+    dms_explicit = sim_diskmargin(Pd, Kd, 0, wgrid)
+    αs = [d.α for d in dms_explicit]
+    @test dm_default.α ≈ minimum(αs) rtol=1e-3
+end
+
 ## MIMO
 
 a = 10