Midpoint integration with for loop

5 views (last 30 days)
Hugo Hernández Hernández
Answered: arushi on 23 Aug 2024
Hi,
Computing a numerical integration with the Midpoint Method I'm struggling with the output of my function, while with Euler I got the expected result, with this method I'm getting a different graphical solution far away than the analytical solution (and the euler's one). What could be wrongly specified?
Thanks for your help.
% Midpoint method
% General information
L = 0.61; % [m]
G = 9.81; % gamma [m/s^2]
% Sampling interval
h = 0.1;
%h = 0.05;
%h = 0.01;
ti = 0; % Initial time
tf = 4; % Final time
t = ti:h:tf; % Time vector
function [md] = midpoint(t,L,G,h)
% Initial Value problem
md = zeros(2, length(t));
md(1,1) = pi/40;
for i = 1:length(t)-1
pmd = md(2,i) + (h/2) + (-G/L)*sin(md(1,i)) + (h/2) * [md(2,i); (-G/L)*sin(md(1,i))];
md(:,i+1) = md(:,i) + h * pmd;
end
end
  1 Comment
Hugo Hernández Hernández
I have updated my code, got a better result but still the same problem, this comes with "pmd" where as I include an array inside another array there comes an error, I know that this is not possible but, I don't find another solution that generates this midpoint condition, and also applies to Runge-Kutta for 4th order.

Sign in to comment.

Answers (1)

arushi
arushi on 23 Aug 2024
Hi Hugo,
It seems there might be a mistake in the implementation of the Midpoint Method in your MATLAB code. The Midpoint Method for solving ordinary differential equations (ODEs) is a type of Runge-Kutta method, and it requires evaluating the function at the midpoint of the interval before taking a step.
  1. Calculate the midpoint value using the current slope.
  2. Evaluate the function at the midpoint.
  3. Use the slope at the midpoint to take a full step.
In your code, the line calculating pmd seems incorrect. You’re adding (h/2) directly to the velocities and accelerations, which doesn’t correspond to any part of the Midpoint Method. Instead, use the midpoint value to calculate the slope for the full step.
Hope this helps.

Categories

Find more on Loops and Conditional Statements in Help Center and File Exchange

Products

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!