Source Code MATLAB metode Runge-Kutta Order 2 dan 4

clc;

clear all;

%fungsi yang akan dicari penyelesaiannya

fu1=inline('9*u1+24*u2+5*cos(t)-(1/3)*sin(t)','t','u1','u2');

fu2=inline('-24*u1-51*u2-9*cos(t)+(1/3)*sin(t)','t','u1','u2');

%penyelesaian secara eksak

fu1_eks=inline('2*exp(-3*t)-exp(-39*t)+(1/3)*cos(t)','t');

fu2_eks=inline('-exp(-3*t)+2*exp(-39*t)-(1/3)*cos(t)','t');

h=0.01;%lebar sub interval

t=0:h:20;%batas interval

N=length(t)-1;

u1=zeros(N+1,1);

u2=zeros(N+1,1);

u3=zeros(N+1,1);

u4=zeros(N+1,1);

u1_eks=zeros(N+1,1);

u2_eks=zeros(N+1,1);

%nilai awal

u1(1)=4/3;

u2(1)=2/3;

u3(1)=4/3;

u4(1)=2/3;

u1_eks(1)=4/3;

u2_eks(1)=2/3;

u3_eks(1)=4/3;

u4_eks(1)=2/3;

error1(1)=0;

error2(1)=0;

for i=1:N

    %Metode Runge Kutta Order 2

    k1=h*fu1(t(i),u1(i),u2(i));

    l1=h*fu2(t(i),u1(i),u2(i));

    k2=h*fu1(t(i)+h,u1(i)+k1,u2(i)+l1);

    l2=h*fu2(t(i)+h,u1(i)+k1,u2(i)+l1);

    %solusi numerik

    u1(i+1)=u1(i)+(k1+k2)/2;

    u2(i+1)=u2(i)+(l1+l2)/2;

    %solusi eksak

    u1_eks(i+1)=fu1_eks(t(i+1));

    u2_eks(i+1)=fu2_eks(t(i+1));

    %error

    error1(i+1)=abs(u1(i+1)-u1_eks(i+1));

    error2(i+1)=abs(u2(i+1)-u2_eks(i+1));

    %Metode Runge Kutta Order 4

    k1t=fu1(t(i),u3(i),u4(i));

    l1t=fu2(t(i),u3(i),u4(i));

    k2t=fu1(t(i)+h/2,u3(i)+k1t*h/2,u4(i)+l1t*h/2);

    l2t=fu2(t(i)+h/2,u3(i)+k1t*h/2,u4(i)+l1t*h/2);

    k3t=fu1(t(i)+h/2,u3(i)+k2t*h/2,u4(i)+l2t*h/2);

    l3t=fu2(t(i)+h/2,u3(i)+k2t*h/2,u4(i)+l2t*h/2);

    k4t=fu1(t(i)+h,u3(i)+k3t*h,u4(i)+l3t*h);

    l4t=fu2(t(i)+h,u3(i)+k3t*h,u4(i)+l3t*h);

    %solusi numerik

    u3(i+1)=u3(i)+(k1t+2*k2t+2*k3t+k4t)*(h/6);

    u4(i+1)=u4(i)+(l1t+2*l2t+2*l3t+l4t)*(h/6);

    %solusi eksak

    u3_eks(i+1)=fu1_eks(t(i+1));

    u4_eks(i+1)=fu2_eks(t(i+1));

    %error

    error3(i+1)=abs(u3(i+1)-u3_eks(i+1));

    error4(i+1)=abs(u4(i+1)-u4_eks(i+1));

end

%menampilkan error iterasi terakhir untuk masing-masing persamaan

fprintf('Error hasil iterasi terakhir u1 metode RK2 adalah %d\n',error1(N));

fprintf('Error hasil iterasi terakhir u1 metode RK4 adalah %d\n',error3(N));

fprintf('Error hasil iterasi terakhir u2 metode RK2 adalah %d\n',error2(N));

fprintf('Error hasil iterasi terakhir u2 metode RK2 adalah %d\n',error4(N));

figure(1)

plot(t,u1_eks,'ko',t,u1,'b*',t,u3,'g--');

legend('u1 eksak','u1 metode RK2','u1 metode RK4');

title('Metode Runge-Kutta untuk u1');xlabel('t');ylabel('u1');

figure(2)

plot(t,u2_eks,'ro',t,u2,'k*',t,u4,'y--');

legend('u2 eksak','u2 metode RK2','u2 metode RK4');

title('Metode Runga-Kutta untuk u2');xlabel('t');ylabel('u2');

figure(3)

plot(t,error1,'go',t,error3,'k*');

xlabel('t');ylabel('galat u1');

legend('Galat u1 metode RK2','Galat u1 metode RK4');

figure(4)

plot(t,error2,'ro',t,error4,'b*');

xlabel('t');ylabel(' galat u2');

legend('Galat u2 metode RK2','Galat u2 metode RK4');