Solving model detection problems

emdsyn Exact synthesis of model detection filters.
amdsyn Approximate synthesis of model detection filters.

FaultDetectionTools.emdsyn — Function

emdsyn(sysm::MDMModel; rdim, MDSelect, nullspace = true, simple = false, minimal = true, 
                       emtest = false, normalize = false, fast = true, 
                       sdeg, smarg, poles, HDesign, MDtol, MDGainTol, MDfreq, 
                       tcond, offset, atol, atol1, atol2, atol3, rtol) 
                       -> (Q::MDFilter, R::MDFilterIF, info)

Solve the exact model detection problem (EMDP) for a given multiple synthesis model sysm::MDMModel. The computed stable and proper filter objects Q and R contain the model detection filter, representing the solution of the EMDP, and its internal form, respectively.

The returned named tuple info, with the components info.tcond, info.degs, info.MDperf, and info.HDesign, contains additional synthesis related information (see below).

For a multiple synthesis model sysm containing N stable component models, the i-th component system sysm.sys[i] has a descriptor realization of the form sysm.sys[i] = (Ai-λEi,Bi,Ci,Di) and the corresponding input-output form is

  yi = Gui(λ)*u + Gdi(λ)*di + Gwi(λ)*wi + Gvi(λ)*vi

with the Laplace- or Z-transformed plant outputs yi, control inputs u, disturbance inputs di, noise inputs wi, and auxiliary inputs vi, and with Gui(λ), Gdi(λ), Gwi(λ) and Gvi(λ), the corresponding transfer function matrices.

The model detection filter object Q, contains in Q.sys the resulting bank of N filters. The i-th filter Q.sys[i] is in a standard state-space form and generates r_i, the i-th component (scalar or vector) of the overall residual vector r := [r_1; r_2; ...; r_N]. The corresponding input-output (implementation) form of the i-th filter is

        r_i = Qyi(λ)*y + Qui(λ)*u   ,

where Qyi(λ) and Qui(λ) are the transfer function matrices from the output and control inputs to the i-th residual component. The dimensions of output and control inputs are contained in the integers Q.ny and Q.mu, respectively.

The model detection filter internal form object R, contains R.sys, the resulting array of N × N internal form of the filters. The (i,j)-th component filter R.sys[i,j] is in a standard state-space form, generates the residual signal r_ij, and corresponds to the input-output form

   r_ij = Ruij(λ)*u + Rdij(λ)*dj + Rwij(λ)*wj + Rvij(λ)*vj ,

where

   | Ruij(λ) Rdij(λ) Rwij(λ) Raij(λ) | := |Qyi(λ) Qui(λ)]*| Guj(λ) Gdj(λ) Gwj(λ) Gaj(λ) |. 
                                                          |   I    0      0      0      |

The solution of the EMDP ensures that for the i-th filter, Ruii(λ) = 0, Rdii(λ) = 0, and [Ruij(λ) Rdij(λ)] $\neq$ 0 for j $\neq$ i.

Various user options can be specified via keyword arguments as follows:

MDSelect = ifilt specifies in the vector ifilt the indices of the desired filters to be designed (default: ifilt = 1:N)

If minimal = true (default), least order filter synthesis is performed to determine each of the component filters Q.sys[i] for i = 1, ...,N and R.sys[i,j] for i, j = 1, ...,N while with minimal = false full order synthesis is performed.

If HDesign = H, then H is an N-dimensional array of full row rank or empty design matrices H = [H_1, ..., H_N], where H_i is the design matrix employed for the synthesis of the i-th component filter (default: HDesign = missing)

rdim = q specifies the vector q, whose i-th component q[i] specifies the number of residual outputs for the i-th component filter Q.sys[i]. If q is a scalar, then a vector rdim with all components equal to q is assumed. The default value of q[i] is chosen as follows: if HDesign = missing or H_i is empty then q[i] = 1, if minimal = true, or q[i] is the number of the nullspace basis vectors used for the synthesis of Q.sys [i], if minimal = false; if H_i specifies a full row rank design matrix, then q[i] is the row dimension of H_i.

If emdtest = false (default), only the input channels are used for model detectability tests. If emdtest = true, extended model detectability tests are performed using both control and disturbance input channels.

MDfreq = ω specifies a vector of real frequency values or a scalar real frequency value for strong model detectability checks (default: MDfreq = missing).

If nullspace = false (default), a full-order observer based nullspace basis is used for the synthesis of the i-th filter whenever the i-th component system has no disturbance inputs. Otherwise, minimal proper nullspace bases are used. If nullspace = true, minimal proper nullspace bases are used for the synthesis of the model detection filters.

If simple = true, simple proper nullspace bases are emplyed for synthesis. The orders of the basis vectors employed for the synthesis of the i-th filter are provided in info.deg[i]. If simple = false (default), then no simple bases are computed.

