第6章：物理建模基础

物理建模合成是一种通过数学模拟真实乐器的物理行为来产生声音的技术。与基于频谱的合成方法（如FM或加法合成）不同，物理建模从声音产生的根源——振动系统的物理特性出发，通过求解描述这些系统的微分方程来生成声音。本章将介绍物理建模的基础理论，从经典的Karplus-Strong算法开始，逐步深入到数字波导理论，最后探讨弦乐器和管乐器的建模技术。

6.1 Karplus-Strong算法

6.1.1 基本原理与数学模型

Karplus-Strong算法（1983年）是物理建模合成的里程碑，它用极其简单的结构实现了逼真的拨弦声音。该算法的核心是一个带有低通滤波器的延迟反馈回路：

     噪声激励
         |
         v
    [延迟线 z^(-L)]<----+
         |              |
         v              |
    [低通滤波器]        |

    [低通滤波器]        |
         |              |

         +-------->-----+
         |
         v
       输出

数学表达式为： $$y(n) = x(n) + g \cdot H(z) \cdot y(n-L)$$ 其中：

• $x(n)$ 是激励信号（初始时为白噪声脉冲）
• $L$ 是延迟长度，决定基频：$f_0 = f_s / L$
• $g$ 是反馈增益（通常 < 1）
• $H(z)$ 是低通滤波器

6.1.2 延迟线与反馈

延迟线长度$L$直接决定了合成音的基频。对于采样率$f_s = 44100$ Hz，要产生A4（440 Hz）音符： $$L = \frac{f_s}{f_0} = \frac{44100}{440} = 100.23$$ 由于$L$必须是整数，实际频率会有偏差。解决方案包括：

1. 使用分数延迟插值
2. 调整采样率
3. 使用全通滤波器进行相位补偿

6.1.3 激励信号设计

激励信号决定了音色的起始特征：

1. 白噪声脉冲：产生明亮的拨弦声 $$x(n) = \begin{cases} \text{rand}(-1, 1) & n < L \\ 0 & n \geq L \end{cases}$$

2. 三角波脉冲：模拟拨弦位置 $$x(n) = \begin{cases} \frac{2n}{p} & 0 \leq n < p \\ 2 - \frac{2n}{p} & p \leq n < L \\ 0 & n \geq L \end{cases}$$ 其中$p$是拨弦位置（0到L之间）

3. 锯齿波脉冲：产生更尖锐的音头

6.1.4 衰减与滤波

低通滤波器$H(z)$模拟弦的频率相关衰减。最简单的是二点平均滤波器： $$H(z) = \frac{1 + z^{-1}}{2}$$ 频率响应： $$|H(e^{j\omega})| = |\cos(\omega/2)|$$ 这个滤波器在高频处衰减更快，符合真实弦的物理特性。更复杂的滤波器设计可以精确控制衰减特性： $$H(z) = \frac{b_0 + b_1z^{-1}}{1 + a_1z^{-1}}$$ 其中系数可以根据目标衰减时间$T_{60}$计算： $$g \cdot |H(e^{j\omega_0})|^{f_s \cdot T_{60} / L} = 0.001$$

6.2 数字波导理论

6.2.1 d'Alembert解与行波分解

一维波动方程： $$\frac{\partial^2 y}{\partial t^2} = c^2 \frac{\partial^2 y}{\partial x^2}$$ d'Alembert通解： $$y(x,t) = y^+(x-ct) + y^-(x+ct)$$ 这表明任何一维波都可以分解为向右传播的波$y^+$和向左传播的波$y^-$的叠加。

6.2.2 波导网络基础

数字波导将连续的波动方程离散化：

    y^+(n)
    ------>----[z^(-1)]------>----
                                  |
                                  v
                              散射结点
                                  |
    y^-(n)                        v
    <------<----[z^(-1)]------<----

离散化后的传播： $$y^+(n+1, m) = y^+(n, m-1)$$ $$y^-(n+1, m) = y^-(n, m+1)$$ 其中$n$是时间索引，$m$是空间索引。

6.2.3 散射结合

当多个波导在一点相遇时，需要散射结合来保证物理约束（如压力连续、流量守恒）：

对于N个波导在一点相遇，压力波的散射方程： $$p_j^- = p_J - p_j^+$$ 其中结点压力： $$p_J = \frac{2\sum_{i=1}^N \Gamma_i p_i^+}{\sum_{i=1}^N \Gamma_i}$$ $\Gamma_i$是第i个波导的特征导纳。

