Ordinary Differential Equation

1. Euler’s Method

% Define the differential equation dy/dx = x + y

f = @(x, y) x + y;

% Initial conditions and step size

x0 = 0;      % initial value of x

y0 = 1;      % initial value of y

h = 0.1;     % step size

xn = 1;      % final value of x

% Create x values from x0 to xn with step size h

x = x0:h:xn;

% Initialize y array with zeros

y = zeros(size(x));

% Assign initial condition

y(1) = y0;

% Euler method iteration

for i = 1:length(x)-1

    % Apply Euler formula: y(i+1) = y(i) + h*f(x,y)

    y(i+1) = y(i) + h * f(x(i), y(i));

end

% Plot the numerical solution

plot(x, y, '-o')

title('Euler Method')

grid on

2. Runge-Kutta 2nd Order (RK2)
% Define the differential equation
f = @(x, y) x + y;

% Initial conditions
x0 = 0;
y0 = 1;
h = 0.1;
xn = 1;

% Generate x values
x = x0:h:xn;

% Initialize y array
y = zeros(size(x));
y(1) = y0;

% RK2 iteration
for i = 1:length(x)-1
    
    % First slope (at beginning of interval)
    k1 = f(x(i), y(i));
    
    % Second slope (at midpoint)
    k2 = f(x(i) + h/2, y(i) + h*k1/2);
    
    % Update y using RK2 formula
    y(i+1) = y(i) + h * k2;
end

% Plot the result
plot(x, y, '-o')
title('RK2 Method')
grid on






3. ODE using ode45 (Single Variable)
% Define the differential equation
f = @(x, y) x + y;

% Solve ODE using ode45 solver
% [x, y] gives solution points and values
[x, y] = ode45(f, [0 1], 1);

% Plot the solution
plot(x, y)
title('ode45 Solution')
grid on

% Simpson's 1/3 Rule for Numerical Integration
clc;
clear;

% Input data points
x = [0 1 2 3 4];              % x values
y = [1 2 4 8 16];             % y values

n = length(x);

% Check if number of intervals is even
if mod(n-1,2) ~= 0
    error('Number of intervals must be even for Simpson''s 1/3 rule');
end

h = x(2) - x(1);              % step size

sum_odd = 0;                  % sum of odd indexed terms
sum_even = 0;                 % sum of even indexed terms

% Loop for summation
for i = 2:n-1
    if mod(i,2) == 0
        sum_odd = sum_odd + y(i);     % odd terms
    else
        sum_even = sum_even + y(i);   % even terms
    end
end

% Simpson's formula
I = (h/3) * (y(1) + y(n) + 4*sum_odd + 2*sum_even);

% Display result
fprintf('The integral using Simpson''s Rule = %.4f\n', I);

% Numerical Differentiation using Finite Difference Methods
clc;
clear;

% Input data points
x = [1 2 3 4 5];
y = [1 4 9 16 25];    % Example: y = x^2

h = x(2) - x(1);      % step size

% Forward Difference (at first point)
dy_forward = (y(2) - y(1)) / h;

% Backward Difference (at last point)
dy_backward = (y(end) - y(end-1)) / h;

% Central Difference (at middle point)
i = 3; % middle point index
dy_central = (y(i+1) - y(i-1)) / (2*h);

% Display results
fprintf('Forward Difference = %.4f\n', dy_forward);
fprintf('Backward Difference = %.4f\n', dy_backward);
fprintf('Central Difference = %.4f\n', dy_central);