# Solution of model-matching problems

Generalized Nehari approximation.`gnehari`

Solution of the least distance problem.`glinfldp`

Approximate solution of the linear rational matrix equation`grasol`

`G(λ)*X(λ) = F(λ)`

.Approximate solution of the linear rational matrix equation`glasol`

`X(λ)*G(λ) = F(λ)`

.

`DescriptorSystems.gnehari`

— Function```
gnehari(sys[, γ]; fast = true, offset = β,
atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, σ1)
```

Compute for the descriptor system `sys = (A-λE,B,C,D)`

with the transfer function matrix `G(λ)`

the optimal or suboptimal Nehari approximation `sysx = (Ax-λEx,Bx,Cx,Dx)`

with the transfer function matrix `X(λ)`

. The optimal Nehari approximation `X(λ)`

satisfies

\[ \| G(\lambda) - X(\lambda) \|_\infty = \| G^{*}_u(\lambda) \|_H := \sigma_1,\]

where ${\small G_u(\lambda)}$ is the proper antistable part of `G(λ)`

. The resulting $σ_1$ is the Hankel-norm of ${\small G^{*}_u(\lambda)}$ (also the L∞-norm of the optimal approximation error). For a given $γ > σ_1$, the suboptimal approximation satisfies

\[ \| G(\lambda) - X(\lambda) \|_\infty \leq \gamma .\]

The Nehari approximation is stable if `G(λ)`

has no poles on the boundary of the stability domain `Cs`

. In the continuous-time, `Cs`

is the set of complex numbers with negative real parts and its boundary is the extended imaginary axis (containig also the poit at infinity), while, in the discrete-time case, `Cs`

is the set complex number with moduli less than one and its boundary is the unit circle centered in the origin.

To assess the presence of poles in `Cs`

, a boundary offset `β`

can be specified via the keyword parameter `offset = β`

. Accordingly, for a continuous-time system, the stable poles in `Cs`

have real parts less than or equal to `β`

, while for a discrete-time system, they have moduli less than or equal to 1+β`. The default value used for`

β`is`

sqrt(ϵ)`, where`

ϵ`is the working machine precision. For a negative values of`

β`(e.g.,`

β = -sqrt(ϵ)`), an extended stability domain corresponding to the closure of`

Cs`is used instead of`