6.2.4 端点条件与反射

不同的端点条件产生不同的反射特性：

1. 刚性端点（完全反射）： $$p^-(n) = p^+(n)$$

2. 开放端点（相位反转）： $$p^-(n) = -p^+(n)$$

3. 有损端点（部分反射）： $$p^-(n) = r \cdot p^+(n)$$ 其中$|r| < 1$

6.3 延迟线与梳状滤波器

6.3.1 FIR梳状滤波器

前馈梳状滤波器的传递函数： $$H(z) = 1 + g \cdot z^{-D}$$ 频率响应： $$|H(e^{j\omega})| = |1 + g \cdot e^{-j\omega D}| = \sqrt{1 + g^2 + 2g\cos(\omega D)}$$ 产生周期为$f_s/D$的谱峰和谱谷。

6.3.2 IIR梳状滤波器

反馈梳状滤波器： $$H(z) = \frac{1}{1 - g \cdot z^{-D}}$$ 这正是Karplus-Strong算法的基础结构。极点位于： $$z_k = g^{1/D} \cdot e^{j2\pi k/D}$$ 当$g$接近1时，产生尖锐的共振峰。

6.3.3 全通滤波器在物理建模中的应用

全通滤波器保持幅度响应不变，只改变相位： $$H_{AP}(z) = \frac{a + z^{-1}}{1 + a \cdot z^{-1}}$$ 用途：

1. 分数延迟实现：通过调整$a$实现非整数延迟
2. 频散模拟：模拟真实介质中的频率相关传播速度
3. 相位均衡：补偿其他滤波器引入的相位失真

6.3.4 分数延迟实现

对于分数延迟$D = D_{int} + d$（其中$0 < d < 1$），可以使用：

1. 线性插值： $$H(z) = (1-d) + d \cdot z^{-1}$$

2. 全通插值（Thiran设计）： $$H(z) = \frac{a_N + a_{N-1}z^{-1} + \cdots + z^{-N}}{1 + a_1z^{-1} + \cdots + a_Nz^{-N}}$$ 其中系数通过最大平坦群延迟条件确定。

3. Lagrange插值：高阶多项式插值，精度更高但计算复杂度也更高。

6.4 弦乐器建模

6.4.1 拨弦模型

拨弦乐器（吉他、竖琴等）的完整模型包含多个组件：

   激励位置     拾音位置
      |           |
      v           v
 [激励]-->[弦模型]-->[共鸣体]-->[输出]
            |
            v
        [阻尼/衰减]

扩展的Karplus-Strong模型：

1. 拨弦位置影响： 在弦长的$\beta$位置（$0 < \beta < 1$）拨弦会抑制第$n$次谐波： $$A_n = \sin(n\pi\beta)$$ 当$\beta = 1/2$（中点）时，所有偶次谐波被抑制。

2. 拾音位置影响： 类似地，在$\alpha$位置拾音： $$H_{pickup}(z) = \frac{1}{2}[z^{-\lfloor\alpha L\rfloor} + z^{-\lfloor(1-\alpha)L\rfloor}]$$

3. 弦的刚度（频散效应）： 真实的弦有刚度，导致高频传播更快： $$H_{stiff}(z) = \frac{1-S}{1+S} \cdot \frac{1+S \cdot z^{-1}}{1 - S \cdot z^{-1}}$$ 其中$S$是刚度参数。

6.4.2 弓弦交互模型

弓弦乐器（小提琴、大提琴等）的物理过程更复杂，涉及stick-slip（粘滑）现象：

Helmholtz运动： 弓弦交互产生的理想运动是锯齿波形，包含两个阶段：

1. 粘着阶段：弦随弓一起运动
2. 滑动阶段：弦快速回弹

数字实现使用弓表查找： $$v_{string} = f_{bow}(v_{bow} - v_{string}, F_{bow})$$ 其中$f_{bow}$是非线性摩擦力函数，典型形式： $$f_{bow}(v_r, F) = F \cdot sign(v_r) \cdot \mu(|v_r|)$$ 摩擦系数$\mu$随相对速度变化： $$\mu(v) = \mu_s \cdot e^{-v/v_s} + \mu_d \cdot (1 - e^{-v/v_d})$$

6.4.3 共鸣体建模

共鸣体（如吉他箱体、小提琴面板）显著影响音色：

