帮我写一个受迫振动的matlab代码

以下是一个简单的受迫振动 Matlab 代码：

% 受迫振动的 Matlab 代码
% 函数：forced_oscillator
% 参数：m - 相对质量，k - 劲度系数，
%       gamma - 粘性阻力系数，F0 - 外力幅值，omega - 外力频率，
%       tspan - 时间跨度，y0 - 初始条件
% 返回值：t - 时间向量，y - 状态向量矩阵

function [t, y] = forced_oscillator(m, k, gamma, F0, omega, tspan, y0)

% 定义受迫振动的微分方程
f = @(t, y) [y(2); -(gamma/m)*y(2)-(k/m)*y(1)+(F0/m)*sin(omega*t)];

% 使用 ODE45 求解微分方程
[t, y] = ode45(f, tspan, y0);

end

使用方法：

1. 将上述代码复制到一个 m 文件中，例如 forced_oscillator.m；
2. 在 Matlab 中调用该函数，例如：

% 参数
m = 1;
k = 1;
gamma = 0.2;
F0 = 10;
omega = 2*pi;
tspan = [0 10];
y0 = [0; 0];

% 求解微分方程
[t, y] = forced_oscillator(m, k, gamma, F0, omega, tspan, y0);

% 绘制振动图像
figure;
plot(t, y(:,1));
title('Forced oscillator');
xlabel('Time (s)');
ylabel('Position (m)');
grid on;

这将绘制出受迫振动的位移随时间的变化图像。注意，这是一个非线性系统，因此可能会出现混沌现象。要进行更全面的研究，请参阅相关文献。