Model Uncertainty and Robustness

Learning objectives

Apply a controller designed for one model to several actual models and delay conditions to find where performance deteriorates.

Applying the nominal controller to other actual models

Starting conditions

$m\in\{.7,1,1.3\}$ kg, $b=3$ N·s/m, $k=2$ N/m. Fix K and N at the values designed for $m=1$.

Procedure

1. Keep the nominal controller fixed and compare the poles and responses of the three actual models.
2. Find the final value and 2% settling time of each response, and identify the slowest condition.

Results to retain

• Three response curves and a comparison table of poles and settling times.

import numpy as np
from scipy.signal import place_poles,StateSpace,lsim
import matplotlib.pyplot as plt
A0=np.array([[0.,1.],[-2.,-3.]])
B0=np.array([[0.],[1.]]);C=np.array([[1.,0.]])
K=place_poles(A0,B0,[-2.,-4.]).gain_matrix;N=8.
t=np.linspace(0,6,1201)
for m in [.7,1.,1.3]:
    A=np.array([[0.,1.],[-2/m,-3/m]]);B=np.array([[0.],[1/m]])
    poles=np.linalg.eigvals(A-B@K)
    _,y,_=lsim(StateSpace(A-B@K,B*N,C,0),np.ones(len(t)),t)
    print('mass',m,'poles',poles);plt.plot(t,y,label=f'm={m}')
plt.xlabel('Time [s]');plt.ylabel('Position');plt.legend();plt.grid()

Record results and self-check

Mass:

Pole:

Final value:

Settling Time:

Points to check

• Did you keep K unchanged rather than redesigning it for each actual model?
• Did you check that the response stays within the 2% band after the settling time?

Changing mass, friction, and stiffness together

Starting conditions

Use the 27 combinations of m∈{.7,1,1.3}, b∈{2,3,4}, and k∈{1,2,3}. Keep $K=\begin{bmatrix}6&3\end{bmatrix}$ and $N=8$ fixed.

Procedure

1. For each combination, find the maximum real part of the eigenvalues of $A-BK$. Record the number of models where this value is negative and where it is 0 or greater.
2. Find the combination with the largest maximum real part. Also calculate $y_{ss}=-C(A-BK)^{-1}BN$ for each combination, and separately find the combination with the largest $|1-y_{ss}|$.
3. For the two selected combinations, tabulate $m,b,k$, maximum pole real part, and steady-state error, and plot unit-step responses on the same time axis.
4. Distinguish the parameter values checked in this calculation from the unchecked values between them. Explain which coefficient conditions in the learning material's characteristic polynomial should be examined to assess stability over the entire range.

Results to retain

• A stability summary for the 27 conditions and responses for the two extreme cases.

from itertools import product
records=[]
for m,b,k in product([.7,1.,1.3],[2.,3.,4.],[1.,2.,3.]):
    A=np.array([[0.,1.],[-k/m,-b/m]]);B=np.array([[0.],[1/m]])
    aa=A-B@K;pole=np.max(np.linalg.eigvals(aa).real)
    final=float((-C@np.linalg.solve(aa,B*N)).item())
    records.append((m,b,k,pole,abs(1-final)))
print('closest pole:',max(records,key=lambda r:r[3]))
print('largest error:',max(records,key=lambda r:r[4]))
# TODO: Add the step responses of both models on the same axes.

Record results and self-check

m,b,k:

Maximum pole real part:

Steady-State Error:

Points to check

• Did you calculate poles from $A-BK$ rather than A?
• Did you distinguish pole calculations for the 27 combinations from analytical justification covering the entire parameter range?

With input delay and saturation

Starting conditions

Euler-integrate the nominal model's state with dt=.002 s steps for 6 s. Input delays are 0,.08,.2 s, and the input limit is ±4.

Procedure

1. Use the practice notebook's function to plot outputs and inputs for the three delays.
2. Halve dt and check that the main responses remain similar.
3. Explain why linear-pole stability alone is insufficient to assess an implementation with delay and saturation.

Results to retain

• Output and input curves, RMS tracking error (after 2 s), input saturation fraction (the full 0–6 s record), and comparison of integration time steps.

def delayed(delay,dt=.002):
    times=np.arange(0,6+dt/2,dt);x=np.zeros(2);commands=[];ys=[];us=[]
    lag=int(round(delay/dt))
    for _ in times:
        command=float(np.clip(N-(K@x).item(),-4,4));commands.append(command)
        u=commands[-1-lag] if len(commands)>lag else 0.
        ys.append(x[0]);us.append(u)
        x=x+dt*(A0@x+B0[:,0]*u)
    return times,np.array(ys),np.array(us)
for delay in [0.,.08,.2]:
    tt,yy,uu=delayed(delay);keep=tt>=2;print(delay,np.sqrt(np.mean((1-yy[keep])**2)),np.mean(np.abs(uu)>=4-1e-8))
    plt.plot(tt,yy,label=str(delay))
plt.xlabel('Time [s]');plt.ylabel('Position');plt.legend();plt.grid()
# TODO: Add input plots and a comparison with the dt=.001 results.

Record results and self-check

Delay:

dt:

RMS tracking error:

Input saturation fraction:

Points to check

• Is delay implemented by retrieving past inputs from the command history?
• Did you avoid confusing time-step effects with model instability?

Official examples and data

• SciPy · Pole placement