MATLAB Examples

Specifying Linearization for Model Components Using System Identification

This example shows how to use System Identification Toolbox to identify a linear system for a model component that does not linearize well and use the identified system to specify its linearization. Note that running this example requires Simscape Power Systems in addition to System Identification Toolbox.


Linearizing Hard Drive Model

Open the simulink model for the hard drive.

model = 'scdpwmharddrive';

In this model, the hard drive plant is driven by a current source. The current source is implemented by a circuit that is driven by a Pulse Width Modulation (PWM) signal so that its output can be adjusted by the duty cycle. For details of the hard drive model, see the example "Digital Servo Control of a Hard-Disk Drive" in Control System Toolbox™ examples.

PWM-driven circuits usually have high frequency switching components, such as the MOSFET transistor in this model, whose average behavior is not well defined. Thus, exact linearization of this type of circuits is problematic. When you linearize the model from duty cycle input to the position error, the result is zero.

SimscapeSolver = [model '/PWM driven converter/Solver Configuration'];
set_param(SimscapeSolver,'UseLocalSolver','off'); %turn off for linearization
io(1) = linio('scdpwmharddrive/Duty cycle',1,'input');
io(2) = linio('scdpwmharddrive/Hard Disk Model',1,'output');
sys = linearize(model,io)
set_param(SimscapeSolver,'UseLocalSolver','on') %turn back on for simulation
sys =
  D = 
                 Duty cycle
   position err           0
Static gain.

Finding a Linear Model for PWM Component

You can use frequency response estimation to obtain the frequency response of the PWM-driven current source and use the result to identify a linear model for it. The current signal has a discrete sample time of 1e-7. Thus, you need to use a fixed sample time sinestream signal. Create a signal that has frequencies between 2K and 200K rad/s.

idinput = frest.createFixedTsSinestream(Ts,{2000,200000});
idinput.Amplitude = 0.1;

You can then define the input and output points for PWM-driven circuit and run the frequency response estimation with the sinestream signal.

pwm_io(1) = linio('scdpwmharddrive/Duty cycle',1,'input');
pwm_io(2) = linio('scdpwmharddrive/PWM driven converter',1,'openoutput');
sysfrd = frestimate(model,pwm_io,idinput);

Using the N4SID command from System Identification Toolbox, you can identify a second-order model using the frequency response data. Then, compare the identified model to the original frequency response data.

sysid = ss(tfest(idfrd(sysfrd),2));

We used frequency response data with frequencies between 2K and 200K rad/s. The identified model has a flat magnitude response for frequencies smaller than 2K. However, our estimation did not have include for those frequencies. Assume that you would like to make sure the response is flat by checking the frequency response for 20 and 200 rad/s. To do so, create another input signal with those frequencies in it.

lowfreq = [20 200];
inputlow = frest.createFixedTsSinestream(Ts,lowfreq)
The sinestream input signal:
      Frequency           : [20 200] (rad/s)
      Amplitude           : 1e-05
      SamplesPerPeriod    : [3141593 314159]
      NumPeriods          : 4
      RampPeriods         : 0
      FreqUnits (rad/s,Hz): rad/s
      SettlingPeriods     : 1
      ApplyFilteringInFRESTIMATE (on/off)    : on
      SimulationOrder (Sequential/OneAtATime): Sequential

In the input signal parameters, we can see that having a very fast sample rate of 1e-7 seconds (10 MHz sampling frequency) for the frequencies 20 and 200 rad/s cause high SamplesPerPeriod values of 3141593 and 314159. Considering that each frequency has 4 periods, frequency response estimation would log output data with around 14 millions samples. This would require a lot of memory and it is quite likely that you might run into memory issues running the estimation.

Obviously, you do not need such a high sampling rate for analyzing 20 and 200 rad/s frequencies. You can use a smaller sampling rate to avoid memory issues:

Tslow = 1e-4;
wslow = 2*pi/Tslow;
inputlow = frest.createFixedTsSinestream(Tslow,wslow./round(wslow./lowfreq));
inputlow.Amplitude = 0.1;

To make the model compatible with the smaller sampling rate, resample the output data point using a rate transition block as in the modified model:

modellow = 'scdpwmharddrive_lowfreq';

You can now run the analysis for the low frequencies and compare it against identification result.

load scdpwmharddrive_lowfreqresults.mat
% sysfrdlow = frestimate(modellow,getlinio(modellow),inputlow);

Specifying the Linearization for PWM Component

As you verified using the frequency response estimation, the low-frequency dynamics of the PWM-driven component are captured well by the identified system. Now you can make linearization use this system as the linearization of the PWM-driven component. To do so, specify block linearization of that subsystem as follows:

pwmblock = 'scdpwmharddrive/PWM driven converter';
rep = struct('Specification','sysid',...
set_param('scdpwmharddrive/Duty cycle','SampleTime','Ts_plant');

Linearizing the model after specifying the linearization of the PWM component gives us a non-zero result:

set_param(SimscapeSolver,'UseLocalSolver','off') %turn off for linearization
sys = linearize(model,io);
set_param(SimscapeSolver,'UseLocalSolver','on') %turn back on for simulation

You might still like to validate this linearization result using frequency response estimation. Doing this as below verifies that our linearization result is quite accurate and all the resonances exist in the actual dynamics of the model.

valinput = frest.Sinestream(sys);
valinput = fselect(valinput,3e3,1e5);
valinput.Amplitude = 0.1;
sysval = frestimate(model,io,valinput);

Close the model: