Avancement avant reu

BE_Beam.m

@@ -55,6 +55,7 @@ eig_zoh2 = eig(c2d(sys, Ts2, 'zoh')) ;
 eig_zoh3 = eig(c2d(sys, Ts3, 'zoh')) ;
 %%
 Question4 = figure ;
+
 subplot(1,3,1) ;
 plot(eig_cont, 'x');
 axis square;
@@ -75,10 +76,10 @@ plot(cos(theta), sin(theta), 'k--', 'HandleVisibility','off','LineWidth',0.25);
 xlim([-1.1, 1.1]);
 ylim([-1.1, 1.1]);
 grid on;
-title({'Valeurs propres en temps discret',' avec la méthode de Tustin'});
 legend("Ts1", "Ts2", "Ts3");
 xlabel('Réel');
 ylabel('Imaginaire');
+title({'Valeurs propres en temps discret',' avec la méthode de Tustin'});
 hold off;
 
 subplot(1,3,3);
@@ -88,7 +89,6 @@ plot(eig_zoh2,'x', 'LineWidth',1);
 plot(eig_zoh3,'x', 'LineWidth',1);
 axis square;
 legend("Ts1", "Ts2", "Ts3");
-% theta = linspace(0, 2*pi, 100);
 plot(cos(theta), sin(theta), 'k--', 'HandleVisibility','off','LineWidth',0.25);
 xlim([-1.1, 1.1]);
 ylim([-1.1, 1.1]);
@@ -135,12 +135,12 @@ end
 [minimum, min_index] = min(TransientTime)
 k = min_index * 0.01 
 
-%% Question 8 
-
 H = inv(Cw_L*inv(-A+B*C*k)*B)
+
+%% Question8
 [y, t, x] = lsim(feedback(sys,k), H*cmd, t) ;
 u = H*cmd'-k*y ;