1. 简化模型（并联二阶滤波器组）： $$H_{body}(z) = \sum_{i=1}^{N} \frac{g_i}{1 - 2r_i\cos(\omega_i)z^{-1} + r_i^2z^{-2}}$$ 每个滤波器对应一个共振模态。

2. 卷积模型： 使用实测的箱体脉冲响应： $$y(n) = \sum_{k=0}^{M-1} h_{body}(k) \cdot x(n-k)$$

3. 模态合成方法： 基于有限元分析的精确物理模型（计算密集）。

6.4.4 参数控制与表现力

物理建模的优势在于参数的物理意义明确：

1. 力度控制： - 拨弦力度 → 激励幅度 - 弓压 → 摩擦力参数

2. 音色变化： - 拨弦/弓弦位置 → 谐波含量 - 阻尼 → 衰减时间

3. 音高控制： - 弦长（延迟线长度） - 张力（需要重新计算波速）

实时表现力映射：

MIDI速度 --> [映射曲线] --> 激励幅度
                         --> 滤波器截止频率
                         --> 衰减时间

调制轮 --> [映射] --> 弓压
                   --> 弓速
                   --> 弓位置

6.5 管乐器建模

6.5.1 单簧管模型

单簧管是最早成功建模的管乐器，其核心是非线性励振器（簧片）和线性共鸣器（管体）的耦合：

    压力 p
      ^
      |
  [簧片模型]<---+
      |         |
   流量 u       |
      v         |
  [管体波导]----+

簧片方程（准静态模型）： $$u = \begin{cases} 0 & p > p_M \\ w \cdot H \cdot \sqrt{\frac{2|p_\Delta - p|}{ρ}} \cdot sign(p_\Delta - p) & p_{thresh} < p < p_M \\ u_{max} & p < p_{thresh} \end{cases}$$ 其中：

• $p_\Delta$ = 嘴部压力 - 管内压力
• $w$ = 簧片开口宽度
• $H$ = 簧片有效高度
• $p_M$ = 簧片完全关闭的压力阈值

管体建模： 使用数字波导，考虑：

1. 音孔影响（分支波导）
2. 锥形管的阻抗变化
3. 端部辐射损耗

6.5.2 铜管乐器模型

铜管乐器的励振机制是嘴唇振动，比簧片更复杂：

嘴唇模型（二质量模型）： $$m_1\ddot{x}_1 + r_1\dot{x}_1 + k_1x_1 = F_{lip1}(p, x_1, x_2)$$ $$m_2\ddot{x}_2 + r_2\dot{x}_2 + k_2x_2 = F_{lip2}(p, x_1, x_2)$$ 简化的单模态模型： $$\ddot{x} + 2\zeta\omega_0\dot{x} + \omega_0^2x = \frac{p}{μ}$$ 其中$\omega_0$是嘴唇的自然频率，可由演奏者调节。

流量方程： $$u = \begin{cases} b \cdot x \cdot \sqrt{\frac{2p}{ρ}} & x > 0 \\ 0 & x \leq 0 \end{cases}$$

6.5.3 音孔与指法建模

音孔改变管的有效长度和辐射特性：

音孔传输矩阵： $$\begin{bmatrix} p_1 \\ u_1 \end{bmatrix} = \begin{bmatrix} 1 & 0 \\ Y_{hole} & 1 \end{bmatrix} \begin{bmatrix} p_2 \\ u_2 \end{bmatrix}$$ 音孔导纳： $$Y_{hole} = \begin{cases} jZ_0 \tan(kl_h) & \text{闭孔} \\ \frac{1}{jZ_h \tan(kl_{eff})} & \text{开孔} \end{cases}$$ 其中$l_{eff}$是考虑端部修正的有效长度。

指法表实现： 预计算不同指法组合的等效管长和共振频率，实时切换。

6.5.4 非线性效应

管乐器中的重要非线性现象：

1. 激波形成： 高声压级时的波形陡化，需要考虑非线性传播： $$\frac{\partial p}{\partial t} + c\frac{\partial p}{\partial x} + \frac{\beta p}{ρc^2}\frac{\partial p}{\partial x} = 0$$

2. 涡流损耗： 在音孔边缘和管端产生，导致额外阻尼。

3. 模态耦合： 强激励下不同振动模态之间的能量传递。

4. 多音现象： 通过控制嘴唇或簧片参数产生和弦效果。

本章小结

