MateLab + Latex
clear; % this examples demostrates how a numerical method can be implemented for % solving a general system of ODEs using vectors and vector-valued functions % Specifically, in this example we solve a system of ODEs of size 3: % y1' = y2+b % y2' = sin(1-exp(y3/a)) % y3' = a/(1+t) % with initial conditions % y1(0) = 0, y2(0) = 1, y3(0) = 0 % where a and b are some parameters % define an anonymous function for the right-hand side % (note that functions are separated by ';'. As the result, function f % returns a column-vector) f = @(t, y, a, b) [ y(2)+b; sin(1-exp(y(3)/a)); a/(1+t); ]; % Alternatively, function f could have been defined in a separate file as: % function out = f(t,y,a,b) % % out = zeros(3,1); % out(1) = y(2)+b; % out(2) = sin(1-exp(y(3)/a)); % out(3) = a/(1+t); % % end % note that f is a function of variables t, y and parameters a and b, % however solve_system_euler requires a function of just two variables t % and y. A new function of just two variables t and y with pre-defined % values for parameters a and b can be created as my_a = 5; my_b = 1.5; my_f = @(t, y) f(t, y, my_a, my_b); % exact solution exact_solution = @(t) [ sin(t)+my_b*t; cos(t); my_a*log(1+t);]; % initial conditions t_start = 0; y_start = [ 0; 1; 0; ]; % final time t_final = 1; % time step dt = 0.05; [t, y] = solve_system_euler(my_f, y_start, t_start, t_final, dt); y_exact = exact_solution(t); plot(t, y_exact(1,:), '-b'); hold on plot(t, y_exact(2,:), '-k'); plot(t, y_exact(3,:), '-r'); plot(t, y(1,:), 'ob'); plot(t, y(2,:), 'ok'); plot(t, y(3,:), 'or'); hold off % you can use latex in legend lgd = legend('$y_1^{\textnormal{exact}}$', '$y_2^{\textnormal{exact}}$', '$y_3^{\textnormal{exact}}$', ... '$y_1^{\textnormal{num}}$', '$y_2^{\textnormal{num}}$', '$y_3^{\textnormal{num}}$'); lgd.Interpreter = 'latex'; lgd.Location = 'northwest'; lgd.FontSize = 16; xlabel('Time'); ylabel('y_1, y_2, y_3'); set(gca,'FontSize',14);