sample.m

clear 
%% Метод Эйлера %% 
% Условие задачи
U = dsolve('Du = (t - u) / (t - 2*u)', 'u(0) = 1'); 
tau0 = 0.03;
b = 0.3; 
% Функция для нахождения аналитического решения
u = @(t)eval(U); 
% Функция для нахождения приближенного решения
f = @(t,u) (t - u) / (t - 2*u);
T = b - 0;
 
%Метод Эйлера для случаев tau, tau/2, tau/4, tau/8, tau/16
% Порядок аппроксимации метода Эйлера
p = 1;           
for k = 0:4
    % Начальное условие
    y = 1;        
    tau = tau0/2^k; 
    N = T/tau; 
    t = 0;
    for m = 2:N+1 
        y = y + tau*f(t,y);
    	  t = t + tau;
end 
% Погрешность относительно точного решения u(b)
    D = abs(u(b)-y); 
if k > 0
       % Апостериорная оценка
        d = 1/(2^p-1)*abs(y - yprev);       
    end    
    yprev = y;
end 

%Метод Эйлера-Коши для случаев tau, tau/2, tau/4, tau/8, tau/16  
p = 2;
for k = 0:4 
    y = 1;        
    tau = tau0/2^k; 
    N = T/tau; 
    t = 0;
for m = 2:N+1 
        % Приближение во вспомогательной точке
        y_ = y + tau*f(t,y);         
  y = y + tau/2*(f(t,y)+f(t+tau,y_));
        t = t + tau;
    end 
D = abs(u(b)-y);    
if k > 0
        d = 1/(2^p-1)*abs(y - yprev);       
    end    
    yprev = y;
end

%Метод Рунге-Кутты для случаев tau, tau/2, tau/4, tau/8, tau/16
p = 4;           
for k = 0:4
    y = 1;        
    tau = tau0/2^k; 
    N = T/tau; 
    t = 0;
    for m = 2:N+1 
        k1 = tau*f(t,y);
        k2 = tau*f(t+0.5*tau, y+0.5*k1);
        k3 = tau*f(t+0.5*tau, y+0.5*k2);
        k4 = tau*f(t+tau, y+k3);
        y = y + 1/6*(k1+2*k2+2*k3+k4); 
        t = t + tau;
    end
    D = abs(u(b)-y);
    if k > 0
        d = 1/(2^p-1)*abs(y - yprev);
    end    
    yprev = y;
end