▲
                           │uw      │dw              │do
                         ┌─┴─┐    ┌─▼─┐            ┌─▼─┐
                         │   │    │   │            │   │
                         │ Wu│    │ Wi│            │ Wd│
                         │   │    │   │            │   │
    ┌────┐       ┌─────┐ └─▲─┘    └───┘   ┌─────┐  └─┬─┘    ┌───┐
 r  │    │       │     │   │      di▼     │     │   d│  y   │   │ e
 ──►│ Wr ├──►+──►│  K  ├───┴────────+────►│  P  ├────▼───┬─►│ We├──►
    │    │   ▲   │     │   u              │     │ Py     │  │   │
    └────┘   │-  └─────┘                  └─────┘        │  └───┘
             │                                           │
             │                                           │  ┌───┐
             │                                           ▼  │   │
             └───────────────────────────────────────────+◄─┤ Wn│◄─
                                                          n │   │ ni
                                                            └───┘

In this example we will show how to manually construct the following transfer function, the $H_\infty$ norm of which is commonly a performance metric for $H_\infty$ synthesis.

\[\begin{bmatrix} W_e S_o W_d \\ W_u K S_o W_d \end{bmatrix}\]

(Note, the function hinfpartition exists to make special instances of this kind of partitioning smoother, this example shows how general transfer matrices can be constructed for more complicated $H_\infty$ designs)

We note that the first transfer function is the transfer function from $d_o \rightarrow e$ and the second transfer function is from $d_o \rightarrow u_w$. It's thus enough for this example to define the blocks $P,W_e,W_d,W_u$. To form the feedback interconnection, we also need to define a summation node after $P$. Note, $K$ is the output of the synthesis procedure, we're thus going to break the loop at $K$ and place the inputs and outputs to $K$ ($u, y$) in the system_mapping, while the external inputs and outputs ($w, z$) go in the performance_mapping.

Using the connect function, we can then form the closed-loop transfer function automatically, without forming the intermediate output sensitivity function $S_o$.

The code below constructs all transfer functions using the same names as in the block diagram above. Before we start, it really helps to have a block diagram in front of you to not get lost. We redraw the diagram focusing on the relevant signals and blocks only

               ▲
               │uw              │do
             ┌─┴─┐            ┌─▼─┐
             │   │            │   │
             │ Wu│            │ Wd│
             │   │            │   │
     ┌─────┐ └─▲─┘   ┌─────┐  └─┬─┘    ┌───┐
     │     │   │     │     │   d│  y   │   │ e
┌───►│  K  ├───┴─────┤  P  ├────▼───┬─►│ We├──►
│    │     │   u     │     │        │  │   │
│    └─────┘         └─────┘        │  └───┘
│                                   │
│                                   │
│                                   │
└───────────────────────────────────┘
using RobustAndOptimalControl, ControlSystemsBase, Plots
P = named_ss(ss(tf(1, [1, 0.02, 1])), :P)
We = named_ss(10P.sys,  :We, u=:y, y=:e)  # Weigh sensitivity by magnitude of plant
Wu = named_ss(makeweight(1e-3, 10, 10), :Wu, u=:Wu, y=:uw) # Above ω=100, we want to limit the control effort
Wd = named_ss(makeweight(10, 1, 1e-3),  :Wd, u=:do, y=:d) # d is low frequency
Wd: NamedStateSpace{Continuous, Float64}
A = 
 -0.01
B = 
 0.25
C = 
 0.39996000000000004
D = 
 0.001

Continuous-time state-space model
With state  names: Wdx
     input  names: do
     output names: d

We also create the summation node after $P$

sumP = sumblock("y = Py + d")
sumblock: NamedStateSpace{Continuous, Float64}
D = 
 1.0  1.0

Continuous-time state-space model
With state  names: 
     input  names: Py d
     output names: y

We also need a splitter to split the external input :u (remember, the loop was broken at $K$) into two signals with unique names.

split_u = splitter(:u, 2)
splitter: NamedStateSpace{Continuous, Float64}
D = 
 1.0
 1.0

Continuous-time state-space model
With state  names: 
     input  names: u
     output names: u1 u2

The next step is to list all the connections

connections = [
    :u1 => :Wu # splitter to input of Wu
    :u2 => :Pu # splitter to input of P
    :Py => :Py # P output to first input of sumblock
    :d => :d   # output of Wd to second input of sumblock
    :y => :y   # output of sumblock to input of We
]

When we specify the external inputs and outputs to connect, we include $y$ and $u$ since they are external from the view of connect:

w1 = [ # External inputs
    :do, :u
]
z1 = [ # External outputs
    :e, :uw, :y
]

We are now ready to form the system we want to minimize the norm of

