matlab多自由度体系振型叠加法采用正则坐标求解无阻尼强迫振动 代码实现

以下是一个简单的 Matlab 代码实现多自由度体系振型叠加法求解无阻尼强迫振动的过程。

首先定义系统的参数，包括质量矩阵、刚度矩阵和外力向量：

% 定义系统参数
M = [1 0; 0 2];  % 质量矩阵
K = [3 -1; -1 2];  % 刚度矩阵
F = [2*sin(2*t); 0];  % 外力向量

接着，由于采用正则坐标求解，需要计算正则坐标变换矩阵：

% 计算正则坐标变换矩阵
[V, D] = eig(K, M);  % 求解广义特征值问题
P = V' * M;  % 正则坐标变换矩阵

然后，根据正则坐标求解每个振型的振幅和相位，将其叠加得到总振动响应：

% 叠加各振型的振动响应
X = zeros(length(F), 1);  % 初始化总振动响应
for i = 1:length(F)
    w = sqrt(D(i, i));  % 求解振型频率
    C = P(i, :) * F / (w^2 - D(i, i));  % 求解振幅
    phi = atan2(w * (P(i, 2) * F(1) - P(i, 1) * F(2)), ...
        w^2 - D(i, i) - w * (P(i, 1)^2 * F(1) + P(i, 2)^2 * F(2)));  % 求解相位
    X = X + C * V(:, i) * sin(w * t + phi);  % 叠加振动响应
end

最后，将得到的总振动响应作图展示：

% 绘制振动响应图像
plot(t, X);
xlabel('Time (s)');
ylabel('Displacement (m)');
legend('m1', 'm2');

完整的代码如下：

% 定义时间范围和时间步长
t = 0:0.01:10;

% 定义系统参数
M = [1 0; 0 2];  % 质量矩阵
K = [3 -1; -1 2];  % 刚度矩阵
F = [2*sin(2*t); 0];  % 外力向量

% 计算正则坐标变换矩阵
[V, D] = eig(K, M);  % 求解广义特征值问题
P = V' * M;  % 正则坐标变换矩阵

% 叠加各振型的振动响应
X = zeros(length(F), 1);  % 初始化总振动响应
for i = 1:length(F)
    w = sqrt(D(i, i));  % 求解振型频率
    C = P(i, :) * F / (w^2 - D(i, i));  % 求解振幅
    phi = atan2(w * (P(i, 2) * F(1) - P(i, 1) * F(2)), ...
        w^2 - D(i, i) - w * (P(i, 1)^2 * F(1) + P(i, 2)^2 * F(2)));  % 求解相位
    X = X + C * V(:, i) * sin(w * t + phi);  % 叠加振动响应
end

% 绘制振动响应图像
plot(t, X);
xlabel('Time (s)');
ylabel('Displacement (m)');
legend('m1', 'm2');