How can I fix my Runge Kutta (4th order) method to solve a 2nd order ODE?
solve x”(t) +δx'(t) + αx(t) + βx(t)^3 = γcos(ωt), x(0)=0 , x'(0)=0
clc
clear;
h=1;
delta = 0.1;
alpha = -1;
beta = 0.25;
gamma = 2.5;
omega = 2;
t(1)=0;
z(1)=0;
x(1)=0;
d=(gamma*(cos(omega*t)))-(beta*(x^3))-(alpha*x)-(delta*z)
g=@(t,x,z) z;
f=@(t,x,z) (gamma*(cos(omega*t)))-(beta*(x^3))-(alpha*x)-(delta*z);
for i=1:50
L1 = h*g(t, x, z);
k1 = h*f(t, x, z);
L2 = h*g(t+h/2, x+k1/2, z+L1/2);
k2 = h*f(t+h/2, x+k1/2, z+L1/2);
L3 = h*g(t+h/2, x+k2/2, z+L2/2);
k3 = h*f(t+h/2, x+k2/2, z+L2/2);
L4 = h*g(t+h, x+k3, z+L3);
k4 = h*f(t+h, x+k3, z+L3);
z = z + (L1+2*L2+2*L3+L4)/6;
x = x + (L1+2*L2+2*L3+L4)/6;
t= t+h;
fprintf(‘i=%8.0f, t=%8.2f, x=%8.6f, z=%8.6fn’,i,t,x,z)
plot(t,z,’-*r’)
hold on
plot(t,x,’-ob’)
endsolve x”(t) +δx'(t) + αx(t) + βx(t)^3 = γcos(ωt), x(0)=0 , x'(0)=0
clc
clear;
h=1;
delta = 0.1;
alpha = -1;
beta = 0.25;
gamma = 2.5;
omega = 2;
t(1)=0;
z(1)=0;
x(1)=0;
d=(gamma*(cos(omega*t)))-(beta*(x^3))-(alpha*x)-(delta*z)
g=@(t,x,z) z;
f=@(t,x,z) (gamma*(cos(omega*t)))-(beta*(x^3))-(alpha*x)-(delta*z);
for i=1:50
L1 = h*g(t, x, z);
k1 = h*f(t, x, z);
L2 = h*g(t+h/2, x+k1/2, z+L1/2);
k2 = h*f(t+h/2, x+k1/2, z+L1/2);
L3 = h*g(t+h/2, x+k2/2, z+L2/2);
k3 = h*f(t+h/2, x+k2/2, z+L2/2);
L4 = h*g(t+h, x+k3, z+L3);
k4 = h*f(t+h, x+k3, z+L3);
z = z + (L1+2*L2+2*L3+L4)/6;
x = x + (L1+2*L2+2*L3+L4)/6;
t= t+h;
fprintf(‘i=%8.0f, t=%8.2f, x=%8.6f, z=%8.6fn’,i,t,x,z)
plot(t,z,’-*r’)
hold on
plot(t,x,’-ob’)
end solve x”(t) +δx'(t) + αx(t) + βx(t)^3 = γcos(ωt), x(0)=0 , x'(0)=0
clc
clear;
h=1;
delta = 0.1;
alpha = -1;
beta = 0.25;
gamma = 2.5;
omega = 2;
t(1)=0;
z(1)=0;
x(1)=0;
d=(gamma*(cos(omega*t)))-(beta*(x^3))-(alpha*x)-(delta*z)
g=@(t,x,z) z;
f=@(t,x,z) (gamma*(cos(omega*t)))-(beta*(x^3))-(alpha*x)-(delta*z);
for i=1:50
L1 = h*g(t, x, z);
k1 = h*f(t, x, z);
L2 = h*g(t+h/2, x+k1/2, z+L1/2);
k2 = h*f(t+h/2, x+k1/2, z+L1/2);
L3 = h*g(t+h/2, x+k2/2, z+L2/2);
k3 = h*f(t+h/2, x+k2/2, z+L2/2);
L4 = h*g(t+h, x+k3, z+L3);
k4 = h*f(t+h, x+k3, z+L3);
z = z + (L1+2*L2+2*L3+L4)/6;
x = x + (L1+2*L2+2*L3+L4)/6;
t= t+h;
fprintf(‘i=%8.0f, t=%8.2f, x=%8.6f, z=%8.6fn’,i,t,x,z)
plot(t,z,’-*r’)
hold on
plot(t,x,’-ob’)
end matlab, differential equations MATLAB Answers — New Questions
