Evaluating a 2nd order ODE using the Runge-Kutta method
47 views (last 30 days)
Show older comments
Emir Alp Karslioglu
on 31 Dec 2020
Commented: Emir Alp Karslioglu
on 1 Jan 2021
Greetings,
I've been working on a 2nd order ODE: y''(t) = -e^(3t)*y'(t) - y(t) + (5-2e^(-3t))*e^(-2t) +1
With initial conditions y(0) = 2 ; y'(0) = -2 where 0<t<1
I need to use a third order Runge Kutta method, which I can code for a 1st order ODE.
However given a 2nd order differential equation, I'm having difficulties implementing the ODE into my Runge Kutta code.
I've tried writing the 2nd order ODE in a linear system of 1st order ODEs but im still stuck.
By the way, i do not want to use ode45 or ode23 commands. Thank you for your feedback!
Here is what I've been working on so far:
function RK3 = func(a,b,n)
a = 0;
b = 1;
n = input('Enter the number of intervals:');
h = (b-a)/n;
y(1) = 2;
y(2) = -2;
for i=1:n
t = (i-1)*h;
k1 = f(t,y(:,i));
k2 = f(t+0.5*h,y(:,i)+0.5*k1*h);
k3 = f(t+h,y(:,i)-k1*h+(2*k2*h));
y(i+1)= y(:,i)+h/6*(k1+4*k2+k3); %approximated value of the ODE using RK3
end
for k=1:n+1
t=a+(k-1)*h;
exact(k)=exp(-2*t)+1; %exact value of the ODE
end
error = max(abs(y-exact));
max(exact)
max(y)
error
end
function fty = f(t,y) %matrix form of the 2nd order ODE
y = [0 1; -1 -exp(-3*t)];
u = [0; 5-2*exp(-3*t)*exp(-2*t)+1];
fty = y + u;
end
0 Comments
Accepted Answer
Alan Stevens
on 31 Dec 2020
More like this:
a = 0;
b = 1;
%n = input('Enter the number of intervals:');
n = 100;
h = (b-a)/n;
y = [2; -2];
t = 0:h:n*h;
for i=1:n
k1 = f(t(i),y(:,i));
k2 = f(t(i)+0.5*h,y(:,i)+0.5*k1*h);
k3 = f(t(i)+h,y(:,i)-k1*h+2*k2*h);
y(:,i+1)= y(:,i)+h/6*(k1+4*k2+k3); %approximated value of the ODE using RK3
end
exact = exp(-2*t)+1;
disp(max(abs(exact-y(1,:))))
plot(t,y(1,:),'r',t,exact,'b--'),grid
xlabel('t'),ylabel('y')
legend('RK3','Exact')
function fty = f(t,y) %matrix form of the 2nd order ODE
Y = [0 1; -1 -exp(-3*t)]*y;
u = [0; (5-2*exp(-3*t))*exp(-2*t)+1];
fty = Y + u;
end
More Answers (0)
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!