# diskmargin

Disk-based stability margins of feedback loops

## Syntax

``[DM,MM] = diskmargin(L)``
``MMIO = diskmargin(P,C)``
``___ = diskmargin(___,sigma)``

## Description

example

````[DM,MM] = diskmargin(L)` computes the disk-based stability margins for the SISO or MIMO negative feedback loop `feedback(L,eye(N))`, where `N` is the number of inputs and outputs in `L`. The `diskmargin` command returns loop-at-a-time stability margins in `DM` and multiloop margins in `MM`. Disk-based margin analysis provides a stronger guarantee of stability than the classical gain and phase margins. For general information about disk margins, see Stability Analysis Using Disk Margins.```

example

````MMIO = diskmargin(P,C)` computes the stability margins when considering independent, concurrent variations at both the plant inputs and plant outputs the negative feedback loop of the following diagram. ```

example

````___ = diskmargin(___,sigma)` specifies an additional skew parameter that biases the modeled gain and phase variation toward gain increase (positive `sigma`) or gain decrease (negative `sigma`). You can use this argument to test the relative sensitivity of stability margins to gain increases versus decreases. You can use this argument with any of the previous syntaxes.```

## Examples

collapse all

`diskmargin` computes both loop-at-a-time and multiloop disk margins. This example illustrates that loop-at-a-time margins can give an overly optimistic assessment of the true robustness of MIMO feedback loops. Margins of individual loops can be sensitive to small perturbations in other loops, and loop-at-a-time margins ignore such loop interactions.

Consider the two-channel MIMO feedback loop of the following illustration.

The plant model `P` is drawn from MIMO Stability Margins for Spinning Satellite and `C` is the static output-feedback gain [1 -2;0 1].

```a = [0 10;-10 0]; b = eye(2); c = [1 10;-10 1]; P = ss(a,b,c,0); C = [1 -2;0 1]; ```

Compute the disk-based margins at the plant output. The negative-feedback open-loop response at the plant output is `Lo = P*C`.

```Lo = P*C; [DMo,MMo] = diskmargin(Lo);```

Examine the loop-at-a-time disk margins returned in the structure array `DM`. Each entry in `DM` contains the stability margins of the corresponding feedback channel.

`DMo(1)`
```ans = struct with fields: GainMargin: [0 Inf] PhaseMargin: [-90 90] DiskMargin: 2 LowerBound: 2 UpperBound: 2 Frequency: Inf WorstPerturbation: [2x2 ss] ```
`DMo(2)`
```ans = struct with fields: GainMargin: [0 Inf] PhaseMargin: [-90 90] DiskMargin: 2 LowerBound: 2 UpperBound: 2 Frequency: 0 WorstPerturbation: [2x2 ss] ```

The loop-at-a-time margins are excellent (infinite gain margin and 90° phase margin). Next examine the multiloop disk margins `MMo`. These consider independent and concurrent gain (phase) variations in both feedback loops. This is a more realistic assessment because plant uncertainty typically affects both channels simultaneously.

`MMo`
```MMo = struct with fields: GainMargin: [0.6839 1.4621] PhaseMargin: [-21.2607 21.2607] DiskMargin: 0.3754 LowerBound: 0.3754 UpperBound: 0.3762 Frequency: 0 WorstPerturbation: [2x2 ss] ```

The multiloop gain and phase margins are much weaker than their loop-at-a-time counterparts. Stability is only guaranteed when the gain in each loop varies by a factor less than 1.46, or when the phase of each loop varies by less than 21°. Use `diskmarginplot` to visualize the gain and phase margins as a function of frequency.

`diskmarginplot(Lo)`

Typically, there is uncertainty in both the actuators (inputs) and sensors (outputs). Therefore, it is a good idea to compute the disk margins at the plant inputs as well as the outputs. Use `Li = C*P` to compute the margins at the plant inputs. For this system, the margins are the same at the plant inputs and outputs.

```Li = C*P; [DMi,MMi] = diskmargin(Li); MMi```
```MMi = struct with fields: GainMargin: [0.6839 1.4621] PhaseMargin: [-21.2607 21.2607] DiskMargin: 0.3754 LowerBound: 0.3754 UpperBound: 0.3762 Frequency: 0 WorstPerturbation: [2x2 ss] ```