G = connect([P,We,Wu,Wd,sumP,split_u], connections; z1, w1)
NamedStateSpace{Continuous, Float64}
A = 
  0.0   1.0    0.0   0.0     0.0                0.0
 -1.0  -0.02   0.0   0.0     0.0                0.0
  0.0   0.0    0.0   1.0     0.0                0.0
  1.0   0.0   -1.0  -0.02    0.0                0.39996000000000004
  0.0   0.0    0.0   0.0   -99.49879346007116   0.0
  0.0   0.0    0.0   0.0     0.0               -0.01
B = 
 0.0     0.0
 0.0     1.0
 0.0     0.0
 0.001   0.0
 0.0    32.0
 0.25    0.0
C = 
 0.0  0.0  10.0  0.0    0.0               0.0
 0.0  0.0   0.0  0.0  -31.09026361897661  0.0
 1.0  0.0   0.0  0.0    0.0               0.39996000000000004
D = 
 0.0     0.0
 0.0    10.0
 0.001   0.0

Continuous-time state-space model
With state  names: Px1##feedback#477 Px2##feedback#477 Wex1##feedback#477 Wex2##feedback#477 Wux##feedback#477 Wdx##feedback#477
     input  names: do u
     output names: e uw y

At this stage, it's good practice to check that the poles, inputs and outputs of $G$ looks correct, it's easy to forget some signal..

Before we perform the $H_\infty$ synthesis, we need to partition this system to indicate which outputs go to and from the controller and which are part of the performance specification

Gsyn = partition(G, u = [:u], y = [:y]) # You can provide either u or w, and either y or z
ExtendedStateSpace{Continuous, Float64}
A = 
  0.0   1.0    0.0   0.0     0.0                0.0
 -1.0  -0.02   0.0   0.0     0.0                0.0
  0.0   0.0    0.0   1.0     0.0                0.0
  1.0   0.0   -1.0  -0.02    0.0                0.39996000000000004
  0.0   0.0    0.0   0.0   -99.49879346007116   0.0
  0.0   0.0    0.0   0.0     0.0               -0.01
B1 = 
 0.0
 0.0
 0.0
 0.001
 0.0
 0.25
B2 = 
  0.0
  1.0
  0.0
  0.0
 32.0
  0.0
C1 = 
 0.0  0.0  10.0  0.0    0.0               0.0
 0.0  0.0   0.0  0.0  -31.09026361897661  0.0
C2 = 
 1.0  0.0  0.0  0.0  0.0  0.39996000000000004
D11 = 
 0.0
 0.0
D12 = 
  0.0
 10.0
D21 = 
 0.001
D22 = 
 0.0
Continuous-time extended state-space model

Let's perform some checks to see that we have a sound problem specification

hinfassumptions(Gsyn)
true

we also check the equivalences we know should hold

system_mapping(Gsyn) == P.sys
true
G[:uw, :u].sys == Wu.sys
false

These will not be identical, the realization might differ, but they should represent the same system

hinfnorm2(G[:e, :do].sys - (We * Wd).sys)[1] < 1e-8
true
hinfnorm2(G[:e, :u].sys - (We * P).sys)[1] < 1e-8
false

Now, let's synthesize!

K, γ = hinfsynthesize(Gsyn, γrel=1.1, transform=false)
γ, K
(0.1972551401220788, StateSpace{Continuous, Float64}
A = 
    0.0                  1.0                    0.0                   0.0                  0.0                   0.0
   -5.029564643626379   -1.4683706824648812    -3.344193594283276    -6.880021065211896    2.283886347486903    -1.7120108840893369
    0.0                  0.0                    0.0                   1.0                  0.0                   0.0
    0.0                  0.0                   -1.0                  -0.02                 0.0                   0.0
 -128.9460685960441    -46.3478618388762     -107.01419501706484   -220.16067408678066   -26.414430340490256   -54.78434829085878
 -250.00000000000006     0.0                    0.0                   0.0                  0.0                -100.00000000000004
B = 
   0.0
   0.0
   0.0
   1.0
   0.0
 250.00000000000006
C = 
 -4.029564643626379  -1.4483706824648812  -3.344193594283276  -6.880021065211896  2.283886347486903  -1.7120108840893369
D = 
 0.0

Continuous-time state-space model)

We may now form the closed-loop system and check the calculated $H_\infty$ norm

Gcl = lft(Gsyn, K)
hinfnorm2(Gcl) # Should be equal to γ
(0.1972551401220788, 2.557529082231186)

The calculated $\gamma$ is quite small, indicating that our performance specifications were easy to satisfy. Let's plot the resulting performance mapping with the lower loop closed around $K$ as well as some interesting sensitivity functions

w = exp10.(LinRange(-3, 3, 300))
bodeplot(Gcl, w, plotphase=false, title=["do -> e" "do -> uw"])
Example block output
S, PS, CS, T = gangoffour(P, K)
specificationplot([S, CS, T], ylims=(1e-3, 1e2))
Example block output