物理建模合成通过数学模拟真实乐器的物理行为来生成声音，提供了前所未有的真实感和表现力。本章涵盖的关键概念：

核心算法

• Karplus-Strong算法：用延迟线+低通滤波器实现弦乐合成，基频$f_0 = f_s/L$
• 数字波导：基于d'Alembert解，将波分解为前向波和后向波
• 散射结合：多个波导相遇时保证物理约束

关键公式

1. 基本KS方程：$y(n) = x(n) + g \cdot H(z) \cdot y(n-L)$
2. 波动方程解：$y(x,t) = y^+(x-ct) + y^-(x+ct)$
3. 二点平均滤波器：$H(z) = \frac{1 + z^{-1}}{2}$
4. 梳状滤波器： - FIR：$H(z) = 1 + g \cdot z^{-D}$ - IIR：$H(z) = \frac{1}{1 - g \cdot z^{-D}}$
5. 全通滤波器：$H_{AP}(z) = \frac{a + z^{-1}}{1 + a \cdot z^{-1}}$

实现要点

• 分数延迟对精确音高控制至关重要
• 激励信号设计决定音色起始特征
• 非线性元素（簧片、弓弦摩擦）产生丰富谐波
• 共鸣体建模显著影响最终音色
• 参数的物理意义使得控制直观自然

练习题

基础题

练习6.1 计算题：若采样率为48000 Hz，要用Karplus-Strong算法生成C4音（261.63 Hz），延迟线长度应该是多少？如果只能使用整数延迟，实际频率是多少？频率误差是多少音分？

练习6.2 分析题：在Karplus-Strong算法中，使用二点平均滤波器$H(z) = (1 + z^{-1})/2$。计算并画出其幅度响应，解释为什么它能模拟弦的物理特性。

练习6.3 推导题：证明FIR梳状滤波器$H(z) = 1 + g \cdot z^{-D}$的频率响应在$\omega = 2\pi k/D$处有峰值（当g>0时）或谷值（当g<0时）。

挑战题

练习6.4 设计题：设计一个改进的Karplus-Strong算法，能够模拟吉他的拨弦位置效应。在弦长1/4处拨弦会抑制哪些谐波？写出激励信号的数学表达式。

练习6.5 分析题：数字波导中，两个不同特征阻抗Z₁和Z₂的波导连接，推导反射系数和透射系数。如果Z₂ = 2Z₁，计算能量反射率。

练习6.6 实现题：描述如何用全通滤波器实现分数延迟D = 100.3样本。给出一阶Thiran全通滤波器的系数计算方法。

练习6.7 开放题：比较Karplus-Strong算法与采样回放在合成拨弦乐器声音时的优缺点。在什么情况下你会选择物理建模而不是采样？

练习6.8 综合题：设计一个简单的单簧管物理模型。给出簧片-管体耦合的差分方程，并解释如何选择参数来控制音色的明亮度。

常见陷阱与错误

1. 延迟线长度的精度问题

错误：直接将浮点延迟长度截断为整数

L = int(fs / f0)  # 音高不准！

正确：使用分数延迟或调整采样率

L_int = round(fs / f0)
frac = (fs / f0) - L_int
使用全通滤波器补偿分数部分

2. 滤波器稳定性

错误：反馈增益过大导致不稳定

g = 1.1  # 系统发散！

正确：确保所有极点在单位圆内

g × |H(z)| < 1 对所有频率

3. 激励信号设计

错误：激励时间过长，产生不自然的音头

激励持续整个音符时长

正确：激励应该是短脉冲

激励长度 ≤ 一个周期（L样本）

4. 非线性函数的不连续性

错误：簧片模型中的硬限幅

if (x < 0) u = 0;  // 产生咔嗒声

正确：使用平滑的非线性函数

u = x × smooth_step(x)  // 平滑过渡

5. 相位响应忽略

错误：只考虑幅度响应设计滤波器 正确：物理建模中相位同样重要，影响波的传播时间

6. 参数突变

错误：直接改变延迟线长度

L = new_L;  // 产生咔嗒声和不连续

正确：使用交叉淡化或参数平滑

L += (new_L - L) × smooth_factor;

7. 数值精度问题

错误：使用float32可能累积误差 正确：关键运算使用double精度，特别是相位累加器

8. 忽视物理约束

错误：设置物理上不可能的参数组合

弦张力为负、簧片质量为零等

正确：参数范围应符合物理规律，建立合理的参数映射