Source Code MATLAB Metode Adam Bashforth Moulton

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');
%fungsi 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:10;%batas interval
N=length(t)-1;
u1=zeros(N+1,1);
u2=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;
u1_eks(1)=4/3;
u2_eks(1)=2/3;
error1(1)=0;
error2(1)=0;
for i=1:N
    % Runge Kutta
    k1t=fu1(t(i),u1(i),u2(i));
    l1t=fu2(t(i),u1(i),u2(i));
    k2t=fu1(t(i)+h/2,u1(i)+k1t*h/2,u2(i)+l1t*h/2);
    l2t=fu2(t(i)+h/2,u1(i)+k1t*h/2,u2(i)+l1t*h/2);
    k3t=fu1(t(i)+h/2,u1(i)+k2t*h/2,u2(i)+l2t*h/2);
    l3t=fu2(t(i)+h/2,u1(i)+k2t*h/2,u2(i)+l2t*h/2);
    k4t=fu1(t(i)+h,u1(i)+k3t*h,u2(i)+l3t*h);
    l4t=fu2(t(i)+h,u1(i)+k3t*h,u2(i)+l3t*h);
 
    u1(i+1)=u1(i)+(k1t+2*k2t+2*k3t+k4t)*(h/6);
    u2(i+1)=u2(i)+(l1t+2*l2t+2*l3t+l4t)*(h/6);
end
for i=4:N
    % Adams-Bashforth -- *predictor*
    u1(i+1) = u1(i) + h/24*(55*fu1(t(i), u1(i),u2(i)) - 59*fu1(t(i-1),...
        u1(i-1),u2(i-1))+37*fu1(t(i-2),u1(i-2),u2(i-2)) - 9*fu1(t(i-3),...
        u1(i-3),u2(i-3)));
    u2(i+1) = u2(i) + h/24*(55*fu2(t(i), u1(i),u2(i)) - 59*fu2(t(i-1),...
        u1(i-1),u2(i-1))+37*fu2(t(i-2),u1(i-2),u2(i-2)) - 9*fu2(t(i-3),...
        u1(i-3),u2(i-3)));
t(i+1) = t(i) + h;
% Adams-Moulton -- *corrector*
u1(i+1) = u1(i) + h/24*(9*fu1(t(i+1),u1(i+1),u2(i+1)) + 19*fu1(t(i),...
        u1(i),u2(i))-5*fu1(t(i-1),u1(i-1),u2(i-1)) + fu1(t(i-2),u1(i-2),...
        u2(i-2)));
    u2(i+1) = u2(i) + h/24*(9*fu2(t(i+1),u1(i+1),u2(i+1)) + 19*fu2(t(i),...
        u1(i),u2(i))-5*fu2(t(i-1),u1(i-1),u2(i-1)) + fu2(t(i-2),u1(i-2),...
        u2(i-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));
end
%menampilkan error iterasi terakhir untuk masing-masing persamaan
fprintf('nilai u1 terakhirnya %d\n',u1(N));
fprintf('nilai u2 terakhirnya %d\n',u2(N));
fprintf('Error hasil iterasi terakhir u1 adalah %d\n',error1(N));
fprintf('Error hasil iterasi terakhir u2 adalah %d\n',error2(N));
figure(1)
plot(t,u1_eks,'ko',t,u1,'y--');
legend('u1 eksak','u1');
title('Metode ABM untuk u1');xlabel('t');ylabel('u1');
figure(2)
plot(t,u2_eks,'ro',t,u2,'k--');
legend('u2 eksak','u2');
title('Metode ABM untuk u2');xlabel('t');ylabel('u2');
figure(3)
plot(t,error1,'go',t,error2,'k*');
xlabel('t');ylabel('galat u1');
legend('Galat u1','Galat u2');

Source code MATLAB metode Euler Eksplisit

%Metode Euler Eksplisit
%dy/dt=f(t,y)
%f(t,y)=-(y+1)*(y+3);
clear all;
dt=0.01;
t=0:dt:2;
N=length(t)-1;
y=zeros(N+1,1);

y(1)=-2;
y_eks(1)=-2;
err(1)=0;
for n=1:N
    y(n+1)=y(n)+dt*(-(y(n)+1)*(y(n)+3));%solusi Euler eksplisit
    y_eks(n+1)=-3+2/(1+exp(-2*t(n+1)));%solusi eksak
    err(n+1)=abs(y(n+1)-y_eks(n+1));%kesalahan mutlak
end;
figure(1)
plot(t,y,'ko',t,y_eks,'r');
figure(2);%hold on
plot(t,err);

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');