15,961,997 members

See more:

I have impulsive differential equation. I want to code them in octave.

**What I have tried:**

clc; clear all; close all; % step sizes t0 = 0; tfinal = 5; h = 0.01; n = ceil((tfinal-t0)/h)+1; %initial conditions s(1) = 0.5; theta = 0.01; N0 = 2; Ta = 0.7; b=N0*(Ta-theta)+theta; x1(1) = -5; x2(1) = 3; t(1) =0; % functions handles Fx1 = @(t,x1,x2) -1/2*x1-x1.^3+x1.*x2.^2-((x1.^2+x2.^2).^1/4).*sign(x1); Fx2 = @(t,x1,x2) -1/2*x2-x2.^3-x2.*x1.^2-((x1.^2+x2.^2).^1/4).*sign(x2); for N = 1:n if (mod(N,N0)==0) s(N+1) = s(N)+b; else s(N+1) = theta; end for i = 1:ceil(s(N+1)) t(i+1) = t(i)+h; ## if (mod(i,N0)==0) ## h = b; ## t(i+1) = t(i)+h; ## else ## h = theta; ## t(i+1) = t(i)+h; ## end %updates x1 and x2 k1x1 = Fx1(t(i), x1(i), x2(i) ); k1x2 = Fx2(t(i), x1(i), x2(i) ); k2x1 = Fx1(t(i)+h/2, x1(i)+h/2*k1x1, x2(i)+h/2*k1x2); k2x2 = Fx2(t(i)+h/2, x1(i)+h/2*k1x1, x2(i)+h/2*k1x2); k3x1 = Fx1(t(i)+h/2, x1(i)+h/2*k2x1, x2(i)+h/2*k2x2); k3x2 = Fx2(t(i)+h/2, x1(i)+h/2*k2x1, x2(i)+h/2*k2x2); k4x1 = Fx1(t(i)+h, x1(i)+h*k3x1, x2(i)+h*k3x2); k4x2 = Fx2(t(i)+h, x1(i)+h*k3x1, x2(i)+h*k3x2); if (i== ceil(s(N+1))) x1(i+1) = 0.5*x1(i)+h/6*(k1x1 + 2*k2x1 + 2*k3x1 +k4x1); x2(i+1) = 0.8*x2(i)+h/6*(k1x2 + 2*k2x2 + 2*k3x2 +k4x2); else x1(i+1) = x1(i)+h/6*(k1x1 + 2*k2x1 + 2*k3x1 +k4x1); x2(i+1) = x2(i)+h/6*(k1x2 + 2*k2x2 + 2*k3x2 +k4x2); end end N = N+1; end %plot plot(t,x1); hold on; plot(t,x2); ##xlim([0, tfinal]); ##ylim([-6, 6]); drawnow();

This content, along with any associated source code and files, is licensed under The Code Project Open License (CPOL)

CodeProject,
20 Bay Street, 11th Floor Toronto, Ontario, Canada M5J 2N8
+1 (416) 849-8900