When you have a script and functions in the same file, the script must be first.
x = linspace(0, 1, 100);
y = rk4(@(x, y) sin(x) + 5*y, 0, 1, 0.01, x);
plot_solution(y, x)
function [y] = rk4(f, x0, y0, h, x)
%RK4 Solves the ODE y'=f(x,y) with initial condition y(x0)=y0 using the
%fourth order Runge-Kutta method.
% Inputs:
% f - Function handle for the right-hand side of the ODE.
% x0 - Initial value of x.
% y0 - Initial value of y.
% h - Step size.
% x - Array of values of x at which to evaluate the solution.
% Outputs:
% y - Array of values of the solution at the specified values of x.
% Initialize the solution.
y = zeros(size(x));
y(1) = y0;
% Loop over the values of x.
for i=2:numel(x)
% Calculate the k1, k2, k3, and k4 coefficients.
k1 = h*f(x(i-1), y(i-1));
k2 = h*f(x(i-1) + h/2, y(i-1) + k1/2);
k3 = h*f(x(i-1) + h/2, y(i-1) + k2/2);
k4 = h*f(x(i-1) + h, y(i-1) + k3);
% Update the solution.
y(i) = y(i-1) + (k1 + 2*k2 + 2*k3 + k4)/6;
end
end
function plot_solution(y, x)
%Plots the solution to the ODE.
% Plot the solution.
figure
plot(x, y)
xlabel('x')
ylabel('y')
title('Solution to the ODE')
end
