Taratura PID, esercizio con codice MATLAB

MATLAB CODE

%-----------------------TARATURA PID-----------------------%

clear; close all;

%% CARICAMENTO DATI ESERCIZIO

%load datisperimentali.mat

s=zpk('s');

G = 0.7/((500*s+1)^4*(100*s+1));

t=0:0.01:10000;

N=1;

figure(1)

step(t,G);

title('Curva di processo dell''impianto originario');

xlabel('tempo');

ylabel('ampiezza');

[valori_y,time]=step(G,t);

W=feedback(G,1);

%% METODO CICLO CHIUSO

% Calcolo il margine di ampiezza ovvero il guadagno critico

[K_cr,Phi_m,omega_pi,omega_c]=margin(G);

% Determino il periodo critico

T_cr = 2*pi/omega_pi;

% Determini i parametri del PID

[Kp_closedloop, Ti_closedloop, Td_closedloop ] = TARATURA_CICLO_CHIUSO(K_cr, T_cr);

% Funzione di trasferimento PID

C_PID_closedloop=Kp_closedloop*(((1+1/(s*Ti_closedloop))+s*Td_closedloop/(1+(Td_closedloop/N)*s)));

% closed loop

W_closedloop = feedback(series(C_PID_closedloop,G),1);

% Grafico con risposta del ciclo chiuso

% figure(2)

% step(W_closedloop); 5

% title('Risposta ciclo chiuso metodo ciclo chiuso');

% xlabel('tempo');

% ylabel('ampiezza');

% Salvo le informazioni della risposta al gradino

parametri_risposta_closedloop=stepinfo(W_closedloop);

Tc_sup_closedloop=parametri_risposta_closedloop.RiseTime/10;

Tc_inf_closedloop=parametri_risposta_closedloop.RiseTime/20;

%% IMPLEMENTO METODO GRAFICO (Valutando la risposta al gradino ottengo i valori di: k,T,tau)

% Estraggo il valore assunto dalla risposa al termine del vettore temporale

% che è uguale

K=valori_y(end);

% Devo trovare il punto di flesso

dY = diff(valori_y);

ddY = diff(dY);

Ind = find(ddY > 0 & min(ddY)); %Indice del vettore in cui si cambia concavità

% Stima del tempo di flesso e della sua ordinata

Tflesso = time(Ind(end));

Yflesso = valori_y(Ind(end));

% Stimo il coefficiente angolare

m = (Yflesso - valori_y(Ind(end)-10))/(Tflesso - time(Ind(end)-10));

% Stimo il tempo di ritardo del modello del I ordine

tau_tan = (m*Tflesso - Yflesso)/m;

% Stimo la costante di tempo del modello del I ordine

Tt_tan = (tau_tan*K)/(abs(Yflesso-m*Tflesso));

% Funzione di trasferimento approssimata mediante metodo della tangente

G_tan = K*exp(-s*tau_tan)/(1+s*Tt_tan);

% Vado a comparare graficamento il modello cosi stimato e la curva di

% processo

% figure(3)

% step(t,G); 6

% title('Risposta sistema traccia');

% xlabel('tempo');

% ylabel('ampiezza');

% hold on;

% step(G_tan);

% title('Risposta sistema approssimato metodo tangente');

% Determino i valori dei parametri del PID

[Kp_tan,Ti_tan,Td_tan] = TARATURA_ZIEG_NIC(K,Tt_tan,tau_tan);

% Funzione di trasferimento del PID

C_PID_tan=Kp_tan*(((1+1/(s*Ti_tan))+s*Td_tan/(1+(Td_tan/N)*s)));

% closed loop

W_tan = feedback(series(C_PID_tan,G_tan),1);

% Grafico con risposta del ciclo chiuso

% figure(4)

% step(W_tan);

% hold on

% step(W);

% title('Risposta ciclo chiuso metodo tangente');

% xlabel('tempo');

% ylabel('ampiezza');

% Salvo le informazioni della risposta al gradino

parametri_risposta_tan=stepinfo(W_tan);

Tc_sup_tan=parametri_risposta_tan.RiseTime/10;

Tc_inf_tan=parametri_risposta_tan.RiseTime/20;

%% CALCOLO DELLE AREE

% Definisco i parametri che approssimano la curva di processo mediante una

% funzione che utilizza il metodo delle aree

