Fourier Analysis on a Damped Harmonic Oscillator

18 次查看（过去 30 天）
Nathan Falkner
Nathan Falkner 2022-4-11
回答： Mathieu NOE 2022-4-14
I need to test if my oscillator model is damped or not, and so I need to derive frequency in order to create a graph that provides me with the correct information. This is my code, I'm going to need to send all of it:
k = 1;
t = 0.01;
ind = 1;
tmax = 100;
v = 5;
gamma = 0.005;
h = 0.01;
t_array(1) = 0.01;
v_array(1) = 5;
x_pos_array(1) = 0;
while (t <= tmax)
a = -(2*gamma*v_array(ind))-(k*x_pos_array(ind));
x_pos_array(ind+1) =x_pos_array(ind) + v_array(ind)*h;
v_array(ind+1) = v_array(ind) + a*h;
t_array(ind+1) = t_array(ind) + h;
t=t_array(ind+1);
ind = ind + 1;
end
F = fft(x_pos_array);
AF = abs(F);
xF = x_pos_array;
for i = 1:10000
D(i) = AF(i)*t_array(i);
i = i + 1;
end
plot (D,x_pos_array)
Only x_pos_array affects the fourier transform, however the majority of it is so intertwined I had to send the whole thing.
Now, when comparing my results to known results (https://en.wikipedia.org/wiki/Harmonic_oscillator, https://au.mathworks.com/help/symbolic/physics-damped-harmonic-oscillator.html), the graph that this code spits out is nowhere similar to the graph that should be there.
Have I written the code relating to the fourier transformation and associated graph correctly? Or do I need to do something else?

请先登录，再回答此问题。

回答（1 个）

Mathieu NOE
Mathieu NOE 2022-4-14
hello
I improved a bit the firts section of your code but I don't get what you are trying to do wit the second portion of the cde
F = fft(x_pos_array);
AF = abs(F);
xF = x_pos_array;
for i = 1:10000
D(i) = AF(i)*t_array(i);
i = i + 1;
end
plot (D,x_pos_array)
improved code :
clc
clearvars
dzeta = 0.05; % damping ratio
% physical constants
m = 10; % mass
k = 1; % stiffness
c = 2*dzeta*sqrt(k*m); % viscous damping
ind = 1;
tmax = 100;
v = 5;
h = 0.01;
t_array(1) = 0.01;
v_array(1) = 5;
x_pos_array(1) = 0;
t = t_array(1);
while (t <= tmax)
a = -1/m*((c*v_array(ind))+(k*x_pos_array(ind))); % acceleration = forces divided by mass
v_array(ind+1) = v_array(ind) + a*h; % update velocity first
x_pos_array(ind+1) = x_pos_array(ind) + v_array(ind+1)*h; % and after position
t_array(ind+1) = t_array(ind) + h;
t=t_array(ind+1);
ind = ind + 1;
end
plot(t_array,x_pos_array);

类别

在 Help Center 和 File Exchange 中查找有关 Fast Fourier Transforms 的更多信息

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by