Evaluating a 2nd order ODE using the Runge-Kutta method

조회 수: 52 (최근 30일)
Emir Alp Karslioglu
Emir Alp Karslioglu 2020년 12월 31일
댓글: Emir Alp Karslioglu 2021년 1월 1일
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

채택된 답변

Alan Stevens
Alan Stevens 2020년 12월 31일
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
  댓글 수: 1
Emir Alp Karslioglu
Emir Alp Karslioglu 2021년 1월 1일
This code works like a charm! Thank you so much for your assistance. I need to work more on matlix manipulation in matlab it seems. Much appreciated!

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by