Fourier Analysis on a Damped Harmonic Oscillator

6 次查看(过去 30 天)
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 CenterFile Exchange 中查找有关 Fourier Analysis and Filtering 的更多信息

产品

Community Treasure Hunt

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

Start Hunting!

Translated by