no oscillation after put in pid controller

Question

I have a problem regarding putting a PID controller in my simulink file.

In my simulink file, i used pid controller to control my process. I used s-function as my process block diagram.

According Ziegler-Nichols method, for the first step we set k equal to the smallest value (0.5) so I put 0.5 in my proportional value . but there is no different between the result with controller and without controller.Even i increase or decrease proportional value.

Why this problem occur? hope someone can help me.Thank you.

my block diagram is look like below picture.refer to this picture

http://s1009.photobucket.com/albums/af218/syarinazulkifeli/?action=view&current=Untitled-1.png

Here is my s-function file:

function [sys,x0,str,ts]= reactor_sfcn(t,x,u,flag)

switch flag

case 0
    [sys,x0,str,ts]=mdlInitializeSizes;
case 1,
    sys = mdlDerivatives(t,x,u);
case 3,
    sys = mdlOutputs(t,x,u);
case 9
    sys =[];
end

function [sys,x0,str,ts] = mdlInitializeSizes()

s = simsizes;
s.NumContStates   = 11;
s.NumDiscStates   = 0;
s.NumOutputs      = 11;
s.NumInputs       = 1;
s.DirFeedthrough  = 0;
s.NumSampleTimes  = 1;

sys = simsizes(s) ;
x0 = [0.0258,0,0,0,0,0,0,0,8.83,303.15,303.15];
str=[] ; 
ts = [0 0];

    function sys = mdlDerivatives (t,x,u)
Tjo = u;
sys = reactor(t,x,Tjo);


        function sys = mdlOutputs(t,x,u)
%                 sys(1)=x(1);
%                 sys(2)=x(2);
%                 sys(3)=x(3);
%                 sys(4)=x(4);
%                 sys(5)=x(5);
%                 sys(6)=x(6);
%                 sys(7)=x(7);
%                 sys(8)=x(8);
%                 sys(9)=x(9);
%                sys(10)=x(10);
%                sys(11)=x(11);
                 sys = x;  

Link of s-function file

function DXDT = reactor(t,x,Tjo)

% -------------------------------------------- %
%       Parameters definition
% -------------------------------------------- %  
I  = x(1);   X = x(2);   P0 = x(3);  
P1 = x(4);  P2 = x(5);   Q0 = x(6);  
Q1 = x(7);  Q2 = x(8);    M = x(9);
Tjo = x(10);  T = x(11); 

% ---------------------------------------
%               Constants
% =======================================

%constant
 M0 = 8.83;
 I0 = 0.0258;
 Qh = 60.44;
 MMWm = 100.3;
 U0 = 55.1;

% Densities
dp = 1200;                
dm = 968-1.225*(T-273.15);               
Rhoc = 998;                            

% volume expansion factor
Fev = (dm - dp)./dp ;                  

% volume fraction
 Fvm = (1 - X)./(1 + Fev*X);           
 Fvp = X*(1-Fev)./(1+Fev*X);

% Total reactant mixture density
 Rho = dm*Fvm + dp*Fvp;            
% Reactor and jacket volume
 Vc = 2;                            
  V = 2;                           
% Reactor dimension      
At =3.1416*0.15*0.113;    
% coolant flow rate
Mc = 0.41/18;                      
%   
 Cpc = 77.22;            
 Cp = 199.13;               
% Average coolant temperature
Tji = 303.15;
Tj =(Tji+Tjo)/2;     

% Overall heat transfer coeff   
 alpha = 0.4;  
 U = U0-alpha*X;        

 delHp = 57800;                        

% ---------------------------------------
%           Rates of reaction
% ======================================= 

% Fujita-Doolittle equation
Tgp = 114 +273.25   ;                      
A = 0.168-8.21e-6*(T-Tgp)^2;
B = 0.03;
g = exp(2.303*Fvm/(A + B*Fvm));          

% Dissociation rate
F = 0.58; 
Kd = (6.32e16*exp(-15.43e3/T));            

% Propagation rate  
Tep = 5.4814e-16*exp(13982/T);               
Kp0 = 2.952e7*exp(-4353/(1.987*T));     
Kp = (Kp0*g)./(g + (Tep*P0.*Kp0));         

% Termination rate  
Tet = (1.1353e-22*exp(17420/T))./I0;      
Kt0 = 5.88e9*exp(-701/(1.987*T));        
Kt = (Kt0*g)./(g + (Tet*P0.*Kt0));        
Ktc = 0;
Ktd = Kt ;                                



