给一个详细完整的在MATLAB中绘制阻尼单自由度系统强迫振动的幅频曲线和相频曲线详细完整代码

在MATLAB中，绘制阻尼单自由度系统的幅频响应和相频响应通常涉及到傅里叶变换和控制系统工具箱。下面是一个详细的步骤以及部分代码示例：

首先，假设我们有一个简谐振子的数学模型，其一般形式为：
\[ m\ddot{x} + c\dot{x} + kx = F_0 \sin(\omega t) \]

其中 \( m \) 是质量，\( c \) 是阻尼系数，\( k \) 是弹性系数，\( F_0 \) 是外部力的幅度，\( \omega \) 是激励频率。

1. 导入必要的库：

import ControlSystems.*

2. 定义模型参数：

m = 1; % 质量 (kg)
c = 0.5; % 阻尼系数 (Ns/m)
k = 2; % 弹性系数 (N/m)
F0 = 1; % 外部力幅度 (N)
w_vec = linspace(0, 10, 1000); % 频率范围从0到10 rad/s，步长自定

3. 创建传递函数模型：

sys = tf([1], [m c k]);

4. 对输入信号进行傅里叶变换得到频域表示：

excitation = input_signal * exp(1i*w_vec'*t); % 输入信号：正弦波
forced_response = lsim(sys, excitation);

`input_signal`需要自行定义，这里假设是一个简单的正弦波。

5. 计算幅频响应（Spectrogram）：

[Hmag, w] = bode(sys, w_vec);
amplitude_response = abs(Hmag);

6. 计算相频响应（Phase plot）：

[Hphase, ~] = bode(sys, w_vec);
phase_response = unwrap(angle(Hphase));

`unwrap`函数用于消除角度测量中的跳跃。

7. 绘制幅频曲线：

plot(w, amplitude_response, 'b', 'LineWidth', 2);
xlabel('Frequency (rad/s)');
ylabel('Amplitude Response');
title('Damped SDOF System Amplitude Frequency Response');

8. 绘制相频曲线：

plot(w, phase_response, 'r', 'LineWidth', 2);
grid on;
hold on;
xlabel('Frequency (rad/s)');
ylabel('Phase Response (degrees)');
legend('Amplitude', 'Phase');
title('Damped SDOF System Phase Frequency Response');

9. 结合幅频和相频图：

subplot(2,1,1);
plot(w, amplitude_response, 'b');
title('Amplitude Frequency Response');
subplot(2,1,2);
plot(w, phase_response, 'r');
title('Phase Frequency Response');