-%%
+
 Question8 = figure ;
 plot(t, u, t, y, t, cmd, t, (Cw_L*x')');
 xlabel("Temps (s)");
@@ -154,35 +154,40 @@ end
 
 %% Question 9
 
-sys_CL = feedback(sys,k);
-% [y, tOut] = lsim(sys_CL, H.*cmd, t);
-% plot(t, Cw(L)*y, t, u)
-eig_cont_CL = eig(sys_CL);
+
+eig_cont_CL = eig(feedback(sys,k));
+eig_zoh1_CL = eig(feedback(c2d(sys, Ts1, 'zoh'), k)) ;
+eig_zoh2_CL = eig(feedback(c2d(sys, Ts2, 'zoh'), k)) ;
+eig_zoh3_CL = eig(feedback(c2d(sys, Ts3, 'zoh'), k)) ;
+
+%%
+
 Question9 = figure;
 subplot(1,2,1) ;
+
 plot(eig_cont_CL, 'x');
+hold on;
 axis square;
 grid on;
 xline(0, 'r--');
 title({'Valeurs propres en temps continu'});
 xlabel('Réel');
 ylabel('Imaginaire');
+hold off;
 
 subplot(1,2,2);
-eig_zoh1 = eig(feedback(c2d(sys, Ts1, 'zoh'), k)) ;
-eig_zoh2 = eig(feedback(c2d(sys, Ts2, 'zoh'), k)) ;
-eig_zoh3 = eig(feedback(c2d(sys, Ts3, 'zoh'), k)) ;
-plot(eig_zoh1,'x', 'LineWidth',1);
+
+plot(eig_zoh1_CL,'x');
 hold on;
-plot(eig_zoh2,'x','LineWidth',1);
-plot(eig_zoh3,'x', 'LineWidth',1);
+plot(eig_zoh2_CL,'x');
+plot(eig_zoh3_CL,'x');
 axis square;
 plot(cos(theta), sin(theta), 'k--', 'HandleVisibility','off','LineWidth',0.25);
 xlim([-1.5, 1.5]);
 ylim([-1.5, 1.5]);
+
 grid on;
 legend("Ts1", "Ts2", "Ts3");
-
 title({'Valeurs propres en temps discret',' avec un bloqueur d''ordre 0'});
 xlabel('Réel');
 ylabel('Imaginaire');
@@ -193,18 +198,31 @@ if SAVE
 end
 
 %%
-sys_CL_Ts = feedback(c2d(sys, 0.037, 'zoh'),k);
-t_d = 0:0.037:10 ;
-cmd = double (t_d >= 1) ;
-[y, t_d, x] = lsim(sys_CL_Ts, H*cmd, t_d) ;
-u_d = H*cmd'-k*y ;
+
+for Ts = Ts2:0.001:Ts3
+    if all(abs(eig(feedback(c2d(sys,Ts,'zoh'),k))) < 1)
+        ;
+    else
+        Ts;
+        break;
+    end
+    
+end
+% Maximum Sampling period
+Ts = Ts-0.001
+%%
+sys_CL_Ts = feedback(c2d(sys, Ts, 'zoh'),k);
+t_d = 0:Ts:10 ;
+cmd_d = double (t_d >= 1) ;
+[y, t_d, x_d] = lsim(sys_CL_Ts, H*cmd_d, t_d) ;
+u_d = H*cmd_d'-k*y ;
 
 Question10 = figure ;
-plot(t_d, u_d, t_d, y, t_d, cmd, t_d, (Cw_L*x')');
+plot(t_d, u_d, t_d, y, t_d, cmd_d, t_d, (Cw_L*x_d')');
 xlabel("Temps (s)");
 ylabel("Amplitude");
 legend("$u(t)$", "$y(t)$", "$w_c(L,t)$", "$w(L,t)$", 'Interpreter', 'latex');
-title({'Simulation du retour de sortie avec pré-gain'},{'avec un bloqueur d''ordre 0'}) ;
+title({'Simulation du retour de sortie avec pré-gain','avec un bloqueur d''ordre 0'}) ;
 
 if SAVE
     exportgraphics(Question10, './latex/Illustrations/Question10.pdf')
@@ -213,6 +231,7 @@ end
 %% Question 11
 
 K_lqr = lqr(sys, eye(size(A)), 1);
+K_lqr
 H_lqr = inv(Cw_L*inv(-A+B*K_lqr)*B) ;
 
 %% Question 12
@@ -231,13 +250,71 @@ if SAVE
     exportgraphics(Question12, './latex/Illustrations/Question12.pdf')
 end
 
+%% Question 13
+
+for Ts = 0.01:0.0001:0.02
+    if all(abs(eig(c2d(sys,Ts,'zoh').A-c2d(sys,Ts,'zoh').B*K_lqr)) < 1)
+        ;
+    else
+        break;
+    end
+end
+Ts = Ts-0.0001
+%%
+sys_d = c2d(sys, Ts, 'zoh'); 
+Ad = sys_d.A; 
+Bd = sys_d.B; 
+Cd = sys_d.C; 
+sys_d_CL = ss(Ad-Bd*K_lqr, Bd*H_lqr, Cd, 0, Ts);
+
+t_d = 0:Ts:10 ;
+cmd_d = double(t_d > 1);
+[y, t_d, x_d] = lsim(sys_d_CL, cmd_d, t_d) ;
+u_d = H*cmd_d'-k*y ;
+
+Question13 = figure ;
+plot(t_d, u_d, t_d, y, t_d, cmd_d, t_d, (Cw_L*x_d')');
+xlabel("Temps (s)");
+ylabel("Amplitude");
+legend("$u(t)$", "$y(t)$", "$w_c(L,t)$", "$w(L,t)$", 'Interpreter', 'latex');
+title({'Simulation du retour d''états avec pré-gain','avec un bloqueur d''ordre 0'}) ;
+
+if SAVE
+    exportgraphics(Question13, './latex/Illustrations/Question13.pdf')
+end
+
+%% Question 14
+
+sys_d = c2d(sys, 0.01, 'zoh'); 
+Ad = sys_d.A; 
+Bd = sys_d.B; 
+Cd = sys_d.C; 
+sys_d_CL = ss(Ad-Bd*K_lqr, Bd*H_lqr, Cd, 0, 0.01);
+
+t_d = 0:0.01:10 ;
+cmd_d = double(t_d > 1);
+[y, t_d, x_d] = lsim(sys_d_CL, cmd_d, t_d) ;
+u_d = H*cmd_d'-k*y ;
+
+Question14 = figure ;
+plot(t_d, u_d, t_d, y, t_d, cmd_d, t_d, (Cw_L*x_d')');
+xlabel("Temps (s)");
+ylabel("Amplitude");
+legend("$u(t)$", "$y(t)$", "$w_c(L,t)$", "$w(L,t)$", 'Interpreter', 'latex');
+title({'Simulation du retour d''états avec pré-gain','avec un bloqueur d''ordre 0'}) ;
+
+if SAVE
+    exportgraphics(Question14, './latex/Illustrations/Question14.pdf')
+end
+
+
 %% Question 15
 Bp = [0; 0; 0; 0; 1/2*L^4-L^3+L; L^3-1/2*L^4; 1/4*L^4-2/3*L^3+1/2*L^2; 1/4*L^4-1/3*L^3];
 
 q0 = -5 * double (t' >= 4);
 sys_perturbations = ss(A-B*K_lqr, [B*H_lqr Bp], C, [0 0]);
 
-%%Question 16
+%% Question 16
 [y, t, x] = lsim(sys_perturbations, [cmd ; q0], t');
 u = H_lqr*cmd'-(K_lqr*x')' ;
 Question16 = figure ;
@@ -259,11 +336,13 @@ B_aug = [B; 0];
 B_p_aug = [Bp; 0];
 B_ref = [zeros(8,1); 1];
 
-%%
+%% Question 18
 K_aug = place(A_aug, B_aug, eig_4)
+
+%% Question 19
 sys_integral = ss(A_aug - B_aug*K_aug, [B_ref B_p_aug], [C, 0], [0 0]);
 
-[y, t, x] = lsim(sys_integral, [cmd ; q], t');
+[y, t, x] = lsim(sys_integral, [cmd ; q0], t');
 u = -(K_aug*x')' ;
 Question19 = figure ;
 plot(t, u, t, y, t, cmd, t, ([Cw_L 0]*x')');
@@ -277,7 +356,7 @@ end
 
 %%
 sys_integral_d = c2d(sys_integral, 0.5, 'zoh');
-t_d = 0:0.5:10; % Define the time vector for discrete simulation
+t_d = 0:0.5:10;
 cmd_d = double (t_d >= 1) ;
 q0_d = -5 * double(t_d >= 4);
 [y_d, t_d, x_d] = lsim(sys_integral_d, [cmd_d ; q0_d], t_d);