offset = β specifies the boundary offset β to assess the stability of poles. Accordingly, for the stability of a continuous-time system all real parts of poles must be at most -β, while for the stability of a discrete-time system all moduli of poles must be at most 1-β. The default value used for β is sqrt(ϵ), where ϵ is the working machine precision.

smarg = α specifies the stability margin which defines the stability domain Cs of poles, as follows: for a continuous-time system, Cs is the set of complex numbers with real parts at most α, while for a discrete-time system, Cs is the set of complex numbers with moduli at most α < 1 (i.e., the interior of a disc of radius α centered in the origin). If smarg is missing, then the employed default values are α = -β for a continuous-time system and α = 1-β for a discrete-time system, where β is the boundary offset specified by the keyword argument offset = β.

sdeg = γ is the prescribed stability degree for the poles of the filters Q and R (default: γ = -0.05 for the real parts of poles for a continuous-time system and γ = 0.95 for the magnitudes of poles for a discrete-time system).

poles = v specifies a complex vector v containing a complex conjugate set of desired poles within the stability domain Cs to be assigned for the filters Q and R (default: poles = missing).

tcond = tcmax specifies the maximum alowed condition number tcmax of the employed non-orthogonal transformations (default: tcmax = 1.e4).

MDtol = tol1 specifies the threshold tol1 for model detectability checks (default: tol1 = 0.0001).

MDGainTol = tol2 specifies the threshold tol2 for strong model detectability checks (default: tol2 = 0.01).

If normalize = false (default), the i-th component filter Q.sys[i] is scaled such that the minimum gain of R.sys[i,j] for j = 1, ..., N, j $\neq$ i, is equal to one. If normalize = true, the standard normalization of component filters is performed to ensure equal gains for R.sys[1,j] and R.sys[j,1]. This option is only possible if ifilt[1] = 1 (see MDSelect).

The rank determinations in the performed reductions are based on rank revealing QR-decompositions with column pivoting if fast = true (default) or the more reliable SVD-decompositions if fast = false.

The keyword arguments atol1, atol2 and rtol, specify, respectively, the absolute tolerance for the nonzero elements of the state-space matrices Ai, Bi, Ci and Di, the absolute tolerance for the nonzero elements of Ei, and the relative tolerance for the nonzero elements of all above matrices. The default relative tolerance is ni*ϵ, where ϵ is the working machine epsilon and ni is the orders of the system sysm.sys[i]. The keyword argument atol3 is an absolute tolerance for observability tests (default: internally determined value). The keyword argument atol can be used to simultaneously set atol1 = atol, atol2 = atol and atol3 = atol.

The resulting named tuple info contains (tcond, degs, MDperf, HDesign), where:

info.tcond is the maximum of the condition numbers of the employed non-orthogonal transformation matrices; a warning is issued if info.tcond >= tcmax;

info.degs is an N-dimensional vector, whose i-th component is an integer vector containing the degrees of the basis vectors of the employed simple nullspace basis for the synthesis of the i-th filter component, if simple = true, and the degrees of the basis vectors of an equivalent polynomial nullspace basis, if simple = false;

info.MDperf is an N × N array containing the resulting distance-mapping performance, representing the peak gains of the associated internal representations of the model detection filters. If the (i,j)-th component filter R.sys[i,j] has the input-output form

 rij = Ruij(λ)*u + Rdij(λ)*dj + Rwij(λ)*wj + Rvij(λ)*vj ,

then, the (i,j)-th performance gain is evaluated as info.MDperf[i,j] = ||Ruij(λ)||∞ if emdtest = false (default) or info.MDperf[i,j] = ||[Ruij(λ) Rdij(λ)]||∞ if emdtest = true. If MDfreq = ω, then each gain info.MDperf[i,j] represents the maximum of 2-norm pointwise gains evaluated for all frequencies in ω.

info.HDesign is an N-dimensional vector, whose i-th component is the design matrix H_i employed for the synthesis of the i-th model detection filter.

Method: The Procedure EMD from [1] is implemented to solve the exact model detection problem. For more details on the least order synthesis of model detection filters see [2].

References:

[1] A. Varga, Solving Fault Diagnosis Problems - Linear Synthesis Techniques. Springer Verlag, 2017; sec. 5.2.

[2] A. Varga, Least order fault and model detection using multi-models. IEEE CDC'09, Shanghai, China, 2009.

FaultDetectionTools.amdsyn — Function

amdsyn(sysm::MDMModel; rdim, MDSelect, nullspace = true, simple = false, minimal = true, 
                       emtest = false, normalize = false, fast = true, 
                       sdeg, smarg, poles, HDesign, MDtol, MDGainTol, MDfreq, 
                       epsreg = 0.1, sdegzer, nonstd = 1, rtolinf= 0.0001, 
                       tcond, offset, atol, atol1, atol2, atol3, rtol) 
                       -> (Q::MDFilter, R::MDFilterIF, info)

Solve the approximate model detection problem (AMDP) for a given multiple synthesis model sysm::MDMModel. The computed stable and proper filter objects Q and R contain the model detection filter, representing the solution of the AMDP, and its internal form, respectively.

