-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathEX11111.m
79 lines (57 loc) · 1.96 KB
/
EX11111.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
%----------------------------------------------------------------------
%
% This program {\it compen.m} demonstrates a dynamic compensator
% design and diplay the simulation result. A finte beam element is
% adopted as a system. For obeserver gain design, the LQR technique
% is used, and the dynamic siumlation is performed using a MATLAB
% control toolbox routine {\it lsim.m}.
%
%---------------------------------------------------------------------
clear
% Provide the system mass and stiffness matrix
M =[0.5571 0 0.0964 -0.3482
0 3.2143 0.3482 -1.2054
0.0964 0.3482 0.2786 -0.5893
-0.3482 -1.2054 -0.5893 1.6071];
K =[2.4178 0 -1.2089 9.0667
0 181.3333 -9.0667 45.3333
-1.2089 -9.0667 1.2089 -9.0667
9.0667 45.3333 -9.0667 90.6667];
F=[1;0;0;0];
% Transform into the first order state space form equation
A=[0*eye(4),eye(4);-inv(M)*K,0*eye(4)];
B=[0*ones(4,1);inv(M)*F];
C=[1 0 0 0 0 0 0 0];
% Use the {\it flqr.m} function to design the observer gain and the
% full state feedback gain
[G,Sc]=felqr(A, B, 1000*eye(8), 0.01);
[L,So]=felqr(A',C', eye(8),0.01);
% Now build the total closed-loop system for both closed-loop system and observer
Atot=[A, -B*G;L'*C, A-L'*C-B*G];
Btot=[B;B];
Ctot=eye(16);
Dtot=eye(16,1);
% Define simulation time, control input, and initial conditions
t=0:0.03:6.0-0.03;
u=zeros(200,1);
x0=zeros(16,1);
x0(1,1)=0.1; x0(2,1)=-0.3; x0(3,1)=0.2;
% Use MATLAB function {\it felresp.m} to simulate the total system
[x,y]=felresp(Atot,Btot,Ctot,Dtot,x0,u,t);
% Plot the result
subplot(221)
plot(t,y(:,1),'-',t,y(:,9),':');
xlabel('Time(sec)');
ylabel('Displacement');
subplot(222)
plot(t,y(:,2),'-',t,y(:,10),':');
xlabel('Time(sec)');
ylabel('Rotation');
subplot(223)
plot(t,y(:,5),'-',t,y(:,13),':');
xlabel('Time(sec)');
ylabel('Linear velocity');
subplot(224)
plot(t,y(:,6),'-',t,y(:,14),':');
xlabel('Time(sec)');
ylabel('Angular velocity');