How can I specify less oschillatory behaviour for PID

How can I speicfy that I do not want a heavily underdamped response as seen in the image attached? The addition of the derivative term makes it difficult to implement in simulink enviroment.

--------------------------------------------------------------First script------------------------------------------------------------

clear all, close all, clc

dt = 0.000005; % this is 10^-6 0.000001

PopSize = 40;

MaxGenerations = 25;

s = tf('s');

G =(1.44e09)/((s*s)+5333*s+9.6e07)

options = optimoptions(@ga,'PopulationSize',PopSize,'MaxGenerations',MaxGenerations,'OutputFcn',@myfun);

[x,fval] = ga(@(K)pidtest(G,dt,K),2,-eye(2),zeros(2,1),[],[],[],[],[],options);

--------------------------------------------------------------Second script--------------------------------------------------------

function J = pidtest(G,dt,parms)

s = tf('s');

K = parms(1)+ parms(2)/s

Loop = series(K,G);

ClosedLoop = feedback(Loop,1);

t = 0:dt:0.005; % this indicates length of time to show

e = 1 - step(ClosedLoop,t);

J = sum(t'.*abs(e)*dt);

step(1*ClosedLoop,t)

h = findobj(gcf,'type','line');

set(h,'linewidth',1);

drawnow

---------------------------------------------------------------end of script------------------------------------------------------

[Mathieu NOE]

hello

have you tried other PID tuning algorithm ? your code also does not include the derivative term of the PID wich accounts for damping your oscillations

[JAMES KEEN]

Hello

I have tried other algorithms. The last I tried was a LQR cost function which has a persistant steady state error of 0.05v but aside from this shows to be great.

The following code seen in the question presents the ITAE cost function which when done as a PI cannot control itself prffoperly. I got rid of the derivative term as I could not port the derivative over to simulink without significant differences in the actucal response.

In short by including the derivative term in the matlab coding, the response would not be the same in the simulink enviroment. I am having issues with making the derivative response be the same as within the matlab coding enviroment.

[Mathieu NOE]

hmm , I don't remember having the same issue regarding the derivative term in simulink models - didi you use the "true" derivative block from the simulink library ? I prefered to use a high pass filter which limits the noise increase.

Never had any issue with that method. i always used fixed step solver / discretized blocks

[JAMES KEEN]

I have tried true derivative block and PID controller block with high pass filter that has N=1e06 as seen in the code.

I cannot use fixed step solver /discretized block I think because the buck converter circuit is continious.

[Mathieu NOE]

so the G transfer function is a linear approximation of the buck converter ?

G =(1.44e09)/((s*s)+5333*s+9.6e07)

[JAMES KEEN]

Yes you are correct in that thought.

[Mathieu NOE]

and the portion of the simulink model with that tf + the PID (including a D term) is ok for you ?

[JAMES KEEN]

The image I have just shown you has gained PID controller values from the following code. It has an offset of 0.05v which I do not like though so I have been trying various codes.

----------------------------1st script-----------------------------------

clear all, close all, clc

dt = 0.000001; % this is 10^-6 0.000001

PopSize = 40; % was 40

MaxGenerations = 25; % was 25

s = tf('s');

G =(1.44e09)/((s*s)+5333*s+9.6e07)

options = optimoptions(@ga,'PopulationSize',PopSize,'MaxGenerations',MaxGenerations,'OutputFcn',@myfun);

[x,fval] = ga(@(K)pidtest(G,dt,K),3,-eye(3),zeros(3,1),[],[],[],[],[],options);

----------------------------------------------------------------------2nd script-----------------------------

function J = pidtest(G,dt,parms)

s = tf('s');

K = parms(1)+ parms(2)/s + parms(3)*s/(0.000001*s+1)%/(1+0.000001*s)% this is 10^-6 0.000001

Loop = series(K,G);

ClosedLoop = feedback(Loop,1);

t = 0:dt:0.005; % this indicates length of time to show

[y,t] = step(ClosedLoop,t);

u = lsim(K,1-y,t);

Q = 1;

R = 0.000005; % 5*10^-6 change these to effect rise, overshoot, settling time? 10^-8 rn

J = dt*sum(Q*(1-y(:)).^2 + R*u(:).^2)

step(5*ClosedLoop,t)

h = findobj(gcf,'type','line');

set(h,'linewidth',1);

drawnow

----------------------------------------end of script------------------------------------------

-----------------------------------------------------end of code used on that simulink model-----------------------------------------------------------------------------------------------------------------------------------------------------------

I have a code which shows to eliminate the steady state error but it does not translate onto the simulink model shown. The new code which gets rid of the steady state error is shown below.

----------------------------1st script-----------------------------------

clear all, close all, clc

dt = 0.000001; % this is 10^-6 0.000001

PopSize = 40; % was 40

MaxGenerations = 25; % was 25

s = tf('s');

G =(1.44e09)/((s*s)+5333*s+9.6e07)

options = optimoptions(@ga,'PopulationSize',PopSize,'MaxGenerations',MaxGenerations,'OutputFcn',@myfun);

[x,fval] = ga(@(K)pidtest1(G,dt,K),3,-eye(3),zeros(3,1),[],[],[],[],[],options);

----------------------------------------------------------------------2nd script-----------------------------

function J = pidtest(G,dt,parms)

s = tf('s');

K = parms(1)+ parms(2)/s + parms(3)*s/(0.000001*s+1) %/(1+0.000001*s)% this is 10^-6 0.000001

Loop = series(K,G);

ClosedLoop = feedback(Loop,1);

t = 0:dt:0.005; % this indicates length of time to show

[y,t] = step(ClosedLoop,t);

S=stepinfo(ClosedLoop);

rt1=S.RiseTime;

rt2=S.SettlingTime;

rt3=S.Overshoot;

%%% J = (dt*sum((1-y(:)).^2))+rt1+rt2+rt3;

u = lsim(K,1-y,t);

Q = 1;

R = 0.000005; % 5*10^-6 change these to effect rise, overshoot, settling time? 10^-8 rn

J = dt*sum(Q*(1-y(:)).^2 + R*u(:).^2)+rt1+rt2+rt3;

step(5*ClosedLoop,t)

h = findobj(gcf,'type','line');

set(h,'linewidth',1);

drawnow

----------------------------------------end of script------------------------------------------

[Mathieu NOE]

hi again

I tried this PID tuning code and found it pretty effective

V-Tiger - File Exchange - MATLAB Central (mathworks.com)

results :

kp = 0.1926

ki = 1.9297e+03

kd = 2.2012e-05

attached the code (demo file modified for your purpose)

% vtiger_demo demonstrates V-Tiger.
    % Parameter setting of Plant G(s), sampling time ts[s], and so on.
    ts=1e-5;    % sampling time [s] (original value was 0.01)
    s=tf('s');  % complex variable of Laplace transform
    z=tf('z',ts);   % shift operator
    p=(1-1/z)/ts;   % differential operator based on backward Euler's rule
%     Gs = 5/(0.01*s^2+0.2*s+10)*exp(-0.1*s);   % plant G(s) to be controlled
%     Gs = 1.44e09/(s^2+5333*s+9.6e07);   % plant G(s) to be controlled
    Gs = 1.44e09/(s^2+5333*s+9.6e07);   % plant G(s) to be controlled
    
    G = c2d(Gs,ts); % G(z) is derived by discretizing G(s) with zero order holder

    % Design initial controller K0 using  The Ziegler-Nichols rule
    [Ku,Pm,Wu,Wcp] = margin(G);    % get Gaim margin Ku at Wu[rad/s]
    Tu = 1/(Wu/2/pi);   % When K=Ku, self-excited vibration with period Tu[s] will occur.
    kp0=0.6*Ku; ki0=kp0/(0.5*Tu);    kd0=kp0*0.125*Tu;  % ZN classical parameters
    K0 = kp0 + ki0/p + kd0*p;   % ZN PID controller K0(z)
    
    % Measuring one-shot data y00(t) and u00(t)
    u00=ones(300,1);    % input u00(t) is a step function
    u00(1)=0;   % initial value must be different from the other values <-- IMPORTANT!
    y00=lsim(G,u00);    % y00 is simulated
    r=u00;  % reference input to feedback system
    
    % Step 1) Make step responses to cyclic, and get frequency data.
    freq.y0jw = fft4step(y00);  % y0(j w) from y00(t)
    freq.u0jw = fft4step(u00);  % u0(j w) from u00(t)
    freq.r0jw = fft4step(r);    % r0(j w) from r(t)
    freq.p = fft4tf(p,length(u00)*2);   % p(j w) from differential operator p
    freq.r = r;   % r(t) is reference input to feedback system
    freq.wST=0.05;   % Error band of settling time for cost function (Settling Time Threshold)
    freq.OVr=5;   % Overshoot [%] for constraints
    freq.GMr=6;   % Gain margin [dB] for constraints:  Regulator 3-10dB, 20-inf deg
    freq.PMr=45;  % Phase margin [deg] for constraints:Servo    10-20dB, 40-60  deg
    
    % Step 2) Optimize PID gains by evaluating overshoot, settling time,
    %         and stability margins using virtual time response.
    [kp,ki,kd] = vtigerPID(freq,[kp0 ki0 kd0]);   % Get optimum PID gains. [kp0 ki0 kd0] is initial value for optimization
    K = kp + ki/p + kd*p;  % PID controller by V-Tiger
    
    disp('-----------------------------------------------------')
    % Verify the controller of V-Tiger and ZN by simulations
    [y,u] = freq2yu(freq,K); % Virtual time responses predicted by V-Tiger when K is used
    Gcl = feedback(ss(G*K),1);  % closed loop transfer function. feedback(a,b)=a/(1+a*b). ss(G) is ss(Gcl);  % State space representation from G
    yt=lsim(Gcl,r); % yt is true y(t) simulated using true plant model G(z)
    yZN=lsim(feedback(ss(G*K0),1),r);% y(t) by K0 (ZN) is simulated using true plant model G(z)
    figure(1),
    p1=plot([yt yZN y-yt y00]);
    hold on
    vtigerPID(freq,[kp ki kd],1);
    hold off; grid; xlabel('sample number k (0.01k [sec])')
    legend(p1,'y (V-Tiger)','y (ZN)','error of true/virtual y','y_{00}','Location','southeast');
    title('PID control result. V-Tiger is better than Ziegler-Nichols rule')
    disp('V-Tiger has completed controller design using y00 instead of the plant model.')
    disp(['Plant G(s) to be discretized with sampling time ts=' num2str(ts) '[s] is as follows:']), 
    [yZN,iZN]=max(yZN); text(iZN,yZN,'\leftarrow y (ZN)','Color','red','FontSize',14);
    [yvt,ivt]=max(yt); text(ivt,yvt+0.06,['y (V-Tiger)';'\downarrow '],'Color','blue','FontSize',14);
    text(length(y)*0.16,y00(end)*0.8,['\uparrow                                      '; ...
                                      'y_{00} (used by V-Tiger instead of model G(z))'],'Color','magenta','FontSize',14);
    text(length(y)*0.25,max(y)*0.58,[ ...
    'V-Tiger optimization is as follows:                       '; ...
    '  Cost function: Settling time (error band is \pm 3%)     '; ...
    '  Constraints: Overshoot<3%, Stability margins> 3dB, 20deg'])
    Gs, disp('');
    disp('Fig.1 shows step resonses. V-Tiger is better than ZN method.'),
    disp('press any key to type code of "vtiger_demo.m".'), 
    pause
    disp('-----------------------------------------------------')
    figure(2), margin(G*K);    text(10,-1000,'Open-loop by V-Tiger')
    dbtype vtiger_demo 1:37
    disp('-----------------------------------------------------')
    
    % one more figure for myself 
     figure(3),plot(yt)
    title('PID control result. V-Tiger is better than Ziegler-Nichols rule')
    xlabel('sample number k (0.01k [sec])')
    ylabel('step response')
    kp
    ki
    kd
    

[Tesfaye Girma]

good question bro

[JAMES KEEN]

Hi I am glad you found a good code to use although I am required to use Genetic algorithm with any cost function for my present project. I appreciate the help neverthless.

[Mathieu NOE]

OK - did you get similar values ?

I believe I saw also some submissions of PID tuning using GA in the FEX section - did you check that ?

[JAMES KEEN]

Hi I get the following PID values when using my LQR cost function.

P= 7.23803323894401

I= 0.617502138949610

D= 0.000315064097110046

filter coefficent = 1e06

I have not checked that submission as I have not seen it

Answer 1

[Tesfaye Girma]

figure(3),plot(yt)

title('PID control result. V-Tiger is better than Ziegler-Nichols rule')

xlabel('sample number k (0.01k [sec])')

ylabel('step response')

gv

gi