Finally, you can also compute the multiloop disk margins for gain or phase variations at both the inputs and outputs of the plant. This approach is the most thorough assessment of stability margins, because it this considers independent and concurrent gain or phase variations in all input and output channels. As expected, of all three measures, this gives the smallest gain and phase margins.

```MMio = diskmargin(P,C); diskmarginplot(MMio.GainMargin)```

Stability is only guaranteed when the gain varies by a less than 2 dB or when the phase varies by less than 13°. However, these variations take place at the inputs and the outputs of P, so the total change in I/O gain or phase is twice that.

By default, `diskmargin` computes a symmetric gain margin, with `gmin = 1/gmax`, and an associated phase margin. In some systems, however, loop stability may be more sensitive to increases or decreases in open-loop gain. Use the skew parameter `sigma` to examine this sensitivity.

Compute the disk margin and associated disk-based gain and phase margins for a SISO transfer function, at three values of `sigma`. Negative `sigma` biases the computation toward gain decrease. Positive `sigma` biases toward gain increase.

```L = tf(25,[1 10 10 10]); DMdec = diskmargin(L,-2); DMbal = diskmargin(L,0); DMinc = diskmargin(L,2); DGMdec = DMdec.GainMargin```
```DGMdec = 1×2 0.4013 1.3745 ```
`DGMbal = DMbal.GainMargin`
```DGMbal = 1×2 0.6273 1.5942 ```
`DGMinc = DMinc.GainMargin `
```DGMinc = 1×2 0.7717 1.7247 ```

Put together, these results show that in the absence of phase variation, stability is maintained for relative gain variations between 0.4 and 1.72. To see how the phase margin depends on these gain variations, plot the stable ranges of gain and phase variations for each `diskmargin` result.

```diskmarginplot([DGMdec;DGMbal;DGMinc]) legend('sigma = -2','sigma = 0','sigma = 2') title('Stable range of gain and phase variations')```

This plot shows that the feedback loop can tolerate larger phase variations when the gain decreases. In other words, the loop stability is more sensitive to gain increase. Although `sigma` = –2 yields a phase margin as large as 30 degrees, this large value assumes a small gain increase of less than 3 dB. However, the plot shows that when the gain increases by 4 dB, the phase margin drops to less than 15 degrees. By contrast, it remains greater than 30 degrees when the gain decreases by 4 dB.

Thus, varying the skew `sigma` can give a fuller picture of sensitivity to gain and phase uncertainty. Unless you are mostly concerned with gain variations in one direction (increase or decrease), it is not recommended to draw conclusions from a single nonzero value of `sigma`. Instead use the default `sigma` = 0 to get unbiased estimates of gain and phase margins. When using nonzero values of `sigma`, use both positive and negative values to compare relative sensitivity to gain increase and decrease.

## Input Arguments

collapse all

Open-loop response, specified as a dynamic system model. `L` can be SISO or MIMO, as long as it has the same number of inputs and outputs. `diskmargin` computes the disk-based stability margins for the negative-feedback closed-loop system `feedback(L,eye(N))`.

To compute the disk margins of the positive feedback system `feedback(L,eye(N),+1)`, use `diskmargin(-L)`.

When you have a plant `P` and a controller `C`, you can compute the disk margins for gain (or phase) variations at the plant inputs or outputs, as in the following diagram.

• To compute margins at the plant outputs, set `L = P*C`.

• To compute margins at the plant inputs, set `L = C*P`.

`L` can be continuous time or discrete time. If `L` is a generalized state-space model (`genss` or `uss`) then `diskmargin` uses the current or nominal value of all control design blocks in `L`.

If `L` is a frequency-response data model (such as `frd`), then `diskmargin` computes the margins at each frequency represented in the model. The function returns the margins at the frequency with the smallest disk margin.

If `L` is a model array, then `diskmargin` computes margins for each model in the array.

Plant, specified as a dynamic system model. `P` can be SISO or MIMO, as long as `P*C` has the same number of inputs and outputs. `diskmargin` computes the disk-based stability margins for a negative-feedback closed-loop system. To compute the disk margins of the system with positive feedback, use `diskmargin(P,-C)`.

`P` can be continuous time or discrete time. If `P` is a generalized state-space model (`genss` or `uss`) then `diskmargin` uses the current or nominal value of all control design blocks in `P`.

If `P` is a frequency-response data model (such as `frd`), then `diskmargin` computes the margins at each frequency represented in the model. The function returns the margins at the frequency with the smallest disk margin.

