Error in Cumtrapz?
26 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
KLETECH MOTORSPORTS
el 24 de Nov. de 2020
Comentada: John
el 19 de Oct. de 2022
Hi, i'm trying to integrate acceleration of a simple pendulum twice to get the displacement, and i've plotted it against time.
i.e, ![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/425598/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/425598/image.png)
On the same graph, i've plotted the displacement versus time
i.e, ![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/425603/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/425603/image.png)
with certain initial conditions as follows:
%this part integrates the acceleration(y1 = d^2(x)/dt^2) twice
% to get both velocity(y2) and displacement(y3)
t=0:0.01:1
w=50;
A=pi/2;
y1=-A*(w.^2)*sin(w*t); %here y1 is diff(x,t,2) , i.e, accln
y2=cumtrapz(t',y1'); %here y2 is diff(x,t) , i.e, velocity
y3=cumtrapz(t',y2); %here y3 is the double integrated acceleration
y3a= A*sin(w*t); %direct displacement expression
plot(t,y3, t,y3a)
legend('y3', 'displ')
There is an error in the blue curve obtained by double integration, as both curves should be similar.
I think it might have to do with an error in trapezoidal integration.
Any idea what the error is and how i can correct it?
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/425608/image.png)
1 comentario
John
el 19 de Oct. de 2022
Should you be taking into consideration the acceleration due to gravity ?
Respuesta aceptada
Walter Roberson
el 24 de Nov. de 2020
Remember after one integration you get a formula with a constant of integration, and when you integrate that, you get the double integral plus the integral of the constant of integration, plus a second constant of integration. And the integral of the first constant of integration is the constant times the variable of integration.
Your plot shows a linear trend. That is the constant of integration multiplied by time.
Your blue plot is therefore one of the mathematically correct solutions.
If it is not the solution you want then you need to adjust the constant of integration.
Note that mathematically, trapz(x, y) with equidistant x, is equivalent to (sum(y) - (y(1)+y(end))/2) * delta x
cumtrapz with constant x difference effectively proceeds step by step, taking the proceeding result and adding half of the previous y value and half of the current y value. Those halfs have an influence.
15 comentarios
Walter Roberson
el 25 de Nov. de 2020
Suppose that you have
A = sin(t)
then you would propose to calculate
V = int(A,t) -> -cos(t)
P = int(V,t) -> -sin(t)
So amplitude 1, frequency 1 cycle per 2*pi seconds, you would propose that the position would be strictly -sin(t)
I would then ask you: If the initial position was (say) 1 light year away from the origin, then how exactly is the object going to travel 1 light year in the infinitesimal time from 0 to dt ?
Your proposed formula does not take into account initial position or initial velocity.
Más respuestas (0)
Ver también
Categorías
Más información sobre Numerical Integration and Differentiation en Help Center y File Exchange.
Productos
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!