solving odes with changing input

13 次查看(过去 30 天)
samee baig
samee baig 2016-8-21
i have a system which need to change its input after reaching certain values and do this until all nine states go to zero.Any help/guidance would be greatly appreciated
my code is as below
function xdot=njc(t,x)
k1=1; k2=1.5; k3=1.5; b1=0.5; b2=0.6; b3=0.75; flag1=1; %flag2=-1;
%Step 1. Drive the system to a nonzero equilibrium such that tan(x4) =/ 0 in finite time using
if flag1==1
uu(1)=0;
uu(2)=-k1*(abs(x(4)-4)).^b1*sign(x(4)-4)-k2*(abs(x(5))).^b2*sign(x(5))-k3*(abs(x(6))).^b3*sign(x(6));
uu=[uu(1);uu(2)]; end %Step 2. Drive the x7 variable to zero while keeping x4 = 4 using
if x(4)==4
flag1=-1;
end
if flag1==-1 uu(1)=-1/tan(4)*(k1*(abs(x(7))).^b1*sign(x(7))+k2*(abs(x(8))).^b2*sign(x(8))+k3*(abs(x(9))).^b3*sign(x(9))); uu(2)=-k1*(abs(x(4)-4)).^b1*sign(x(4)-4)-k2*(abs(x(5))).^b2*sign(x(5))-k3*(abs(x(6))).^b3*sign(x(6)); uu=[uu(1);uu(2)]; end
% % Step 3. Drive the x4 variable to zero while keeping x7 = 0 using
if x(7)==0
flag1=0;
end
%
if flag1==0
uu(1)=-k1*(abs(x(7))).^b1*sign(x(7))-k2*(abs(x(8))).^b2*sign(x(8))-k3*(abs(x(9))).^b3*sign(x(9));
uu(2)=-k1*(abs(x(4))).^b1*sign(x(4))-k2*(abs(x(5))).^b2*sign(x(5))-k3*(abs(x(6))).^b3*sign(x(6));
uu=[uu(1);uu(2)];
end
% % Step 4. Drive the x1 variable to zero while keeping x4 = 0 using
if x(4)==0 flag1=2; end
if flag1==2 uu(1)=-k1*(abs(x(1))).^b1*sign(x(1))-k2*(abs(x(2))).^b2*sign(x(2))-k3*(abs(x(3))).^b3*sign(x(3)); uu(2)=-k1*(abs(x(4))).^b1*sign(x(4))-k2*(abs(x(5))).^b2*sign(x(5))-k3*(abs(x(6))).^b3*sign(x(6)); uu=[uu(1);uu(2)]; end
xdot(1)=x(2); xdot(2)=x(3); xdot(3)=uu(1);
xdot(4)=x(5); xdot(5)=x(6); xdot(6)=uu(2);
xdot(7)=x(8); xdot(8)=x(9); xdot(9)=uu(1)*tan(x(4));
xdot=xdot'; end
clear all close all clc
Tspan =linspace(0,60,250);
IC=[1;0;0;pi/4;0;0;-1;0;0]; %opts = odeset('InitialStep',1e-5,'events',@g); %opts = odeset('InitialStep',1e-15); %options = odeset('RelTol',1e-3,'AbsTol',[1e-3 1e-3 1e-3 1e-3 1e-3 1e-3 1e-3 1e-3 1e-3]); [T,Y] = ode45(@njc, Tspan,IC)
x1 = Y(:,1); x2 = Y(:,2); x3 = Y(:,3); x4 = Y(:,4); x5 = Y(:,5); x6 = Y(:,6); x7 = Y(:,7); x8 = Y(:,8); x9 = Y(:,9);
% if x4==4 % i=3 % end
figure(1) plot(T,x1,T,x2,T,x3,'linewidth',2) legend('x1','x2','x3',3); xlabel('t(s)') ylabel('x_1,x_2,x_3') grid on
figure(2) plot(T,x4,T,x5,T,x6,'linewidth',2) legend('x4','x5','x6',3); xlabel('t(s)') ylabel('x_4,x_5,x_6') grid on
figure(3) plot(T,x7,'+',T,x8,'o',T,x9,'linewidth',2) legend('x7','x8','x9',3); xlabel('t(s)') ylabel('x_7,x_8,x_9') grid on

回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Programming 的更多信息

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by