% -------------------------------------------- %
%       ODE's
% -------------------------------------------- % 
 dIdt = -Kd*I - ((Fev*I.*P0.*(1 - X)*Kp)./(1 + Fev*X)); 
 dXdt = Kp*(1 - X).*P0;     
 dP0dt = (-Fev*P0.*P0.*(1 - X)./(1 + Fev*X)).*Kp + 2*F*Kd*I - Kt*P0.*P0; 
 dP1dt = (-Fev*P0.*P1.*(1 - X)./(1 + Fev*X)).*Kp + 2*F*Kd*I - Kt*P0.*P1 + (Kp*M0*(1 -   X)./(1 + Fev*X)).*P0;
dP2dt = (-Fev*P0.*P2.*(1 - X)./(1 + Fev*X)).*Kp + 2*F*Kd*I - Kt*P0.*P2 + (Kp*M0*(1 - X)./(1 + Fev*X)).*(2*P1 + P0);
dQ0dt = (-Fev*P0.*Q0.*(1 - X)./(1 + Fev*X)).*Kp + Ktd*P0.*P0 + 0.5*Ktc*P0.*P0;
dQ1dt = (-Fev*P0.*Q1.*(1 - X)./(1 + Fev*X)).*Kp + Ktd*P0.*P1 + Ktc*P0.*P1;
dQ2dt = (-Fev*P0.*Q2.*(1 - X)./(1 + Fev*X)).*Kp + Ktd*P0.*P2 + Ktc*(P1.*P0 + P1.^2);
dMdt = (-Kp*P0*M0*(1 - X)./(1 + Fev*X)).*((Fev*(1 - X)./(1 + Fev*X)) + 1);
  Rm = (-delHp)*dMdt;                  
dTjodt = (Mc*Cpc*(Tji-Tjo)+ U*At*(T-Tj))/(Vc*Rhoc*Cpc/18);
dTdt = (Qh+(Rm*V*MMWm)-(U*At*(T-Tj)))/(V*Rho*Cp);


DXDT =[dIdt;dXdt;dP0dt;dP1dt;dP2dt;dQ0dt;dQ1dt;dQ2dt;dMdt;dTjodt;dTdt];

Answer 1

You change the Kp term but there is no change...did you give the system a step command? Is the actual and desired signal error being fed back to the PID controller? If these are satisfied and the system is still insensitive to Kp you may need to go through your model block by block to find the issue. Also, why are you modelling your plant with an S function? Implementing these can be tricky. I would much rather see the plant model in diagram form, or at least in embedded m.

For the step input I would recommend this link or google it. You need to set the Kp fairly high and then give it a step command input to destabilize the system and then vary Kp until it's stable and measure the period of oscillation.

I put together a toy model for you to test, I am inclined as well to beleive your s-function is at fault here. If you want to post/send the sfunction, I'd be glad to look at it. However, if you want to learn Ziegler-Nichols, start with a very simple model and follow the steps in the Ziegler-Nichols process. Here's the output for some of my data with a Kp = 200 and a plant of 1/(s+1):

So you can see the oscillation above. Kp = 200 is way too high, you'd have to reduce it quite a bit for Ziegler-Nichols, but I just wanted to give you an example.

EDIT I downloaded your readtor_sfcn and reactor function into two files named reactor_sfcn.m and reactor.m Right away, I can tell why you never see any change in the input. In the reactor_sfcn mdlDerivatives, you pass the input to the system u as the third parameter to the reactor function. In the reactor function I can see that the third parameter, Tjo is never read, but overwritten with x(10). So, if Tjo is supposed to be the input parameter it shouldn't also be a state. That will be up to you to solve that problem, it is specific to your implementation of the plant. I've got a test model I've used to look at your files, I'll try to put it somewhere you can get access to it soon, but it likely won't be as useful as you analyzing what the plant is supposed to be doing. Hope this helps!

Thanks!

Answer 2

First I woult check if the s-function for the process is working as expected. Disconnect the PID controller and connect a step, ramp or a sin block. This could give you a hint if the block for the process is OK and how large the static coefficient is.

A second hint: the Ziegler-Nichols method says one should increase the Kp until the system gets marginaly stable. Have you tested only one value for the proportional gain? Depending on the plant 0.5 can be really small value and way under the levels you need to control the system.