Controller, specified as a dynamic system model. `C` can be SISO or MIMO, as long as `P*C` has the same number of inputs and outputs. `diskmargin` computes the disk-based stability margins for a negative-feedback closed-loop system. To compute the disk margins of the system with positive feedback, use `diskmargin(P,-C)`.

`C` can be continuous time or discrete time. If `C` is a generalized state-space model (`genss` or `uss`) then `diskmargin` uses the current or nominal value of all control design blocks in `C`.

If `C` is a frequency-response data model (such as `frd`), then `diskmargin` computes the margins at each frequency represented in the model. The function returns the margins at the frequency with the smallest disk margin.

Skew of uncertainty region used to compute the stability margins, specified as a real scalar value. This parameter biases the uncertainty used to model gain and phase variations toward gain increase or gain decrease.

• The default `sigma` = 0 uses a balanced model of gain variation in a range `[gmin,gmax]`, with ```gmin = 1/gmax```.

• Positive `sigma` uses a model with more gain increase than decrease (`gmax > 1/gmin`).

• Negative `sigma` uses a model with more gain decrease than increase (`gmin < 1/gmax`).

Use the default `sigma` = 0 to get unbiased estimates of gain and phase margins. You can test relative sensitivity to gain increase and decrease by comparing the margins obtained with both positive and negative `sigma` values. For an example, see Sensitivity of Disk-Based Margins to Gain Increase and Decrease. For more detailed information about how the choice of `sigma` affects the margin computation, see Stability Analysis Using Disk Margins.

## Output Arguments

collapse all

Disk margins for each feedback channel with all other loops closed, returned as a structure for SISO feedback loops, or an N-by-1 structure array for a MIMO loop with N feedback channels. The fields of `DM(i)` are:

FieldValue
`GainMargin`Disk-based gain margins of the corresponding feedback channel, returned as a vector of the form `[gmin,gmax]`. These values express in absolute units the amount by which the loop gain in that channel can decrease or increase while preserving stability. For example, if ```DM(i).GainMargin = [0.8,1.25]``` then the gain of the ith loop can be multiplied by any factor between 0.8 and 1.25 without causing instability. When `sigma` = 0, `gmin = 1/gmax`. If the open-loop gain can change sign without loss of stability, `gmin` can be less than zero for large enough negative `sigma`. If the nominal closed-loop system is unstable, then `DM(i).GainMargin = [1 1]`.
`PhaseMargin`Disk-based phase margin of the corresponding feedback channel, returned as a vector of the form `[-pm,pm]` in degrees. These values express the amount by which the loop phase in that channel can decrease or increase while preserving stability. If the closed-loop system is unstable, then `DM(i).PhaseMargin = [0 0]`.
`DiskMargin`Maximum ɑ compatible with closed-loop stability for the corresponding feedback channel. ɑ parameterizes the uncertainty in the loop response (see Algorithms). If the closed-loop system is unstable, then ```DM(i).DiskMargin = 0```.
`LowerBound`Lower bound on disk margin. This value is the same as `DiskMargin`.
`UpperBound`Upper bound on disk margin. This value represents an upper limit on the actual disk margin of the system. In other words, the disk margin is guaranteed to be no worse than `LowerBound` and no better than `UpperBound`.
`Frequency`Frequency at which the weakest margin occurs for the corresponding loop channel. This value is in rad/`TimeUnit`, where `TimeUnit` is the `TimeUnit` property of `L`.
`WorstPerturbation`

Smallest gain and phase variation that drives the feedback loop unstable, returned as a state-space (`ss`) model with N inputs and outputs, where N is the number of inputs and outputs in `L`. The system `F(s) = WorstPerturbation` is such that the following feedback loop is marginally stable, with a pole on the stability boundary at the frequency `DM(i).Frequency`.

This state-space model is a diagonal perturbation of the form `F(s) = diag(f1(s),...,fN(s))`. Each `fj(s)` is a real-parameter dynamic system that realizes the worst-case complex gain and phase variation applied to each channel of the feedback loop. For the loop-at-a-time margin of the kth feedback loop, only the kth entry `fk(s)` of `DM(k).WorstPerturbation` differs from unity.

For more information on interpreting `WorstPerturbation`, see Disk Margin and Smallest Destabilizing Perturbation

When analyzing a linear approximation of a nonlinear system, it can be useful to inject `WorstPerturbation` into the nonlinear simulation to further analyze the destabilizing affect of this worst-case gain and phase variation. For an example, see Robust MIMO Controller for Two-Loop Autopilot.

