Runge-kutta third order method:
%rk3:runge kutta of thirdorder
clc;
clear all;
close all;
% y' = y-x ode condition
f = @(x,y) y-x;
fex = @(x) exp(x)+x+1; % exact solution
a=0;
b= 3.2;
n =16;
h=(b-a)/n;
y(1) =2; %initial value
i = 0;
for x= a:h:b
i = i+1;
K1 = f(x,y(i)); %initializing solution
K2 = f(x+h*0.5,y(i)+h*K1*0.5);
K3 = f(x+h, y(i)-h*K1 +2*K2*h);
y(i+1) =y(i)+h*(1/6)*(K1 +4*K2+K3);
g(i) = fex(x);
xx(i) = x;
Error(i) = abs(g(i) - y(i)); %error obtain
end
%plot result
plot(xx,y(1:n+1),'k',xx,g,'y')
legend('RK3','Exact solution')
xlabel('x')
ylabel('y')
title('RK3 vs exact solution')
I am getting wrong error value.. please check my code

댓글 수: 4

KSSV
KSSV 2019년 5월 15일
Wrong error value????
NOTE: initialze the variable inside the loop.
SHIVANI TIWARI
SHIVANI TIWARI 2019년 5월 15일
am not getting the right value of absolute error.
i have taken same problem 2. but am not getting the right value of absolute error.
and plz also my program for butcher 5 for solving 1st order ode
clc;
clear all;
close all;
% y' = y-x ode condition
f = @(x,y) y-x;
fex = @(x) exp(x)+x+1; % exact solution
a=0;
b= 3.2;
n =16;
h=(b-a)/n;
y(1) =2; %initial value
i = 0;
for x= a:h:b
i = i+1;
K1 = f(x,y(i));
K2 = f(x+h*(1/4),y(i)+h*(1/4)*K1);
K3 = f(x+h*0.5, y(i)+ h*K1*(1/8)+h*K2*(1/8));
K4 = f(x+h*0.5,y(i)- h*K2*0.5+h*K3);
K5 = f(x+3*h*(1/4),y(i)+ (3/16)*h*K1+(9/16)*h*K4);
K6 = f(x+h,y(i)- 3*h*K1*(1/7)+2*h*K2*(1/7) + (12/7)*h*K3 -(12/7)*h*K4 + (8/7)*h*K5);
y(i+1) =y(i)+h*(1/90)*(7*K1 +32*K3+12*K4+ 32*K5 +7*K6);
g(i) = fex(x);
xx(i) = x;
Error(i) = abs(g(i) - y(i)); %error obtain
%plot result
plot(xx,y(1:n+1),'r',xx,g,'g','LineWidth',2)
legend('butcher5','Exact solution')
xlabel('x')
ylabel('y')
title('butcher5 vs exact solution')
Ana Paula Cervantes Martínez
Ana Paula Cervantes Martínez 2021년 4월 30일
excuse me, how do you determine those formulas for finding the values for K's?
SHIVAM
SHIVAM 2023년 11월 6일
thank u so much

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

답변 (0개)

카테고리

도움말 센터File Exchange에서 Numerical Integration and Differential Equations에 대해 자세히 알아보기

질문:

2019년 5월 15일

댓글:

2023년 11월 6일

Community Treasure Hunt

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

Start Hunting!

Translated by