拟合衰减振动模型，估算阻尼比和阻尼系数

拟合衰减振动模型，估算阻尼比和阻尼系数

引用

CSDN

1.

https://blog.csdn.net/flyfish1986/article/details/140323857

在工程和物理学中，衰减振动模型是分析阻尼系统振动特性的重要工具。本文将介绍如何通过Python代码拟合衰减振动数据，从而估算系统的阻尼比和阻尼系数。

衰减振动模型

在自由振动系统中，阻尼振动可以用以下公式描述：

$$
x(t) = x_0 e^{-\zeta \omega_n t} \cos(\omega_d t + \phi)
$$

其中：

• $x(t)$ 是时间 $t$ 时的位移（Displacement at time $t$）
• $x_0$ 是初始位移（Initial displacement）
• $\zeta$ 是阻尼比（Damping ratio）
• $\omega_n$ 是无阻尼固有频率（Undamped natural frequency）
• $t$ 是时间（Time）
• $\omega_d$ 是阻尼振动频率（Damped natural frequency），其公式为：
  $$
  \omega_d = \omega_n \sqrt{1 - \zeta^2}
  $$
• $\phi$ 是相位角（Phase angle）

公式说明

1. 初始位移 $x_0$ (Initial displacement):

• 系统开始自由振动时的位移。

2. 阻尼比 $\zeta$ (Damping ratio):

• 衡量系统阻尼程度的无量纲参数。
• $\zeta$ 越大，阻尼越强。

3. 无阻尼固有频率 $\omega_n$ (Undamped natural frequency):

• 系统在没有阻尼情况下的固有振动频率。
• 单位：弧度每秒（radians per second, rad/s）。

4. 阻尼振动频率 $\omega_d$ (Damped natural frequency):

• 系统在有阻尼情况下的实际振动频率。
• 计算公式：$\omega_d = \omega_n \sqrt{1 - \zeta^2}$
• 当 $\zeta < 1$ 时，系统会振荡，$\omega_d$ 表示振荡频率。

5. 相位角 $\phi$ (Phase angle):

• 描述系统初始位移和初始速度之间的关系。

衰减振动的描述

在阻尼系统中，振动会逐渐衰减，幅度随着时间指数下降。振动系统的运动可以分解为以下几个部分：

1. 指数衰减部分 $e^{-\zeta \omega_n t}$:

• 描述振幅随时间的衰减。

2. 余弦振荡部分 $\cos(\omega_d t + \phi)$:

• 描述系统的振荡行为，频率为 $\omega_d$，初始相位为 $\phi$。

代码实现

import numpy as np
import scipy.optimize as opt
import matplotlib.pyplot as plt

# 定义衰减振动函数
def damped_vibration(t, x0, zeta, omega_n, phi):
    omega_d = omega_n * np.sqrt(1 - zeta**2)
    return x0 * np.exp(-zeta * omega_n * t) * np.cos(omega_d * t + phi)

# 合成实验数据
t_data = np.linspace(0, 10, 1000)  # 时间数据
x0_true = 1.0  # 初始位移
zeta_true = 0.1  # 真正的阻尼比
omega_n_true = 2 * np.pi  # 真正的无阻尼固有频率
phi_true = 0  # 真正的相位角
x_data = damped_vibration(t_data, x0_true, zeta_true, omega_n_true, phi_true)  # 生成理想数据

# 加入噪声
x_data_noisy = x_data + 0.05 * np.random.normal(size=t_data.shape)  # 加入随机噪声

# 拟合数据
initial_guess = [x0_true, zeta_true, omega_n_true, phi_true]  # 初始猜测
params, params_covariance = opt.curve_fit(damped_vibration, t_data, x_data_noisy, p0=initial_guess)  # 拟合模型

# 提取拟合参数
x0_est, zeta_est, omega_n_est, phi_est = params

# 假设质量为1 kg
m = 1.0
gamma_est = 2 * zeta_est * omega_n_est * m

# 输出结果
print(f"Estimated Initial Displacement (x0): {x0_est}")
print(f"Estimated Damping Ratio (zeta): {zeta_est}")
print(f"Estimated Natural Frequency (omega_n): {omega_n_est}")
print(f"Estimated Damping Coefficient (gamma): {gamma_est} N·s/m")

# 绘图比较
plt.figure(figsize=(10, 6))
plt.plot(t_data, x_data_noisy, label='Noisy Data', alpha=0.5)
plt.plot(t_data, damped_vibration(t_data, *params), label='Fitted Curve', color='red')
plt.xlabel('Time (s)')
plt.ylabel('Displacement (m)')
plt.legend()
plt.show()

结果解读

• t_data: 时间数据，生成一个从0到10秒的1000个时间点的数组。
• x0_true: 初始位移（True initial displacement），设定为1.0米。
• zeta_true: 真正的阻尼比（True damping ratio），设定为0.1。
• omega_n_true: 真正的无阻尼固有频率（True undamped natural frequency），设定为 $2\pi$（即每秒一个完整的振动周期）。
• phi_true: 真正的相位角（True phase angle），设定为0。

拟合结果如下：

Estimated Initial Displacement (x0): 0.9969144391835894
Estimated Damping Ratio (zeta): 0.09975574475118448
Estimated Natural Frequency (omega_n): 6.2712971489170615
Estimated Damping Coefficient (gamma): 1.2511958352924026 N·s/m

每个结果的单位及其读法：

1. Estimated Initial Displacement (x0):

• 单位: 米 (meters, m)
• 读法: “Estimated Initial Displacement is value meters”

2. Estimated Damping Ratio (zeta):

• 单位: 无单位（damping ratio 是一个无量纲参数）
• 读法: “Estimated Damping Ratio is value ”

3. Estimated Natural Frequency (omega_n):

• 单位: 弧度每秒 (radians per second, rad/s)
• 读法: “Estimated Natural Frequency is value radians per second”

4. Estimated Damping Coefficient (gamma):

• 单位: 牛顿·秒每米 (Newton-seconds per meter, N·s/m)
• 读法: “Estimated Damping Coefficient is value Newton-seconds per meter”