When `L = P*C` is the open-loop response of a system comprising a controller and plant with unit negative feedback in each channel, `DM` contains the stability margins for variations at the plant outputs. To compute the stability margins for variations at the plant inputs, use `L = C*P`. To compute the stability margins for simultaneous, independent variations at both the plant inputs and outputs, use ```MMIO = diskmargin(P,C)```.

When `L` is a model array, `DM` has additional dimensions corresponding to the array dimensions of `L`. For instance, if `L` is a 1-by-3 array of two-input, two-output models, then `DM` is a 2-by-3 structure array. `DM(j,k)` contains the margins for the jth feedback channel of the kth model in the array.

Multiloop disk margins, returned as a structure. The gain (or phase) margins quantify how much gain variation (or phase variation) the system can tolerate in all feedback channels at once while remaining stable. Thus, `MM` is a single structure regardless of the number of feedback channels in the system. (For SISO systems, `MM` = `DM`.) The fields of `MM` are:

FieldValue
`GainMargin`Multiloop disk-based gain margins, returned as a vector of the form `[gmin,gmax]`. These values express in absolute units the amount by which the loop gain can vary in all channels independently and concurrently while preserving stability. For example, if ```MM.GainMargin = [0.8,1.25]``` then the gain of all loops can be multiplied by any factor between 0.8 and 1.25 without causing instability. When `sigma` = 0, `gmin = 1/gmax`.
`PhaseMargin`Multiloop disk-based phase margin, returned as a vector of the form `[-pm,pm]` in degrees. These values express the amount by which the loop phase can vary in all channels independently and concurrently while preserving stability.
`DiskMargin`Maximum ɑ compatible with closed-loop stability. ɑ parameterizes the uncertainty in the loop response (see Algorithms).
`LowerBound`Lower bound on disk margin. This value is the same as `DiskMargin`.
`UpperBound`Upper bound on disk margin. This value represents an upper limit on the actual disk margin of the system. In other words, the disk margin is guaranteed to be no worse than `LowerBound` and no better than `UpperBound`.
`Frequency`Frequency at which the weakest margin occurs. This value is in rad/`TimeUnit`, where `TimeUnit` is the `TimeUnit` property of `L`.
`WorstPerturbation`

Smallest gain and phase variation that drives the feedback loop unstable, returned as a state-space (`ss`) model with N inputs and outputs, where N is the number of inputs and outputs in `L`. The system `F(s) = WorstPerturbation` is such that the following feedback loop is marginally stable, with a pole on the stability boundary at `MM.Frequency`.

This state-space model is a diagonal perturbation of the form `F(s) = diag(f1(s),...,fN(s))`. Each `fj(s)` is a real-parameter dynamic system that realizes the worst-case complex gain and phase variation applied to each channel of the feedback loop.

For more information on interpreting `WorstPerturbation`, see Disk Margin and Smallest Destabilizing Perturbation

When analyzing a linear approximation of a nonlinear system, it can be useful to inject `WorstPerturbation` into the nonlinear simulation to further analyze the destabilizing affect of this worst-case gain and phase variation. For an example, see Robust MIMO Controller for Two-Loop Autopilot.

When `L = P*C` is the open-loop response of a system comprising a controller and plant with unit negative feedback in each channel, `MM` contains the stability margins for variations at the plant outputs. To compute the stability margins for variations at the plant inputs, use `L = C*P`. To compute the stability margins for simultaneous, independent variations at both the plant inputs and outputs, use ```MMIO = diskmargin(P,C)```.

When `L` is a model array, `MM` is a structure array with one entry for each model in `L`.

Disk margins for independent variations applied simultaneously at input and output channels of the plant `P`, returned as a structure having the same fields as `MM`.

For variations applied simultaneously at inputs and outputs, the `WorstPerturbation` field is itself a structure with fields `Input` and `Output`. Each of these fields contains a state-space model such that for ```Fi(s) = MMIO.WorstPerturbation.Input``` and ```Fo(s) = MMIO.WorstPerturbation.Output```, the system of the following diagram is marginally unstable, with a pole on the stability boundary at the frequency `MMIO.Frequency`.

These state-space models `Input` and `Output` are diagonal perturbations of the form `F(s) = diag(f1(s),...,fN(s))`. Each `fj(s)` is a real-parameter dynamic system that realizes the worst-case complex gain and phase variation applied to each channel of the feedback loop.

