type stability fixes and tests for freqresp (#642)

* use getproperty instead of getfield for LQG

* merge

* merge

* type stability fixes and tests for freqresp

Author: baggepinnen

src/freqresp.jl

@@ -39,8 +39,9 @@ _freq(w, te::Discrete) = cis(w*te.Ts)
 @autovec () function freqresp(sys::AbstractStateSpace, w_vec::AbstractVector{W}) where W <: Real
     ny, nu = size(sys)
     T = promote_type(Complex{real(eltype(sys.A))}, Complex{W})
+    PDT = PermutedDimsArray{T,3,(3,1,2),(2,3,1),Array{T,3}}
     if sys.nx == 0 # Only D-matrix
-        return PermutedDimsArray(repeat(T.(sys.D), 1, 1, length(w_vec)), (3,1,2))
+        return PDT(repeat(T.(sys.D), 1, 1, length(w_vec)))
     end
     local F, Q
     try
@@ -68,7 +69,7 @@ _freq(w, te::Discrete) = cis(w*te.Ts)
         ldiv!(A, Bc, shift = w) # B += (A - w*I)\B # solve (A-wI)X = B, storing result in B
         mul!(Ri, C, Bc, -1, 1) # use of 5-arg mul to subtract from D already in Ri. - rather than + since (A - w*I) instead of (w*I - A)
     end
-    PermutedDimsArray(R, (3,1,2)) # PermutedDimsArray doesn't allocate to perform the permutation
+    PDT(R) # PermutedDimsArray doesn't allocate to perform the permutation
 end
 
 """
@@ -82,8 +83,9 @@ freqresp_nohess
     ny, nu = size(sys)
     nx = sys.nx
     T = promote_type(Complex{real(eltype(sys.A))}, Complex{W})
+    PDT = PermutedDimsArray{T,3,(3,1,2),(2,3,1),Array{T,3}}
     if nx == 0 # Only D-matrix
-        return PermutedDimsArray(repeat(T.(sys.D), 1, 1, length(w_vec)), (3,1,2))
+        return PDT(repeat(T.(sys.D), 1, 1, length(w_vec)))
     end
     A,B,C,D = ssdata(sys)
     te = sys.timeevol
@@ -100,7 +102,7 @@ freqresp_nohess
         Bc = Ac \ B # Bc = (A - w*I)\B # avoid inplace to handle sparse matrices etc.
         mul!(Ri, C, Bc, -1, 1) # use of 5-arg mul to subtract from D already in Ri. - rather than + since (A - w*I) instead of (w*I - A)
     end
-    PermutedDimsArray(R, (3,1,2)) # PermutedDimsArray doesn't allocate to perform the permutation
+    PDT(R) # PermutedDimsArray doesn't allocate to perform the permutation
 end
 
 

test/test_freqresp.jl

@@ -10,6 +10,7 @@ G = ss([-5 0 0 0; 0 -1 -2.5 0; 0 4 0 0; 0 0 0 -6], [2 0; 0 1; 0 0; 0 2],
 @test evalfr(H, -5) == [0.0 -0.5; Inf 2.0]
 @test evalfr(H, -1) == [0.0 -0.3; 0.0 0.4]
 @test evalfr(H, 0) ≈ [0.0 0.0; 0.2 1/3]
+@inferred evalfr(H, 0)
 
 @test (@test_logs (:warn, "Got exception SingularException(4), returning Inf") evalfr(G, -6)) == [Inf Inf; Inf Inf]
 @test (@test_logs (:warn, "Got exception SingularException(1), returning Inf") evalfr(G, -5)) == [Inf Inf; Inf Inf]
@@ -33,6 +34,9 @@ resp1 = ones(ComplexF64, length(w), 1, 1)
 @test freqresp(sys1s, w) == resp1
 @test freqresp(G1, w) == resp1
 @test freqresp(H1, w) == resp1
+@inferred freqresp(sys1, w)
+@inferred freqresp(G1, w)
+@inferred freqresp(H1, w)
 
 for G in [2, 2I, 2I(2), randn(2,2)]
     @test dcgain(G) == G
@@ -58,6 +62,11 @@ resp2 = reshape((im*w .+ 2)./(im*w  .+ 1), length(w), 1, 1)
 @test freqresp(G2, w) == resp2
 @test freqresp(H2, w) == resp2
 
+@inferred freqresp(sys2, w)
+@inferred freqresp(sys2s, w)
+@inferred freqresp(G2, w)
+@inferred freqresp(H2, w)
+
 
 ## Complex-coefficient system
 sys3 = ss(-1+im, 1, (1-im), (1-im))