[K_aree, T_aree, tau_aree] = METODO_AREE(time,valori_y);

% Funzione di trasferimento sistema approssimato metodo delle aree 7

G_aree = K_aree*exp(-s*tau_aree)/(1+s*T_aree);

% Vado a comparare graficamento il modello stimato e la curva di

% processo originale

% figure(5)

% step(t,G);

% title('Risposta sistema originale');

% hold on

% step(t,G_aree);

% xlabel('tempo');

% ylabel('ampiezza');

% step(G_aree);

% title('Risposta sistema approssimato metodo aree');

% Parametri del PID metodo aree

[Kp_aree,Ti_aree,Td_aree] = TARATURA_ZIEG_NIC(K_aree,T_aree,tau_aree);

% Funzione trasferimento PID metodo aree

C_PID_aree=Kp_aree*(1+1/(s*Ti_aree)+s*Td_aree/(1+(Td_aree/N)*s));

% Closed loop

W_aree = feedback(series(C_PID_aree,G_aree),1);

% figure(6)

% step(W_aree);

% hold on

% step(W);

% title('Risposta ciclo chiuso metodo aree');

% xlabel('tempo');

% ylabel('ampiezza');

% Salvo i dati della risposta al gradino

parametri_risposta_aree=stepinfo(W_aree);

Tc_sup_aree=parametri_risposta_aree.RiseTime/10;

Tc_inf_aree=parametri_risposta_aree.RiseTime/2; 8

% % METODO APPROSSIMAZIONE RITARDO PARI AL 5PERCENTO DEL GUADAGNO

% Cerco l'indice del vettore per cui la risposta raggiunge il 5% del suo

% valore finale

indice_k=1;

while valori_y(indice_k) <= 0.05*K

tau_5percent=time(indice_k);

indice_k=indice_k+1;

end

% Cerco l'indice del vettore per cui la risposta raggiunge il 63% del suo

% valore finale

while valori_y(indice_k) <= 0.63*K

Tt_5percent=time(indice_k);

indice_k=indice_k+1;

end

Tt_5percent=Tt_5percent-tau_5percent;

% Funzione di trasferimento sistema approssimato metodo 5%

G_5percent = K*exp(-s*tau_5percent)/(1+s*Tt_5percent);

%Vado a comparare graficamento il modello cosi stimato e la curva di

%processo

% figure(7)

% step(G);

% title('Risposta sistema originario');

% xlabel('tempo');

% ylabel('ampiezza');

% hold on;

% step(G_5percent);

% title('Risposta sistema approssimato metodo 5%');

% Determino i parametri del PID

[Kp_5percent,Ti_5percent,Td_5percent] = TARATURA_ZIEG_NIC(K,Tt_5percent,tau_5percent);

% Funzione di trasferimento PID

C_PID_5percent=Kp_5percent*(1+1/(s*Ti_5percent)+s*Td_5percent/(1+(Td_5percent/N)*s)); 9

% Close loop

W_5percent = feedback(series(C_PID_5percent,G_5percent),1);

% figure(8)

% step(W_5percent);

% hold on

% step(W);

% title('Risposta ciclo chiuso 5%');

% xlabel('tempo');

% ylabel('ampiezza');

%% CONFRONTO APPROSSIMAZIONI E RISPOSTE AL GRADINO

figure(2)

subplot(2,2,1)

step(G)

hold on

step(G_tan)

title('Approssimazione metodo tangente')

xlabel('tempo');

ylabel('ampiezza');

subplot(2,2,2)

step(G)

hold on

step(G_5percent)

title('Approssimazione metodo 5%')

xlabel('tempo');

ylabel('ampiezza');

subplot(2,2,3) 10

step(G)

hold on

step(G_aree)

title('Approssimazione metodo aree')

xlabel('tempo');

ylabel('ampiezza');

subplot(2,2,4)

step(G)

hold on

step(G_tan)

step(G_5percent)

step(G_aree)

title('Confronto approssimazioni')

xlabel('tempo');

ylabel('ampiezza');

%Grafico tutte le risposte

figure(3)

subplot(2,2,1)

step(W)

hold on

step(W_tan)

title('Risposta metodo tangente')

xlabel('tempo');

ylabel('ampiezza');

subplot(2,2,2)

step(W)

hold on

step(W_5percent) 11