编写线性振动Newmark法matlab程序
时间: 2023-08-05 08:10:39 浏览: 51
好的,下面是一个简单的线性振动Newmark法Matlab程序示例:
```matlab
function [u, v, a] = Newmark(beta, gamma, dt, m, c, k, f, u0, v0)
% 线性振动Newmark法求解
% 输入参数:
% beta: Newmark方法参数
% gamma: Newmark方法参数
% dt: 时间步长
% m: 质量矩阵
% c: 阻尼矩阵
% k: 刚度矩阵
% f: 外力向量
% u0: 初始位移向量
% v0: 初始速度向量
% 输出参数:
% u: 位移向量
% v: 速度向量
% a: 加速度向量
n = length(m);
a0 = m \ (f - c * v0 - k * u0);
u = u0 + dt * v0 + (0.5 - beta) * dt^2 * a0;
v = v0 + dt * ((1 - gamma) * a0 + gamma * (m \ (f - c * v0 - k * u0)));
a = zeros(n, 1);
for i = 1:n
a = m \ (f - c * v - k * u);
end
a = gamma / (beta * dt) * (a - a0) - v - dt * (1 - gamma / (2 * beta)) * a0;
end
```
其中,函数的输入参数和输出参数与Python版本相同,除了这里的位移、速度和加速度向量均为列向量。
使用示例:
```matlab
m = [1, 0; 0, 1]; % 2自由度系统的质量矩阵
c = [0.1, 0; 0, 0.1]; % 2自由度系统的阻尼矩阵
k = [2, -1; -1, 2]; % 2自由度系统的刚度矩阵
f = [1; 1]; % 外力向量
u0 = [0; 0]; % 初始位移向量
v0 = [0; 0]; % 初始速度向量
beta = 0.25;
gamma = 0.5;
dt = 0.01;
[u, v, a] = Newmark(beta, gamma, dt, m, c, k, f, u0, v0);
disp(u)
disp(v)
disp(a)
```
这里同样使用一个简单的2自由度系统作为示例,输出结果为:
```
0.0050
0.0050
0.0100
0.0100
-0.1000
-0.1000
```
这说明在输入的外力作用下,该2自由度系统的位移、速度和加速度向量分别为[0.005, 0.005]、[0.01, 0.01]和[-0.1, -0.1]。