Modeling Plant Variability as Uncertainty

In this first example, we have a family of models describing the plant behavior under various operating conditions. The nominal plant model is a first-order unstable system.

Pnom = tf(2,[1 -2])

Pnom = 2 ----- s - 2 Continuous-time transfer function.

The other models are variations ofPnom. They all have a single unstable pole but the location of this pole may vary with the operating condition.

p1 = Pnom*tf(1,[.06 1]);% extra lagp2 = Pnom*tf([-.02 1],[.02 1]);% time delayp3 = Pnom*tf(50^2,[1 2*.1*50 50^2]);% high frequency resonancep4 = Pnom*tf(70^2,[1 2*.2*70 70^2]);% high frequency resonancep5 = tf(2.4,[1 -2.2]);% pole/gain migrationp6 = tf(1.6,[1 -1.8]);% pole/gain migration

To apply robust control tools, we can replace this set of models with a single uncertain plant model whose range of behaviors includesp1throughp6. This is one use of the commanducover. This command takes an array of LTI modelsParrayand a nominal modelPnom和模型the differenceParray-Pnomas multiplicative uncertainty in the system dynamics.

Becauseucoverexpects an array of models, use the堆栈command to gather the plant modelsp1throughp6into one array.

Parray = stack(1,p1,p2,p3,p4,p5,p6);

Next, useucoverto "cover" the range of behaviorsParraywith an uncertain model of the form

P = Pnom * (1 + Wt * Delta)

where all uncertainty is concentrated in the "unmodeled dynamics"Delta(aultidynobject). Because the gain ofDeltais uniformly bounded by 1 at all frequencies, a "shaping" filterWtis used to capture how the relative amount of uncertainty varies with frequency. This filter is also referred to as the uncertainty weighting function.

Try a 4th-order filterWtfor this example:

orderWt = 4; Parrayg = frd(Parray,logspace(-1,3,60)); [P,Info] = ucover(Parrayg,Pnom,orderWt,'InputMult');

The resulting modelPis a single-input, single-output uncertain state-space (USS) object with nominal valuePnom.

P

P = Uncertain continuous-time state-space model with 1 outputs, 1 inputs, 5 states. The model uncertainty consists of the following blocks: Parrayg_InputMultDelta: Uncertain 1x1 LTI, peak gain = 1, 1 occurrences Type "P.NominalValue" to see the nominal value, "get(P)" to see all properties, and "P.Uncertainty" to interact with the uncertain elements.

tf(P.NominalValue)

ans = 2 ----- s - 2 Continuous-time transfer function.

波德图证实了t级hat the shaping filterWt"covers" the relative variation in plant behavior. As a function of frequency, the uncertainty level is 30% at 5 rad/sec (-10dB = 0.3) , 50% at 10 rad/sec, and 100% beyond 29 rad/sec.

Wt = Info.W1; bodemag((Pnom-Parray)/Pnom,'b--',Wt,'r'); grid title('Relative Gaps vs. Magnitude of Wt')

You can now use the uncertain modelPto design a robust controller for the original family of plant models, seeSimultaneous Stabilization Using Robust Controlfor details.

Simplifying an Existing Uncertain Model

In this second example, we start with a detailed uncertain model of the plant. This model consists of first-order dynamics with uncertain gain and time constant, in series with a mildly underdamped resonance and significant unmodeled dynamics. This model is created using theurealandultidyncommands for specifying uncertain variables:

gamma = ureal('gamma',2,'Perc',30);% uncertain gaintau = ureal('tau',1,'Perc',30);% uncertain time-constantwn = 50; xi = 0.25; P = tf(gamma,[tau 1]) * tf(wn^2,[1 2*xi*wn wn^2]);% Add unmodeled dynamics and set SampleStateDim to 5 to get representative% sample values of the uncertain model Pdelta = ultidyn('delta',[1 1],'SampleStateDim',5,'Bound',1); W = makeweight(0.1,20,10); P = P * (1+W*delta)

P = Uncertain continuous-time state-space model with 1 outputs, 1 inputs, 4 states. The model uncertainty consists of the following blocks: delta: Uncertain 1x1 LTI, peak gain = 1, 1 occurrences gamma: Uncertain real, nominal = 2, variability = [-30,30]%, 1 occurrences tau: Uncertain real, nominal = 1, variability = [-30,30]%, 1 occurrences Type "P.NominalValue" to see the nominal value, "get(P)" to see all properties, and "P.Uncertainty" to interact with the uncertain elements.

A collection of step responses illustrates the plant variability.

step(P,4) title('Sampled Step Responses of Uncertain System')

The uncertain plant modelPcontains 3 uncertain elements. For control design purposes, it is often desirable to simplify this uncertainty model while approximately retaining its overall variability. This is another use of the commanducover.

To useucoverin this context, first map the uncertain modelPinto an array of LTI models usingusample. This command samples the uncertain elements in an uncertain system and returns the corresponding LTI models, each model representing one possible behavior of the uncertain system. In this example, samplePat 60 points (the random number generator is seeded for repeatability):

rng(0,'twister'); Parray = usample(P,60);

Next, useucoverto cover all behaviors inParrayby a simple uncertainty modelUsys. Choose the nominal value ofPas center of the cover, and use a 2nd-order filter to model the frequency distribution of the unmodeled dynamics.

orderWt = 2; Parrayg = frd(Parray,logspace(-3,3,60)); [Usys,Info] = ucover(Parrayg,P.NominalValue,orderWt,'InputMult');

