Matlab mass pendulum simulation

I have tried to simulate a mass pendulum on my own to train my matlab skills. I have coded a function that shows no errors, but however it does not do what it is supposed to do.
I wanted to iterate it for every time step. My idea was to calculate the current forces and the resulting acceleration. The velocity and the position s would be the sum of all previous accelerations. But somehow its ignoring the spring and the friction.
function pendel (m,te,k,c,s0)
t = linspace(1,te,100*te);
g = 9.81;
s = zeros(1,length(t));
v = zeros(1,length(t));
a = zeros(1,length(t));
for i = 1:length(t)
f1 = -m*g;
f2 = -k*s(i);
f3 = -v(i)*c;
a(i) = (f1+f2+f3)/m;
v(i) = sum(a(1:i));
s(i) = sum(v(1:i))+s0;
end
subplot(3,1,1);
plot (t,s)
title('Ortsdiagramm')
subplot(3,1,2);
plot(t,v)
title('Geschwindigkeitsdiagramm')
subplot(3,1,3);
plot(t,a)
title('Beschleunigungsdiagramm')
end

回答(1 个)

It appears to ignore the spring and friction forces because s(i) and v(i) are always 0. That is, s(i) and v(i) are calculated after they are used, each time through the loop, so the s(i) and v(i) used in the force calculations are always the 0 used to initialize s and v. One way to fix this is to use s(i-1) and v(i-1) to calculate a(i), v(i) and s(i).
Also note that you need to include the time-step in calculating velocity from acceleration or position from velocity.
(I don't claim that the following is physically accurate; I used your equations as they were, and just changed time/index-related stuff as described above, but it seems more physically plausible.)
pendel(1,10,1,1,0);
function pendel (m,te,k,c,s0)
t = linspace(1,te,100*te);
g = 9.81;
s = zeros(1,length(t));
v = zeros(1,length(t));
a = zeros(1,length(t));
dt = diff(t);
s(1) = s0;
for i = 2:length(t)
f1 = -m*g;
f2 = -k*s(i-1);
f3 = -v(i-1)*c;
a(i) = (f1+f2+f3)/m;
v(i) = v(i-1)+a(i)*dt(i-1);
s(i) = s(i-1)+v(i)*dt(i-1);
end
subplot(3,1,1);
plot (t,s)
title('Ortsdiagramm')
subplot(3,1,2);
plot(t,v)
title('Geschwindigkeitsdiagramm')
subplot(3,1,3);
plot(t,a)
title('Beschleunigungsdiagramm')
end

类别

帮助中心File Exchange 中查找有关 Programming 的更多信息

产品

版本

R2021a

回答:

2022-1-8

Community Treasure Hunt

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

Start Hunting!

Translated by