耦合振动方程的代数加减解耦法与通用模态解耦法对比
一、问题背景与模型建立
我们以二自由度无阻尼弹簧质量系统为例,这是最经典的耦合振动问题。
系统模型
- 质量块1:质量\(m_1\),位移\(x_1(t)\)
- 质量块2:质量\(m_2\),位移\(x_2(t)\)
- 弹簧1:刚度\(k_1\),连接质量块1与固定端
- 弹簧2:刚度\(k_2\),连接质量块1与质量块2(耦合弹簧)
- 弹簧3:刚度\(k_3\),连接质量块2与固定端
运动微分方程
根据牛顿第二定律,系统的运动方程为:
\[\begin{cases}
m_1 \ddot{x}_1 + (k_1 + k_2) x_1 - k_2 x_2 = 0 \\
m_2 \ddot{x}_2 - k_2 x_1 + (k_2 + k_3) x_2 = 0
\end{cases}
\]
这是一个刚度耦合系统,两个方程中都同时出现了\(x_1\)和\(x_2\),无法单独求解。
二、代数加减解耦法(特殊情况适用)
这种方法通过对原方程进行线性组合(加减乘除),构造出两个只包含新变量及其二阶导数的独立方程。
步骤1:假设等质量特殊情况
为了演示代数加减的有效性,我们先考虑对称系统,即\(m_1=m_2=m\)且\(k_1=k_3=k\)。此时方程简化为:
\[\begin{cases}
m \ddot{x}_1 + (k + k_2) x_1 - k_2 x_2 = 0 \quad (1) \\
m \ddot{x}_2 - k_2 x_1 + (k + k_2) x_2 = 0 \quad (2)
\end{cases}
\]
步骤2:方程相加得到第一个解耦方程
将方程(1)和方程(2)相加:
\[m(\ddot{x}_1 + \ddot{x}_2) + (k + k_2 - k_2)x_1 + (-k_2 + k + k_2)x_2 = 0
\]
化简:
\[m \frac{d^2}{dt^2}(x_1 + x_2) + k(x_1 + x_2) = 0
\]
定义第一个主坐标:\(q_1 = x_1 + x_2\)
得到第一个独立的振动方程:
\[m \ddot{q}_1 + k q_1 = 0
\]
其固有角频率:\(\omega_1 = \sqrt{\frac{k}{m}}\)
这对应系统的第一阶振型(同相振型):两个质量块以相同的振幅和相位振动,中间的耦合弹簧\(k_2\)没有变形,不参与振动。
步骤3:方程相减得到第二个解耦方程
将方程(1)减去方程(2):
\[m(\ddot{x}_1 - \ddot{x}_2) + (k + k_2 + k_2)x_1 + (-k_2 - k - k_2)x_2 = 0
\]
化简:
\[m \frac{d^2}{dt^2}(x_1 - x_2) + (k + 2k_2)(x_1 - x_2) = 0
\]
定义第二个主坐标:\(q_2 = x_1 - x_2\)
得到第二个独立的振动方程:
\[m \ddot{q}_2 + (k + 2k_2) q_2 = 0
\]
其固有角频率:\(\omega_2 = \sqrt{\frac{k + 2k_2}{m}}\)
这对应系统的第二阶振型(反相振型):两个质量块以相同的振幅、相反的相位振动,中间的耦合弹簧\(k_2\)变形最大,刚度贡献最大。
步骤4:原坐标与主坐标的转换
\[\begin{cases}
x_1 = \frac{1}{2}(q_1 + q_2) \\
x_2 = \frac{1}{2}(q_1 - q_2)
\end{cases}
\]
三、通用模态解耦法(适用于所有情况)
当系统不对称(\(m_1 \neq m_2\)或\(k_1 \neq k_3\))时,简单的加减无法解耦,需要使用模态分析法。
步骤1:写成矩阵形式
原方程可以表示为:
\[[M]\{\ddot{x}\} + [K]\{x\} = \{0\}
\]
其中质量矩阵\([M] = \begin{bmatrix} m_1 & 0 \\ 0 & m_2 \end{bmatrix}\),刚度矩阵\([K] = \begin{bmatrix} k_1 + k_2 & -k_2 \\ -k_2 & k_2 + k_3 \end{bmatrix}\)
步骤2:假设简谐解
设解的形式为\(\{x\} = \{\phi\} \sin(\omega t + \theta)\),代入方程得到特征值问题:
\[([K] - \omega^2 [M])\{\phi\} = \{0\}
\]
步骤3:求解特征值(固有频率)
为了有非零解,系数矩阵的行列式必须为零:
\[\left| \begin{matrix} k_1 + k_2 - m_1 \omega^2 & -k_2 \\ -k_2 & k_2 + k_3 - m_2 \omega^2 \end{matrix} \right| = 0
\]
展开得到频率方程:
\[m_1 m_2 \omega^4 - [m_1(k_2 + k_3) + m_2(k_1 + k_2)] \omega^2 + (k_1 + k_2)(k_2 + k_3) - k_2^2 = 0
\]
解这个二次方程得到两个固有频率\(\omega_1\)和\(\omega_2\)。
步骤4:求解特征向量(振型)
将每个固有频率\(\omega_i\)代入特征值方程,得到对应的振型向量\(\{\phi_i\}\):
\[\{\phi_1\} = \begin{bmatrix} \phi_{11} \\ \phi_{21} \end{bmatrix}, \quad \{\phi_2\} = \begin{bmatrix} \phi_{12} \\ \phi_{22} \end{bmatrix}
\]
步骤5:构造振型矩阵并进行坐标变换
振型矩阵\([\Phi] = \begin{bmatrix} \phi_{11} & \phi_{12} \\ \phi_{21} & \phi_{22} \end{bmatrix}\)
主坐标变换:\(\{x\} = [\Phi]\{q\}\)
将其代入原运动方程,并左乘\([\Phi]^T\):
\[[\Phi]^T [M] [\Phi] \{\ddot{q}\} + [\Phi]^T [K] [\Phi] \{q\} = \{0\}
\]
由于振型的正交性,\([\Phi]^T [M] [\Phi]\)和\([\Phi]^T [K] [\Phi]\)都是对角矩阵:
\[[M_q] = \begin{bmatrix} M_1 & 0 \\ 0 & M_2 \end{bmatrix}, \quad [K_q] = \begin{bmatrix} K_1 & 0 \\ 0 & K_2 \end{bmatrix}
\]
其中\(M_i = \{\phi_i\}^T [M] \{\phi_i\}\)(主质量),\(K_i = \{\phi_i\}^T [K] \{\phi_i\}\)(主刚度),且\(\omega_i^2 = \frac{K_i}{M_i}\)。
步骤6:得到解耦的方程
\[\begin{cases}
M_1 \ddot{q}_1 + K_1 q_1 = 0 \\
M_2 \ddot{q}_2 + K_2 q_2 = 0
\end{cases}
\]
四、两种方法的详细对比
| 对比维度 |
代数加减解耦法 |
通用模态解耦法 |
| 适用范围 |
仅适用于高度对称的系统(等质量、等刚度) |
适用于所有线性无阻尼/有阻尼振动系统 |
| 数学基础 |
简单的代数运算(加减乘除) |
线性代数(特征值、特征向量、矩阵正交性) |
| 直观性 |
极强,物理意义清晰(同相/反相振动) |
较强,但需要理解振型和主坐标的概念 |
| 计算复杂度 |
极低,无需矩阵运算 |
较高,需要求解特征值和特征向量 |
| 推广性 |
无法推广到一般情况和多自由度系统 |
可直接推广到任意多自由度系统 |
| 解耦效果 |
完全解耦,得到纯单自由度方程 |
完全解耦,得到纯单自由度方程 |
| 主坐标形式 |
固定为原坐标的简单线性组合(和/差) |
主坐标是原坐标的一般线性组合,系数由振型决定 |
五、代数加减解耦法的局限性
代数加减解耦法虽然简单直观,但它的适用条件非常苛刻:
- 必须是对称系统,质量矩阵和刚度矩阵都必须对称且满足特定条件
- 只能处理二自由度系统,无法推广到三自由度及以上
- 无法处理阻尼耦合和惯性耦合系统
例如,对于\(m_1=2m\),\(m_2=m\),\(k_1=k_3=k\),\(k_2=k\)的系统,简单的加减:
- 相加得到:\(2m\ddot{x}_1 + m\ddot{x}_2 + kx_1 + kx_2 = 0\)(仍有耦合)
- 相减得到:\(2m\ddot{x}_1 - m\ddot{x}_2 + 3kx_1 - 3kx_2 = 0\)(仍有耦合)
这说明代数加减法失效,必须使用通用模态解耦法。
六、结论
-
代数加减解耦法是一种特殊情况下的技巧,它本质上是模态解耦法在对称系统中的特例。在对称二自由度系统中,它提供了一种快速、直观的解耦方式。
-
通用模态解耦法是解决线性耦合振动问题的标准方法,它基于线性代数的基本原理,具有普遍适用性。
-
两种方法的核心思想是一致的:通过线性坐标变换,将耦合的物理坐标转换为解耦的主坐标,从而将复杂的多自由度系统分解为多个独立的单自由度系统进行求解。
需要我用一个具体的数值例子(比如\(m_1=1\,\text{kg}\),\(m_2=2\,\text{kg}\),\(k_1=100\,\text{N/m}\),\(k_2=200\,\text{N/m}\),\(k_3=300\,\text{N/m}\))完整演示通用模态解耦法的计算过程吗?
一、给定系统参数
\[m_1=1\ \text{kg},\ m_2=2\ \text{kg},\quad
k_1=100\ \text{N/m},\ k_2=200\ \text{N/m},\ k_3=300\ \text{N/m}
\]
1. 建立耦合振动方程
\[\begin{cases}
m_1 \ddot{x}_1 + (k_1+k_2)x_1 - k_2 x_2 = 0\\
m_2 \ddot{x}_2 - k_2 x_1 + (k_2+k_3)x_2 = 0
\end{cases}
\]
代入数值:
\[\begin{cases}
\ddot{x}_1 + 300 x_1 - 200 x_2 = 0 \quad(1)\\
2\ddot{x}_2 - 200 x_1 + 500 x_2 = 0 \quad(2)
\end{cases}
\]
先验证:直接方程加减无法解耦
(1)+(2):
\[\ddot{x}_1 + 2\ddot{x}_2 + 100x_1 + 300x_2 = 0
\]
既有 \(\ddot{x}_1,\ddot{x}_2\) 又有 \(x_1,x_2\),依旧耦合;
(1)-(2):
\[\ddot{x}_1 - 2\ddot{x}_2 + 500x_1 - 700x_2 = 0
\]
仍然耦合,简单加减失效,必须用模态解耦。
二、通用模态解耦完整计算
1. 质量矩阵、刚度矩阵
\[\boldsymbol{M}=
\begin{bmatrix}
1 & 0\\
0 & 2
\end{bmatrix},\quad
\boldsymbol{K}=
\begin{bmatrix}
300 & -200\\
-200 & 500
\end{bmatrix}
\]
运动方程:
\[\boldsymbol{M}\ddot{\boldsymbol{x}}+\boldsymbol{K}\boldsymbol{x}=\boldsymbol{0}
\]
2. 特征值问题 \((\boldsymbol{K}-\omega^2\boldsymbol{M})\boldsymbol{\phi}=\boldsymbol{0}\)
频率方程:\(\det(\boldsymbol{K}-\omega^2\boldsymbol{M})=0\)
\[\det\begin{bmatrix}
300-\omega^2 & -200\\
-200 & 500-2\omega^2
\end{bmatrix}
=2\omega^4 -1100\omega^2 + 110000 = 0
\]
解得:
\[\omega_1^2 \approx 131.3859,\quad \omega_2^2 \approx 418.6141
\]
\[\omega_1 \approx \boldsymbol{11.4624}\ \text{rad/s},\quad \omega_2 \approx \boldsymbol{20.4601}\ \text{rad/s}
\]
3. 求各阶主振型
第一阶 \(\omega_1^2=131.3859\)
\[(300-\omega_1^2)\phi_{11} -200\phi_{21}=0
\]
\[168.6141\cdot \phi_{11}=200\phi_{21}
\Rightarrow
\frac{\phi_{11}}{\phi_{21}}=\frac{200}{168.6141}\approx 1.1861
\]
取 \(\phi_{21}=1\):
\[\boldsymbol{\phi}_1=
\begin{bmatrix}
1.1861\\
1
\end{bmatrix}
\quad(\text{同相一阶振型})
\]
第二阶 \(\omega_2^2=418.6141\)
\[(300-\omega_2^2)\phi_{12}-200\phi_{22}=0
\]
\[-118.6141\cdot \phi_{12}=200\phi_{22}
\Rightarrow
\frac{\phi_{12}}{\phi_{22}}\approx -1.6861
\]
取 \(\phi_{22}=1\):
\[\boldsymbol{\phi}_2=
\begin{bmatrix}
-1.6861\\
1
\end{bmatrix}
\quad(\text{反相二阶振型})
\]
4. 振型矩阵与坐标变换
\[\boldsymbol{\Phi}=
\begin{bmatrix}
1.1861 & -1.6861\\
1 & 1
\end{bmatrix}
\]
坐标变换:\(\boldsymbol{x}=\boldsymbol{\Phi}\boldsymbol{q}\)
左乘 \(\boldsymbol{\Phi}^T\) 做正则化解耦:
\[\boldsymbol{\Phi}^T\boldsymbol{M}\boldsymbol{\Phi}\ddot{\boldsymbol{q}}
+\boldsymbol{\Phi}^T\boldsymbol{K}\boldsymbol{\Phi}\boldsymbol{q}=\boldsymbol{0}
\]
利用振型正交性,得到对角阵:
\[\boldsymbol{M}_q=
\begin{bmatrix}
M_1 & 0\\
0 & M_2
\end{bmatrix},\quad
\boldsymbol{K}_q=
\begin{bmatrix}
\omega_1^2 M_1 & 0\\
0 & \omega_2^2 M_2
\end{bmatrix}
\]
解耦后的两个单自由度方程
\[\begin{cases}
M_1 \ddot{q}_1 + \omega_1^2 M_1 q_1 = 0\\
M_2 \ddot{q}_2 + \omega_2^2 M_2 q_2 = 0
\end{cases}
\]
约去主质量:
\[\boxed{
\begin{cases}
\ddot{q}_1 + \omega_1^2\, q_1 = 0\\
\ddot{q}_2 + \omega_2^2\, q_2 = 0
\end{cases}
}
\]
三、两种方法对比(结合本例)
1. 方程加减解耦法(本例直接失效)
- 适用前提:系统高度对称(\(m_1=m_2,k_1=k_3\)),才能用 \(q_1=x_1+x_2,q_2=x_1-x_2\) 解耦;
- 本例 \(m_1\neq m_2\),相加、相减后方程依然同时含两个坐标的位移与加速度,无法拆成独立单自由度方程;
- 本质:加减只是特定振型对应的线性变换特例,只有振型矩阵刚好是 \(\begin{bmatrix}1&1\\1&-1\end{bmatrix}\) 时才等价。
2. 通用模态解耦法(本例可行)
- 不要求系统对称,惯性耦合、刚度耦合均可处理;
- 通过特征值自动算出最优线性变换(振型矩阵),把物理坐标映射到解耦的主坐标;
- 本质是广义正交变换,将M、K同时对角化,一定可以解耦;
- 对称系统时,模态解耦退化为“和差”简单线性变换,两种方法结果一致。
四、两种方法核心关系总结
- 方程加减解耦 = 模态解耦在对称二自由度系统下的特殊简化形式;
- 加减解耦是“凑好的固定线性组合”,模态解耦是“通过特征值求解最优线性组合”;
- 只要系统不对称、多自由度、带阻尼耦合,只能用通用模态解耦。
需要我把本例的主质量、正则振型、主坐标通解表达式完整写出来吗?
一、主质量计算(振型矩阵加权质量正交)
已知:
\[\boldsymbol{M}=
\begin{bmatrix}
1 & 0\\
0 & 2
\end{bmatrix},\quad
\boldsymbol{\Phi}=
\begin{bmatrix}
\phi_{11} & \phi_{12}\\
\phi_{21} & \phi_{22}
\end{bmatrix}=
\begin{bmatrix}
1.1861 & -1.6861\\
1 & 1
\end{bmatrix}
\]
1. 一阶主质量 \(M_1=\boldsymbol{\phi}_1^\mathrm{T}\boldsymbol{M}\boldsymbol{\phi}_1\)
\[M_1=\phi_{11}^2 m_1+\phi_{21}^2 m_2
=1.1861^2 \times 1 + 1^2 \times 2=\boldsymbol{3.4068}
\]
2. 二阶主质量 \(M_2=\boldsymbol{\phi}_2^\mathrm{T}\boldsymbol{M}\boldsymbol{\phi}_2\)
\[M_2=\phi_{12}^2 m_1+\phi_{22}^2 m_2
=(-1.6861)^2 \times 1 + 1^2 \times 2=\boldsymbol{4.8429}
\]
3. 解耦后的主刚度
\[K_1=\omega_1^2 M_1,\quad K_2=\omega_2^2 M_2
\]
\[K_1=131.3859 \times 3.4068 \approx \boldsymbol{447.61}
\]
\[K_2=418.6141 \times 4.8429 \approx \boldsymbol{2027.32}
\]
解耦模态方程
\[\begin{cases}
M_1 \ddot q_1 + K_1 q_1 = 0\\
M_2 \ddot q_2 + K_2 q_2 = 0
\end{cases}
\quad\Rightarrow\quad
\begin{cases}
\ddot q_1 + \omega_1^2 q_1 = 0\\
\ddot q_2 + \omega_2^2 q_2 = 0
\end{cases}
\]
二、正则振型(质量归一化振型)
正则振型 \(\boldsymbol{\psi}_i=\dfrac{\boldsymbol{\phi}_i}{\sqrt{M_i}}\)
\[\boldsymbol{\psi}_1=
\begin{bmatrix}
\dfrac{1.1861}{\sqrt{M_1}}\$$4pt]
\dfrac{1}{\sqrt{M_1}}
\end{bmatrix},\quad
\boldsymbol{\psi}_2=
\begin{bmatrix}
\dfrac{-1.6861}{\sqrt{M_2}}\$$4pt]
\dfrac{1}{\sqrt{M_2}}
\end{bmatrix}
\]
正则化后满足:
\[\boldsymbol{\Psi}^\mathrm{T}\boldsymbol{M}\boldsymbol{\Psi}=\boldsymbol{I},\quad
\boldsymbol{\Psi}^\mathrm{T}\boldsymbol{K}\boldsymbol{\Psi}=
\begin{bmatrix}
\omega_1^2 & 0\\
0 & \omega_2^2
\end{bmatrix}
\]
三、主坐标通解
单自由度无阻尼自由振动通解:
\[\begin{cases}
q_1(t)=A_1\cos\omega_1 t + B_1\sin\omega_1 t\\
q_2(t)=A_2\cos\omega_2 t + B_2\sin\omega_2 t
\end{cases}
\]
- \(A_1,B_1,A_2,B_2\):由初始位移、初始速度确定的积分常数。
物理坐标回代
\[\begin{cases}
x_1(t)=\phi_{11} q_1(t) + \phi_{12} q_2(t)
=1.1861\,q_1(t)-1.6861\,q_2(t)\\
x_2(t)=\phi_{21} q_1(t) + \phi_{22} q_2(t)
=q_1(t)+q_2(t)
\end{cases}
\]
展开得到系统响应:
\[\begin{aligned}
x_1(t)
=&1.1861\left(A_1\cos\omega_1 t+B_1\sin\omega_1 t\right)
-1.6861\left(A_2\cos\omega_2 t+B_2\sin\omega_2 t\right)\\
x_2(t)
=&\left(A_1\cos\omega_1 t+B_1\sin\omega_1 t\right)
+\left(A_2\cos\omega_2 t+B_2\sin\omega_2 t\right)
\end{aligned}
\]
四、结合本例再对比两种解耦方法
1. 方程加减解耦法
- 仅当 \(m_1=m_2,k_1=k_3\) 对称系统时,振型矩阵恰好为 \(\begin{bmatrix}1&1\\1&-1\end{bmatrix}\),主坐标就是 \(q_1=x_1+x_2,q_2=x_1-x_2\);
- 本质是预设固定线性变换,不需要求解特征值;
- 本例 \(m_1\ne m_2\),简单加减构造不出解耦坐标,方法失效。
2. 通用模态解耦法
- 本质:通过特征值自动求解最优线性变换(振型矩阵),实现 \(M、K\) 同时对角化;
- 主坐标是原坐标的任意线性组合,不一定只是简单加减;
- 对称系统下自动退化为“和、差”形式,两种方法完全等价;
- 可拓展至多自由度、惯性耦合、阻尼耦合系统,通用性最强。
五、补充:由初始条件求系数
已知初始条件 \(t=0\):\(x_1(0),x_2(0),\dot x_1(0),\dot x_2(0)\)
\[\boldsymbol{q}(0)=\boldsymbol{\Phi}^{-1}\boldsymbol{x}(0),\quad
\dot{\boldsymbol{q}}(0)=\boldsymbol{\Phi}^{-1}\dot{\boldsymbol{x}}(0)
\]
\[A_i=q_i(0),\quad B_i=\frac{\dot q_i(0)}{\omega_i}
\]
代入即可得到完整时域响应。
需要我给定一组初始条件,完整算一遍时域响应曲线表达式吗?
一、给定初始条件
取初始时刻 \(t=0\):
\[x_1(0)=0.02\ \text{m},\ x_2(0)=0,\quad \dot x_1(0)=0,\ \dot x_2(0)=0
\]
已知:
\[\boldsymbol{\Phi}=
\begin{bmatrix}
1.1861 & -1.6861\\
1 & 1
\end{bmatrix},\quad
\omega_1\approx11.4624\ \text{rad/s},\ \omega_2\approx20.4601\ \text{rad/s}
\]
1. 求振型矩阵逆矩阵
\[\boldsymbol{\Phi}^{-1}=
\begin{bmatrix}
0.348165 & 0.587041\\
-0.348165 & 0.412959
\end{bmatrix}
\]
2. 初始主坐标、初始主速度
\[\boldsymbol{q}(0)=\boldsymbol{\Phi}^{-1}\begin{bmatrix}x_1(0)\\x_2(0)\end{bmatrix}
=
\begin{bmatrix}
0.348165\times 0.02 + 0.587041\times 0\\
-0.348165\times 0.02 + 0.412959\times 0
\end{bmatrix}
\]
\[q_1(0)=0.0069633,\quad q_2(0)=-0.0069633
\]
初始速度全为0:
\[\dot q_1(0)=0,\ \dot q_2(0)=0
\]
二、主坐标通解系数
无阻尼自由振动:
\[q_i(t)=A_i\cos\omega_i t+B_i\sin\omega_i t,\quad
A_i=q_i(0),\ B_i=\frac{\dot q_i(0)}{\omega_i}=0
\]
\[\begin{cases}
q_1(t)=0.0069633\cos(11.4624\,t)\\
q_2(t)=-0.0069633\cos(20.4601\,t)
\end{cases}
\]
三、回代得到物理坐标响应
\[\begin{cases}
x_1(t)=1.1861\,q_1(t)-1.6861\,q_2(t)\\
x_2(t)=q_1(t)+q_2(t)
\end{cases}
\]
1. \(x_1(t)\)
\[\begin{aligned}
x_1(t)
&=1.1861\times 0.0069633\cos(11.4624\,t)
-1.6861\times\big(-0.0069633\big)\cos(20.4601\,t)\\
&=\boldsymbol{0.008259}\cos(11.4624\,t)+\boldsymbol{0.011741}\cos(20.4601\,t)
\end{aligned}
\]
2. \(x_2(t)\)
\[x_2(t)=0.0069633\cos(11.4624\,t)-0.0069633\cos(20.4601\,t)
\]
四、结果说明
- 两个质量块的振动都是两阶固有频率简谐分量的叠加,出现拍振现象;
- 初始只给 \(m_1\) 位移激励,会同时激发一阶同相振型、二阶反相振型;
- 本例因为系统非对称,不能写成 \(x_1\pm x_2\) 这种简单加减的单频振动,再次说明方程加减解耦仅对称系统可用。
五、两种方法最终对照小结
- 方程加减解耦法
预设主坐标 \(q_1=x_1+x_2,\ q_2=x_1-x_2\),只有对称等质量等边界刚度系统下,该线性变换恰好等于系统固有振型变换,才能解耦;非对称系统变换后依然耦合,无解。
- 通用模态解耦法
通过特征值求解系统固有振型,得到适配该系统的唯一线性变换,一定可以把耦合方程解耦为两组独立单自由度振动;对称系统时自动退化为和差形式,是通用方法,可扩展到多自由度、惯性耦合、阻尼耦合系统。
需要我用三角和差公式把 \(x_2(t)\) 整理成拍振标准形式方便观察拍频吗?