MATLAB: Is it possible to scale the LTI models for greater accuracy in Control Systems Toolbox 6.0 (R14)

Control System Toolbox

I have created a stable transfer function in Control System toolbox 6.0(R14). However, when I try to plot the step and impulse responses for this transfer function using the STEP and IMPULSE functions, the response does not appear to be completely stable.

I would like to introduce some scaling to this transfer function so that the generated step and impulse responses are completely stable.

Best Answer

The ability to scale the LTI models is available with the help of the PRESCALE function in Control System Toolbox 8.2 (R2008b). Note that the LTI model must be in state space form before it is passed in to the PRESCALE function.

For previous versions of the Control System Toolbox, refer to an example below in which a transfer function is scaled using the appropriate coefficients.

num = [8.76e+020 4.06464e+027];
den = [1 6.8e+007 3.09e+014 4.38e+020 2.03e+026];
pzmap(num,den) % Poles are on the left half of s-plane. System is stable.
step(num,den); %Response appears to be tending to a steady state value.   
s = 1e6; % Introduction of the scaling coefficients. Replacing s with 1e6 in the                              %following step & divide numerator and denominator by 1e25
v = s .^ [4 3 2 1 0];
den1 = ( den .* v )/ 1e25;
num1 = ( num .* [ss 1]) / 1e25;
step(num1,den1)% Now, the response is much more stable.

In the stable response plot, we notice that time has also been scaled by a factor of 1e6. For a better understanding of the significance of scaling, we could also inspect the condition number for the system matrix A.

The condition number of a matrix measures the sensitivity of the solution of a system of linear equations to errors in the data. It gives an indication of the accuracy of the results from matrix inversion and the linear equation solution. Condition numbers with values closer to 1 indicate well conditioned matrices.

This inspection can be done for each of the above two cases as follows.

[a1 b1 c1 d1]=tf2ss(num,den);     % not scaled
[a2 b2 c2 d2]=tf2ss(num1,den1);   % scaled
cond(a1) % Condition number of the non scaled matrix a1
cond(a2) % Condition number of the scaled matrix a2

Related Solutions

MATLAB: Does relative error model reduction using the BSTSCHMR and BSTSCHML functions in the Robust Control Toolbox 2.0.10 (R14) produce wrong results

This phenomenon is a result of ill-conditioned system which causes numerical inaccuracies in the model reduction algorithm.

This is often a result of a nearly singular D matrix. In this case, adjustments to this D matrix may result in a better solution. However, this may distort the systems high frequency dynamics. This procedure is shown in the following code (variables "sys" and num_reduced_states must be defined in the workspace):

% With a system with a nearly singular D matrix
% Assume sys.d = 1e-5
sys1 = sys;
sys1.d = 1e-2; %this will disregard the higher frequency dynamics
[sysr1,aug,hsv] = bstschml(sys1,1,num_reduced_states);
sysr1a = sysr1;
sysr1a.d = sys.d; % put D back to original value
% Now compare the frequency response

bode(sys, sysr1, sysr1a);

As an alternative, if you have the Mu-Analysis and Synthesis Toolbox, you can first try rebalancing the system using the SRELBAL function before performing this model reduction. This would be done as follows:

[a,b,c,d] = ssdata(sys);
sysmat = pck(a,b,c,d);
[sysmatrb,sigrel] = srelbal(sysmat);
[a2,b2,c2,d2] = unpck(sysmatrb);
sys2 = ss(a2,b2,c2,d2); 
[sysr2,aug,hsv] = bstschml(sys,1,num_reduced_states);
% Now compare the frequency response
bode(sys,sysr2);

MATLAB: Do I not get a perfect fit while estimating a frequency response data with a higher natural frequency using the System Identification Toolbox 7.3 (R2009a)

This is an expected behavior while using System Identification Toolbox 7.3 (R2009a). The reason for this behavior is as follows:

IDPROC imposes a lower limit on the value of Tw = 1/wn = 1e-3 as a default. For 4 kHz, wn is 2.5133e4 which leads to Tw value of the order of 1e-5.

Hence one is unlikely to get anything good using default settings. One must configure the parameter value bounds before estimating, as in:

f = 4000;
wn = 2*pi*f; %rad/sec
zeta = 0.25;
tf('s');
num = [1];
sec= 2*wn*zeta;
third = wn^2;
den = [1 sec third];
csys = tf(num,den); % continuous transfer function
Ts = 0.00001; % Ts = 0.00001 sec;
dsys = c2d(csys,Ts); 
[P,PHA,W] = bode(dsys); % bode of discrete sys, W is in rad/sec, P is in Db, PHA is in degrees
zfr = P.*exp(i*PHA*pi/180);
gfr = idfrd(zfr,W,Ts);
figure(1);
bode(gfr), legend('gfr')
mproc = idproc('p2u');
mproc.Tw.min = 1e-6; % lower the min bound to be smaller than expected value
m = pem(gfr, mproc);

In general, it is advisable to follow this multi-step process rather than directly doing "pem(data, 'p2u')". One should carefully inspect the min and max values of each parameter and preferably also set the initial guess values in the 'value' field of parameter (wherever possible).

Best Answer

Related Solutions

MATLAB: Does relative error model reduction using the BSTSCHMR and BSTSCHML functions in the Robust Control Toolbox 2.0.10 (R14) produce wrong results

MATLAB: Do I not get a perfect fit while estimating a frequency response data with a higher natural frequency using the System Identification Toolbox 7.3 (R2009a)

Related Question