The returned named tuple info, with the components info.tcond, info.degs, info.MDperf, info.MDgap, info.HDesign and info.nonstandard, contains additional synthesis related information (see below).

  yi = Gui(λ)*u + Gdi(λ)*di + Gwi(λ)*wi + Gvi(λ)*vi

        r_i = Qyi(λ)*y + Qui(λ)*u ,

where Qyi(λ) and Qui(λ) are the transfer function matrices from the measured outputs and control inputs to the i-th residual component. The dimensions of measured outputs and control inputs are contained in the integers Q.ny and Q.mu, respectively.

The model detection filter internal form object R, contains R.sys, the resulting array of N × N internal form filters. The (i,j)-th component filter R.sys[i,j] is in a standard state-space form, generates the residual signal r_ij, and corresponds to the input-output form

   r_ij = Ruij(λ)*u + Rdij(λ)*dj + Rwij(λ)*wj + Rvij(λ)*vj ,

where

   | Ruij(λ) Rdij(λ) Rwij(λ) Raij(λ) | := |Qyi(λ) Qui(λ)]*| Guj(λ) Gdj(λ) Gwj(λ) Gaj(λ) |. 
                                                          |   I    0      0      0      |

The solution of the AMDP ensures that for the i-th filter, Ruii(λ) = 0, Rdii(λ) = 0, [Ruij(λ) Rdij(λ)] $\neq$ 0 for j $\neq$ i, and, additionally, the maximization of the achieved gaps (see description of info.MDgap).

Various user options can be specified via keyword arguments as follows:

MDSelect = ifilt specifies in the vector ifilt the indices of the desired filters to be designed (default: ifilt = 1:N)

MDfreq = ω specifies a vector of real frequency values or a scalar real frequency value for strong model detectability checks (default: MDfreq = missing).

tcond = tcmax specifies the maximum alowed condition number tcmax of the employed non-orthogonal transformations (default: tcmax = 1.e4).

MDtol = tol1 specifies the threshold tol1 for model detectability checks (default: tol1 = 0.0001).

MDGainTol = tol2 specifies the threshold tol2 for strong model detectability checks (default: tol2 = 0.01).

epsreg = ϵ specifies the value of the regularization parameter ϵ (default: ϵ = 0.1)

sdegzer = δ specifies the prescribed stability degree δ for zeros shifting (default: δ = −0.05 for a continuous-time system sysf.sys and δ = 0.95 for a discrete-time system sysf.sys).

nonstd = job specifies the option to handle nonstandard optimization problems as follows:

  job = 1 – use the quasi-co-outer–co-inner factorization (default);
  job = 2 – use the modified co-outer–co-inner factorization with the
            regularization parameter `ϵ`;
  job = 3 – use the Wiener-Hopf type co-outer–co-inner factorization;
  job = 4 – use the Wiener-Hopf type co-outer-co-inner factorization with
            zero shifting of the non-minimum phase factor using the
            stabilization parameter `δ`;
  job = 5 – use the Wiener-Hopf type co-outer-co-inner factorization with
            the regularization of the non-minimum phase factor using the
            regularization parameter `ϵ`.

The resulting named tuple info contains (tcond, degs, MDperf, MDgap, HDesign, nonstandard), where:

info.tcond is the maximum of the condition numbers of the employed non-orthogonal transformation matrices; a warning is issued if info.tcond >= tcmax;

 rij = Ruij(λ)*u + Rdij(λ)*dj + Rwij(λ)*wj + Rvij(λ)*vj ,

info.MDgap is an N-dimensional vector of achieved noise gaps, where the i-th gap is computed as info.MDgap[i] = β[i]/γ[i], where β[i] = min(||Rij(λ)||∞ for i $\neq$ j) and γ[i] = ||Rwii(λ)||∞. Rij(λ) is defined as Rij(λ) = Ruij(λ) if emtest = false (default) or Rij(λ) = [Ruij(λ) Rdij(λ)] if emtest = true. If MDfreq = ω, where ω is a given vector of real frequency values, then each gain β[i] represents the minimum of the maximum of 2-norm pointwise gains evaluated for all frequencies in ω.

info.HDesign is an N-dimensional vector, whose i-th component is the design matrix H_i employed for the synthesis of the i-th model detection filter.

info.nonstandard is an N-dimensional vector, whose i-th component is true if the synthesis of the i-th filter involved a nonstandard optimization problem and is false otherwise.

Method: An extension of the Procedure AMD from [1] is implemented to solve the approximate model detection problem. This procedure relies on the nullspace-based synthesis method of model detection filters proposed in [2].

References:

[1] A. Varga, Solving Fault Diagnosis Problems - Linear Synthesis Techniques. Springer Verlag, 2017; sec. 5.2.

[2] A. Varga, Least order fault and model detection using multi-models. IEEE CDC'09, Shanghai, China, 2009.