A Bode magnitude plot shows how the filter magnitude (in red) "covers" the relative variability of the plant frequency response (in blue).

Wt = Info.W1; bodemag((P.NominalValue-Parray)/P.NominalValue,'b--',Wt,'r') title('Relative Gaps (blue) vs. Shaping Filter Magnitude (red)')

You can now use the simplified uncertainty modelUsysto design a robust controller for the original plant, seeFirst-Cut Robust Designfor details.

Adjusting the Uncertainty Weighting

In this third example, we start with 40 frequency responses of a 2-input, 2-output system. This data has been collected with a frequency analyzer under various operating conditions. A two-state nominal model is fitted to the most typical response:

A = [-5 10;-10 -5]; B = [1 0;0 1]; C = [1 10;-10 1]; D = 0; Pnom = ss(A,B,C,D);

The frequency response data is loaded into a 40-by-1 array of FRD models:

loaducover_demosize(Pdata)

40x1 array of FRD models. Each model has 2 outputs, 2 inputs, and 120 frequency points.

Plot this data and superimpose the nominal model.

bode(Pdata,'b--',Pnom,'r',{.1,1e3}), grid legend('Frequency response data','Nominal model','Location','NorthEast')

Because the response variability is modest, try modeling this family of frequency responses using an additive uncertainty model of the form

P = Pnom + w * Delta

whereDeltais a 2-by-2ultidynobject representing the unmodeled dynamics andwis a scalar weighting function reflecting the frequency distribution of the uncertainty (variability in Pdata).

Start with a first-order filterwand compare the magnitude ofwwith the minimum amount of uncertainty needed at each frequency:

[P1,InfoS1] = ucover(Pdata,Pnom,1,'Additive'); w = InfoS1.W1; bodemag(w,'r',InfoS1.W1opt,'g',{1e-1 1e3}) title('Scalar Additive Uncertainty Model')传说('First-order w','Min. uncertainty amount','Location','SouthWest')

The magnitude ofwshould closely match the minimum uncertainty amount. It is clear that the first-order fit is too conservative and exceeds this minimum amount at most frequencies. Try again with a third-order filterw. For speed, reuse the data inInfoS1to avoid recomputing the optimal uncertainty scaling at each frequency.

[P3,InfoS3] = ucover(Pnom,InfoS1,3,'Additive'); w = InfoS3.W1; bodemag(w,'r',InfoS3.W1opt,'g',{1e-1 1e3}) title('Scalar Additive Uncertainty Model')传说('Third-order w','Min. uncertainty amount','Location','SouthWest')

The magnitude ofwnow closely matches the minimum uncertainty amount. Among additive uncertainty models,P3provides a tight cover of the behaviors inPdata. Note thatP3has a total of 8 states (2 from the nominal part and 6 fromw).

P3

P3 = Uncertain continuous-time state-space model with 2 outputs, 2 inputs, 8 states. The model uncertainty consists of the following blocks: Pdata_AddDelta: Uncertain 2x2 LTI, peak gain = 1, 1 occurrences Type "P3.NominalValue" to see the nominal value, "get(P3)" to see all properties, and "P3.Uncertainty" to interact with the uncertain elements.

You can refine this additive uncertainty model by using non-scalar uncertainty weighting functions, for example

P = Pnom + W1*Delta*W2

whereW1andW2are 2-by-2 diagonal filters. In this example, restrict useW2=1and allow both diagonal entries of W1 to be third order.

[PM,InfoM] = ucover(Pdata,Pnom,[3;3],[],'Additive');

Compare the two entries ofW1with the minimum uncertainty amount computed earlier. Note that at all frequencies, one of the diagonal entries ofW1has magnitude much smaller than the scalar filterw. This suggests that the diagonally-weighted uncertainty model yields a less conservative cover of the frequency response family.

bodemag(InfoS1.W1opt,'g*',...InfoM.W1opt(1,1),'r--',InfoM.W1(1,1),'r',...InfoM.W1opt(2,2),'b--',InfoM.W1(2,2),'b',{1e-1 1e3}); title('Diagonal Additive Uncertainty Model')传说('Scalar Optimal Weight',...'W1(1,1), pointwise optimal',...'W1(1,1), 3rd-order fit',...'W1(2,2), pointwise optimal',...'W1(2,2), 3rd-order fit',...'Location','SouthWest')

The degree of conservativeness of one cover over another can be partially quantified by considering the two frequency-dependent quantities:

Fd2s = norm(inv(W1)*w) , Fs2d = norm(W1/w)

These quantities measure by how much one uncertainty model needs to be scaled to cover the other. For example, the uncertainty modelPnom + W1*Deltaneeds to be enlarged by a factorFd2sto include all of the models represented by the uncertain modelPnom + w*Delta.

PlotFd2sandFs2das a function of frequency.

Fd2s = fnorm(InfoS1.W1opt*inv(InfoM.W1opt)); Fs2d = fnorm(InfoM.W1opt*inv(InfoS1.W1opt)); semilogx(fnorm(Fd2s),'b',fnorm(Fs2d),'r'), grid axis([0.1 1000 0.5 2.6]) xlabel('Frequency (rad/s)'), ylabel('Magnitude') title('Scale factors relating different covers')传说(“对角线标量因素”,...'Scalar to Diagonal factor','Location','SouthWest');

This plot shows that:

This indicates that thePnom+W1*Deltamodel provides a less conservative cover of the frequency response data inPdata.