Severity: Warning
Message: fopen(/home/polpe/.phpsession/ci_session3a98b9fc6bf83d5afb3150e8831fd89c9d61ce23): failed to open stream: No space left on device
Filename: drivers/Session_files_driver.php
Line Number: 159
Backtrace:
File: /home/polpe/public_html/application/controllers/Main.php
Line: 17
Function: library
File: /home/polpe/public_html/index.php
Line: 315
Function: require_once
Teljes Matlab script kiegészítő függvényekkel.
file: ipend.m author: Peter Polcz <ppolcz@gmail.com>
Created on 2017.03.27. Monday, 17:42:37 Reviewed on 2017. October 30. [Published] Reviewed on 2017. November 08. [Corrected and published]
try c = evalin('caller','persist'); catch; c = []; end
persist = pcz_persist(mfilename('fullpath'), c); clear c;
persist.backup();
%clear persist
global PUBLISH
PUBLISH = 1;
│ - Persistence for `ipend1` reused (inherited) [run ID: 48, 036]Persistence for `2017.11.09. Thursday, 11:25:23` │ - Script `ipend1` backuped
% model parameters
M = 0.5;
m = 0.2;
l = 1;
g = 9.8;
b = 0;
A = [
0 1 0 0
0 -(4*b)/(4*M + m) -(3*g*m)/(4*M + m) 0
0 0 0 1
0 -(3*b)/(l*(4*M + m)) -(3*g*(M + m))/(l*(4*M + m)) 0
];
B = [
0
4/(m+4*M)
0
3/(m+4*M)/l
];
C = [
1 0 0 0
0 0 1 0
];
D = [ 0 ; 0 ];
sys = ss(A,B,C,D);
H = tf(sys);
H_r = H(1)
H_phi = H(2)
H_r = 1.818 s^2 + 13.36 ----------------- s^4 + 9.355 s^2 Continuous-time transfer function. H_phi = 1.364 ----------- s^2 + 9.355 Continuous-time transfer function.
T = 5;
[y,t] = impulse(H,T);
ipend_simulate_Pi(t,y)
T = 8;
figure('Position', [ 226 318 1276 392 ], 'Color', [1 1 1])
subplot(121)
step(H_r,T)
subplot(122)
step(H_phi,T), hold on
DC_gain = dcgain(H_phi);
plot([0 T],[1 1]*DC_gain)
Simulate
[y,t] = step(H,T);
ipend_simulate_Pi(t,y)
Is the system stable, exponentially stable or unstable
eig(A)
ans = 0.0000 + 0.0000i 0.0000 + 0.0000i 0.0000 + 3.0585i 0.0000 - 3.0585i
Plot only for the transfer function from F to omega. The frequency unit set to be in Hz. Determine the own (or resonance) frequency (f_r) of the system.
figure
bopts = bodeoptions;
bopts.FreqUnits = 'Hz';
% bopts.FreqScale = 'linear';
% bopts.MagUnits = 'abs';
bodeplot(H_phi, bopts), grid on
[gpeak,fpeak] = getPeakGain(H_phi);
fr = fpeak / (2*pi)
fr = 0.4868
figure(6), nyquist(H_phi)
T = 20;
Ts = 0.01;
t = 0:Ts:T;
f = fr;
Amp = 1;
u = Amp*sin(2*pi*f*t);
y = lsim(H,u,t);
ipend_simulate_Pi(t,y,u);
% After a while we shall notice, that the system's motion is quite
% unusual, why?
T = 10;
Ts = 0.01;
t = 0:Ts:T;
f = 4;
Amp = 20;
u = Amp*sin(2*pi*f*t)-0.1;
y = lsim(H,u,t);
ipend_simulate_Pi(t,y,u);
T = 20;
Ts = 0.01;
t = 0:Ts:T;
f = 0.1;
Amp = 1;
u = Amp*sin(2*pi*f*t + pi/2);
y = lsim(H,u,t);
ipend_simulate_Pi(t,y,u);
x0 = [0 0 pi/6 0]';
u = 0;
T = 5;
[t,x] = ode45(@(t,x) A*x + B*u, [0 T], x0);
ipend_simulate_Pi(t,x)
C4 = [ B A*B A*A*B A*A*A*B ]
C4_ = ctrb(A,B);
rank(C4)
C4 = 0 1.8182 0 -3.6446 1.8182 0 -3.6446 0 0 1.3636 0 -12.7562 1.3636 0 -12.7562 0 ans = 4
C = [
1 0 0 0
0 0 1 0
];
O4 = obsv(A,C)
O4_ = ctrb(A',C')'
rank(O4)
% Measure only phi
C = [ 0 0 1 0 ];
O4 = obsv(A,C);
S = [ orth(O4) null(O4) ];
T = inv(S);
A_ = T*A/T
B_ = T*B
C_ = C/T
O4 = 1.0000 0 0 0 0 0 1.0000 0 0 1.0000 0 0 0 0 0 1.0000 0 0 -2.6727 0 0 0 -9.3545 0 0 0 0 -2.6727 0 0 0 -9.3545 O4_ = 1.0000 0 0 0 0 0 1.0000 0 0 1.0000 0 0 0 0 0 1.0000 0 0 -2.6727 0 0 0 -9.3545 0 0 0 0 -2.6727 0 0 0 -9.3545 ans = 4 A_ = 0 -1.0000 0 0 9.3545 0 0 0 0 0 0 -1.0000 3.6519 0 0 0 B_ = 0 1.3714 0 1.9640 C_ = -0.9943 0 0 0
$$\dot x = f(x) + g(x) u$$
syms t x v phi omega
% model parameters
M = 0.5;
m = 0.2;
l = 1;
g = 9.8;
b = 0;
q = 4*(M+m) - 3*m*cos(phi)^2;
f_sym = [
v
(4*m*l*sin(phi)*omega^2 - 1.5*m*g*sin(2*phi) -4*b*v) / q
omega
3*(-m*l*sin(2*phi)*omega^2 / 2 + (M+m)*g*sin(phi) + b*cos(phi)*v) / (l*q)
];
g_sym = [
0
4*l
0
-3*cos(phi)
] / (l*q);
f_ode = matlabFunction(f_sym, 'vars', {t , [ x ; v ; phi ; omega ]});
g_ode = matlabFunction(g_sym, 'vars', {t , [ x ; v ; phi ; omega ]});
[tt,xx] = ode45(f_ode, [0,10], [0 0 0.1 0]');
ipend_simulate_0(tt,xx)
[tt,xx] = ode45(f_ode, [0,10], [0 0 0.1 2]');
ipend_simulate_0(tt,xx)
Ampl = 20;
f = 4;
u = @(t) Ampl*sin(2*pi*f*t);
[tt,xx] = ode45(@(t,x) f_ode(t,x) + g_ode(t,x)*u(t), [0,5], [0 0 pi 0]');
ipend_simulate_0(tt,xx,u)
Ampl = 1;
f = fr;
u = @(t) Ampl*sin(2*pi*f*t);
[tt,xx] = ode45(@(t,x) f_ode(t,x) + g_ode(t,x)*u(t), [0,30], [0 0 pi 0]');
ipend_simulate_0(tt,xx,u)
Ampl = 1;
f = 0.1;
u = @(t) Ampl*sin(2*pi*f*t + pi/2);
[tt,xx] = ode45(@(t,x) f_ode(t,x) + g_ode(t,x)*u(t), [0,30], [0 0 pi 0]');
ipend_simulate_0(tt,xx,u)
$$\dot x = f(x) + g(x) u$$
syms t x v phi omega
% model parameters
M = 0.5;
m = 0.2;
l = 1;
g = 9.8;
b = 10;
q = 4*(M+m) - 3*m*cos(phi)^2;
f_sym = [
v
(4*m*l*sin(phi)*omega^2 - 1.5*m*g*sin(2*phi) -4*b*v) / q
omega
3*(-m*l*sin(2*phi)*omega^2 / 2 + (M+m)*g*sin(phi) + b*cos(phi)*v) / (l*q)
];
g_sym = [
0
4*l
0
-3*cos(phi)
] / (l*q);
f_ode = matlabFunction(f_sym, 'vars', {t , [ x ; v ; phi ; omega ]});
[t,x] = ode45(f_ode, [0,10], [0 0 pi/12 2]');
ipend_simulate_0(t,x)
model parameters
M = 0.5;
m = 0.2;
l = 1;
g = 9.8;
b = 0;
A = [
0 1 0 0
0 -(4*b)/(4*M + m) -(3*g*m)/(4*M + m) 0
0 0 0 1
0 -(3*b)/(l*(4*M + m)) -(3*g*(M + m))/(l*(4*M + m)) 0
];
B = [
0
4/(m+4*M)
0
3/(m+4*M)/l
];
C = [
1 0 0 0
0 0 1 0
];
D = [ 0 ; 0 ];
sys = ss(A,B,C,D);
H = tf(sys);
H_r = H(1);
H_phi = H(2);
Only if $b \not= 0$!
open_system ipend_PID_demo
print -dpng fig/ipend_PID.png -sipend_PID_demo
C_r = pidtune(H_r,'pid');
% pidTuner(H_r,C_r)
At home: derive the resulting transfer function
Control_Loop = [
C_r*H_r / (1 + C_r*H_r)
C_r*H_phi / (1 + C_r*H_r)
];
T = 20;
t = linspace(0,T,T/0.1);
u = t*0 + 5;
y = lsim(Control_Loop,u,t);
ipend_simulate_Pi(t,y)