## Tips

• `diskmargin` assumes negative feedback. To compute the disk margins of a positive feedback system, use `diskmargin(-L)` or `diskmargin(P,-C)`.

• To compute disk margins for a system modeled in Simulink®, first linearize the model to obtain the open-loop response at a particular operating point. Then, use `diskmargin` to compute stability margins for the linearized system. For more information, see Stability Margins of a Simulink Model.

• To compute classical gain and phase margins, use `allmargin`.

• You can visualize disk margins using `diskmarginplot`.

## Algorithms

collapse all

`diskmargin` computes gain and phase margins by applying a disk-based uncertainty model to represent gain and phase variations, and then finding the largest such disk for which the closed-loop system is stable.

### Gain and Phase Uncertainty Model

For SISO L, the uncertainty model for disk-margin analysis incorporates a multiplicative complex uncertainty F into the loop transfer function as follows:

`$F=\frac{1+\alpha \left[\left(1-\sigma \right)/2\right]\delta }{1-\alpha \left[\left(1+\sigma \right)/2\right]\delta }.$`

Here,

• δ is a gain-bounded dynamic uncertainty, normalized so that it always varies within the unit disk (|δ| < 1).

• α sets the amount of gain and phase variation modeled by F. For fixed σ, the parameter ɑ controls the size of the disk. For α = 0, the multiplicative factor is 1, corresponding to the nominal L.

• σ, called the skew, biases the modeled uncertainty toward gain increase or gain decrease. (For details about the effect of skew on the uncertainty model, see Stability Analysis Using Disk Margins.)

For MIMO systems, the model allows the uncertainty to vary independently in each channel:

`${F}_{j}=\frac{1+\alpha \left[\left(1-\sigma \right)/2\right]{\delta }_{j}}{1-\alpha \left[\left(1+\sigma \right)/2\right]{\delta }_{j}}.$`

The model replaces the MIMO open-loop response L with L*F, where

`$F=\left(\begin{array}{ccc}{F}_{1}& 0& 0\\ 0& \ddots & 0\\ 0& 0& {F}_{N}\end{array}\right).$`

### Disk-Margin Computation

For a given `sigma`, the disk margin is the largest ɑ for which the closed-loop system `feedback(L*F,1)` (or `feedback(L*F,eye(N))` for MIMO systems) is stable for all values of F. To find this value, `diskmargin` solves a robust stability problem: Find the largest α such that the closed-loop system is stable for all F in the uncertainty disk Δ(α,σ) described by

`$\Delta \left(\alpha ,\sigma \right)=\left\{F=\frac{1+\alpha \left[\left(1-\sigma \right)/2\right]\delta }{1-\alpha \left[\left(1+\sigma \right)/2\right]\delta }\text{\hspace{0.17em}}\text{\hspace{0.17em}}:\text{\hspace{0.17em}}\text{\hspace{0.17em}}|\delta |<1\right\}.$`

In the SISO case, the robust stability analysis leads to

`${\alpha }_{max}={‖\frac{1}{S+\left(\sigma -1\right)/2}‖}_{\infty },$`

where S is the sensitivity function (1 + L)–1 .

In the MIMO case, the robust stability analysis leads to

`${\alpha }_{max}=\frac{1}{{\mu }_{\Delta }\left(S+\frac{\left(\sigma -1\right)I}{2}\right)}.$`

Here, μΔ is the structured singular value (`mussv`) for the diagonal structure

`$\Delta =\left(\begin{array}{ccc}{\delta }_{1}& 0& 0\\ 0& \ddots & 0\\ 0& 0& {\delta }_{N}\end{array}\right),$`

and δj is the normalized uncertainty for each Fj.

For more details about the margin computation, see [2].

## References

[1] Blight, James D., R. Lane Dailey, and Dagfinn Gangsaas. “Practical Control Law Design for Aircraft Using Multivariable Techniques.” International Journal of Control 59, no. 1 (January 1994): 93–137. https://doi.org/10.1080/00207179408923071.

[2] Seiler, Peter, Andrew Packard, and Pascal Gahinet. “An Introduction to Disk Margins [Lecture Notes].” IEEE Control Systems Magazine 40, no. 5 (October 2020): 78–95.

## Version History

Introduced in R2018b

expand all

Behavior changed in R2020a