Cs`.

The keyword arguments `atol1`

, `atol2`

, and `rtol`

, specify, respectively, the absolute tolerance for the nonzero elements of `A`

, the absolute tolerance for the nonzero elements of `E`

, and the relative tolerance for the nonzero elements of `A`

and `E`

. The default relative tolerance is `n*ϵ`

, where `ϵ`

is the working machine epsilon and `n`

is the order of the system `sys`

. The keyword argument `atol`

can be used to simultaneously set `atol1 = atol`

and `atol2 = atol`

.

The separation of the finite and infinite eigenvalues is performed using rank decisions based on rank revealing QR-decompositions with column pivoting if `fast = true`

or the more reliable SVD-decompositions if `fast = false`

.

*Method:* Let `G(λ) = Gs(λ) + Gu(λ)`

be an additive decomposition of `G(λ)`

such that `Gs(λ)`

has only poles in the closure of `Cs`

(or its closure if `β < 0`

) and `Gu(λ)`

is the antistable part having only unstable poles. The Hankel-norm approximation methods of [1] and [2], with extensions for descriptor systems, are used to compute a stable Nehari approximation `Gn(λ)`

of the unstable part `Gu(λ)`

and the resulting solution is computed as `X(λ) = Gs(λ) + Gn(λ)`

.

*References:*

[1] K. Glover. All optimal Hankel-norm approximations of linear multivariable systems and their L∞ error bounds, Int. J. Control, vol. 39, pp. 1115-1193, 1984.

[2] M. G. Safonov, R. Y. Chiang, and D. J. N. Limebeer. Optimal Hankel model reduction for nonminimal systems. IEEE Trans. Automat. Control, vol. 35, pp. 496–502, 1990.

`DescriptorSystems.glinfldp`

— Function```
glinfldp(sys1, sys2, [, γ]; nehari = false, reltol = 0.0001, fast = true, offset = β,
atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, mindist)
```

Determine for the descriptor systems `sys1 = (A1-λE1,B1,C1,D1)`

and `sys2 = (A2-λE2,B2,C2,D2)`

with the transfer function matrices $G_1(λ)$ and $G_2(λ)$, respectively, the descriptor system `sysx`

with the transfer function matrix $X(λ)$ such that $X(λ)$ is the solution of the 2-block `L∞`

*least distance problem* (`LDP`

)

\[ \text{mindist} := \min \|G_1(λ)-X(λ) \mid G_2(λ) \|_\infty \]

`mindist`

is the achieved minimum distance corresponding to the optimal solution. If `sys2 = []`

, an 1-block `LDP`

is solved. `sys1`

and `sys2`

must not have poles on the boundary of the stability domain `Cs`

. It is assumed that `sys2`

has no poles on the boundary of the stability domain (see below). The resulting solution is stable, provided `sys1`

has no poles on the boundary of the stability domain as well.

If ${\small γ > \|G_2(λ)\|_\infty}$ is a desired sub-optimality degree, then the `γ`

-suboptimal `LDP`

\[ \text{mindist} := \|G_1(λ)-X(λ) \mid G_2(λ) \|_\infty < γ\]

is solved and `mindist`

is the achieved suboptimal distance.

The call with

```
glinfldp(sys[, m2[, γ]]; nehari = false, fast = true, offset = β,
atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, mindist)
```

uses the compound descriptor system `sys = (A-λE,[B1 B2],C,[D1 D2])`

, where `B2`

has `m2`

columns, to define the descriptor systems `sys1 = (A-λE,B1,C,D1)`

and `sys2 = (A-λE,B2,C,D2)`

(i.e., `A1-λE1 = A2-λE2 = A-λE`

and `C1 = C2 = C`

). If `m2 = 0`

, an 1-block `LDP`

is solved. `sys2`

must not have poles on the boundary of the stability domain `Cs`

.

If `nehari = true`

, the optimal or suboptimal Nehari approximation is used to solve the `LDP`

. If `nehari = false`

(default), the optimal solution is computed using the `γ`

-iteration [1].

The keyword argument `reltol`

specifies the relative tolerance for the desired accuracy of `γ`

-iteration. The iterations are performed until the current estimations of maximum $γ_u$ and minimum $γ_l$ of the optimal distance $γ_o$, $γ_l \leq γ_o \leq γ_u$, satisfies

\[ γ_u-γ_l \leq \text{reltol} * \text{gap} ,\]

where `gap`

is the original gap (internally determined).

To assess the presence of poles on the boundary of the stability domain `Cs`

, a boundary offset `β`

can be specified via the keyword parameter `offset = β`

. Accordingly, for a continuous-time setting, the boundary of `Cs`

contains the complex numbers with real parts within the interval `[-β,β]`

, while for a discete-time setting, the boundary of `Cs`

contains the complex numbers with moduli within the interval `[1-β,1+β]`

. The default value used for `β`

is `sqrt(ϵ)`

, where `ϵ`

is the working machine precision.

The keyword arguments `atol1`

, `atol2`

, and `rtol`

, specify, respectively, the absolute tolerance for the nonzero elements of `A`

, the absolute tolerance for the nonzero elements of `E`

, and the relative tolerance for the nonzero elements of `A`

and `E`

. The default relative tolerance is `n*ϵ`

, where `ϵ`

is the working machine epsilon and `n`

is the order of the system `sys`

. The keyword argument `atol`

can be used to simultaneously set `atol1 = atol`

and `atol2 = atol`

.

The rank decisions in the underlying pencil manipulation algorithms are based on rank revealing QR-decompositions with column pivoting if `fast = true`

or the more reliable SVD-decompositions if `fast = false`

.

*Method:* The approach of [1] is used for the solution of the 2-block least distance problem. Necessary conditions for solvability of the LDP with a stable solution is that `sys1`

and `sys2`

have no poles on the boundary of the stability domain `Cs`

.

*References:* [1] C.-C. Chu, J. C. Doyle, and E. B. Lee The general distance problem in H∞ optimal control theory, Int. J. Control, vol 44, pp. 565-596, 1986.

`DescriptorSystems.grasol`

— Function```
grasol(sysg, sysf[, γ]; L2sol = false, nehari = false, reltol = 0.0001, mindeg = false, poles, sdeg,
fast = true, offset = β, atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, info)
```

Determine for the descriptor systems `sysg = (Ag-λEg,Bg,Cg,Dg)`

and `sysf = (Af-λEf,Bf,Cf,Df)`

with the transfer function matrices `G(λ)`

and `F(λ)`

, respectively, the descriptor system `sysx`

with the transfer function matrix `X(λ)`

such that `X(λ)`

is the approximate solution of the linear rational equation `G(λ)X(λ) = F(λ)`

, which achieves the minimum error norm ${\small \text{mindist} := \min \|G(λ)X(λ) - F(λ)\|}$. The resulting `X(λ)`

has all poles stable or lying on the boundary of the stability domain `Cs`

. If `L2sol = false`

(default) then the `L∞`

-norm optimal solution is computed, while if `L2sol = true`

the `L2`

-norm optimal solution is computed. `sysg`

and `sysf`

must not have poles on the boundary of the stability domain `Cs`

.

If `γ > 0`

is a desired sub-optimality degree, then the `γ`

-suboptimal model-matching problem

\[ \text{mindist} := \|G(λ)X(λ) - F(λ) \| < γ\]

is solved and `mindist`

is the achieved suboptimal distance.

The call with

```
grasol(sysgf[, mf[, γ]]; L2sol = false, nehari = false, reltol = 0.0001, mindeg = false, poles, sdeg,
fast = true, offset = β, atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, info)
```

uses the compound descriptor system `sysgf = (A-λE,[Bg Bf],C,[Dg Df])`

, where `Bf`

and `Df`

have `mf`

columns, to define the descriptor systems `sysg = (A-λE,Bg,C,Dg)`

and `sysf = (A-λE,Bf,C,Df)`

(i.e., `Ag-λEg = Af-λEf = A-λE`

and `Cg = Cf = C`

). `sysgf`

must not have poles on the boundary of the stability domain `Cs`

.

If `nehari = true`

, the optimal or suboptimal Nehari approximation is used to compute a `L∞`

-suboptimal solution of the underlying *least-distance problem* (`LDP`

). If `nehari = false`

(default), the `L∞`

-optimal solution is computed using the `γ`

-iteration in the underlying `LDP`

[2].

If `mindeg = true`

, a minimum order solution is determined (if possible), while if `mindeg = false`

(default) a particular solution of non-minimal order is determined.

The resulting named tuple `info`

contains additional information: `info.nrank`

is the normal rank of `G(λ)`

, `info.nr`

is the number of freely assignable poles of the solution `X(λ)`

, `info.mindist`

is the achieved approximation error norm and `info.nonstandard`

is `true`

for a non-standard problem, with `G(λ)`

having zeros on the boundary of the stability domain, and `false`

for a standard problem, when `G(λ)`

has no zeros on the boundary of the stability domain.

The keyword argument `reltol`

specifies the relative tolerance for the desired accuracy of the `γ`

-iteration employed to solve the underlying least-distance problem. The iterations are performed until the current estimations of maximum $γ_u$ and minimum $γ_l$ of the optimal distance satisfies ${\small γ_u-γ_l \leq \text{reltol} * \text{gap}}$, where `gap`

is the initial estimation of the error gap.

To assess the presence of poles on the boundary of the stability domain `Cs`

, a boundary offset `β`

can be specified via the keyword parameter `offset = β`

. Accordingly, for a continuous-time setting, the boundary of `Cs`

contains the complex numbers with real parts within the interval `[-β,β]`

, while for a discrete-time setting, the boundary of `Cs`

contains the complex numbers with moduli within the interval `[1-β,1+β]`

. The default value used for `β`

is `sqrt(ϵ)`

, where `ϵ`

is the working machine precision.

The vector `poles`

specified as a keyword argument, can be used to specify the desired poles of `sysx`

alternatively to or jointly with enforcing a desired stability degree `sdeg`

of poles.

The keyword arguments `atol1`

, `atol2`

, and `rtol`

, specify, respectively, the absolute tolerance for the nonzero elements of `Ag`

, `Af`

, `A`

, `Bg`

, `Bf`

, `Cg`

, `Cf`

, `Dg`

, `Df`

, the absolute tolerance for the nonzero elements of `Eg`

, `Ef`

, and the relative tolerance for the nonzero elements of all above matrices. The default relative tolerance is `n*ϵ`

, where `ϵ`

is the working machine epsilon and `n`

is the order of the system `sysg`

. The keyword argument `atol`

can be used to simultaneously set `atol1 = atol`

and `atol2 = atol`

.

The rank decisions in the underlying pencil manipulation algorithms are based on rank revealing QR-decompositions with column pivoting if `fast = true`

or the more reliable SVD-decompositions if `fast = false`

.

*Method:* An extension of the approach of [1] to descriptor systems is used.

*References:*

[1] B. A. Francis. A Course in H-infinity Theory, Vol. 88 of Lecture Notes in Control and Information Sciences, Springer-Verlag, New York, 1987.

[2] C.-C. Chu, J. C. Doyle, and E. B. Lee. The general distance problem in H∞ optimal control theory, Int. J. Control, vol 44, pp. 565-596, 1986.

`DescriptorSystems.glasol`

— Function```
glasol(sysg, sysf[, γ]; L2sol = false, nehari = false, reltol = 0.0001, mindeg = false, poles, sdeg,
fast = true, offset = β, atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, info)
```

Determine for the descriptor systems `sysg = (Ag-λEg,Bg,Cg,Dg)`

and `sysf = (Af-λEf,Bf,Cf,Df)`

with the transfer function matrices `G(λ)`

and `F(λ)`

, respectively, the descriptor system `sysx`

with the transfer function matrix `X(λ)`

such that `X(λ)`

is the approximate solution of the linear rational equation `X(λ)G(λ) = F(λ)`

, which achieves the minimum error norm ${\small \text{mindist} := \min \|X(λ)G(λ) - F(λ)\|}$. The resulting `X(λ)`

has all poles stable or lying on the boundary of the stability domain `Cs`

. If `L2sol = false`

(default) then the `L∞`

-norm optimal solution is computed, while if `L2sol = true`

the `L2`

-norm optimal solution is computed. `sysg`

and `sysf`

must not have poles on the boundary of the stability domain `Cs`

.

If `γ > 0`

is a desired sub-optimality degree, then the `γ`

-suboptimal model-matching problem

\[ \text{mindist} := \|X(λ)G(λ) - F(λ) \| < γ\]

is solved and `mindist`

is the achieved suboptimal distance.

The call with

```
glasol(sysgf[, pf[, γ]]; L2sol = false, nehari = false, reltol = 0.0001, mindeg = false, poles, sdeg,
fast = true, offset = β, atol = 0, atol1 = atol, atol2 = atol, rtol = n*ϵ) -> (sysx, info)
```

uses the compound descriptor system `sysgf = (A-λE,B,[Cg; Cf],[Dg; Df])`

, where `Cf`

and `Df`

have `pf`

rows, to define the descriptor systems `sysg = (A-λE,B,Cg,Dg)`

and `sysf = (A-λE,B,Cf,Df)`

(i.e., `Ag-λEg = Af-λEf = A-λE`

and `Bg = Bf = B`

). `sysgf`

must not have poles on the boundary of the stability domain `Cs`

.

If `nehari = true`

, the optimal or suboptimal Nehari approximation is used to solve the underlying *least-distance problem* (`LDP`

). If `nehari = false`

(default), the optimal solution is computed using the `γ`

-iteration in the underlying `LDP`

[2].

If `mindeg = true`

, a minimum order solution is determined (if possible), while if `mindeg = false`

(default) a particular solution of non-minimal order is determined.

The resulting named tuple `info`

contains additional information: `info.nrank`

is the normal rank of `G(λ)`

, `info.nl`

is the number of freely assignable poles of the solution `X(λ)`

, `info.mindist`

is the achieved approximation error norm and `info.nonstandard`

is `true`

for a non-standard problem, with `G(λ)`

having zeros on the boundary of the stability domain, and `false`

for a standard problem, when `G(λ)`

has no zeros on the boundary of the stability domain.

The keyword argument `reltol`

specifies the relative tolerance for the desired accuracy of the `γ`

-iteration employed to solve the underlying least-distance problem. The iterations are performed until the current estimations of maximum $γ_u$ and minimum $γ_l$ of the optimal distance satisfies ${\small γ_u-γ_l < \text{reltol}* \text{gap}}$, where `gap`

is the initial estimation of the error gap.

To assess the presence of poles on the boundary of the stability domain `Cs`

, a boundary offset `β`

can be specified via the keyword parameter `offset = β`

. Accordingly, for a continuous-time setting, the boundary of `Cs`

contains the complex numbers with real parts within the interval `[-β,β]`

, while for a discrete-time setting, the boundary of `Cs`

contains the complex numbers with moduli within the interval `[1-β,1+β]`

. The default value used for `β`

is `sqrt(ϵ)`

, where `ϵ`

is the working machine precision.

The vector `poles`

specified as a keyword argument, can be used to specify the desired poles of `sysx`

alternatively to or jointly with enforcing a desired stability degree `sdeg`

of poles.

The keyword arguments `atol1`

, `atol2`

, and `rtol`

, specify, respectively, the absolute tolerance for the nonzero elements of `Ag`

, `Af`

, `A`

, `Bg`

, `Bf`

, `Cg`

, `Cf`

, `Dg`

, `Df`

, the absolute tolerance for the nonzero elements of `Eg`

, `Ef`

, and the relative tolerance for the nonzero elements of all above matrices. The default relative tolerance is `n*ϵ`

, where `ϵ`

is the working machine epsilon and `n`

is the order of the system `sysg`

. The keyword argument `atol`

can be used to simultaneously set `atol1 = atol`

and `atol2 = atol`

.

The rank decisions in the underlying pencil manipulation algorithms are based on rank revealing QR-decompositions with column pivoting if `fast = true`

or the more reliable SVD-decompositions if `fast = false`

.

*Method:* An extension of the approach of [1] to descriptor systems is used.

*References:*

[1] B. A. Francis. A Course in H-infinity Theory, Vol. 88 of Lecture Notes in Control and Information Sciences, Springer-Verlag, New York, 1987.

[2] C.-C. Chu, J. C. Doyle, and E. B. Lee. The general distance problem in H∞ optimal control theory, Int. J. Control, vol 44, pp. 565-596, 1986.