波动方程的求解

波动方程是描述波传播现象的基本偏微分方程,在物理、工程等领域有着广泛应用,例如弦的振动、声波、电磁波的传播等。理解波动方程的解法及其物理意义,是数学物理方程学习的核心内容之一。下面先讨论无限长弦的行波,再讨论有限长弦的驻波、高维柯西问题,以及能量估计。

达朗贝尔公式

达朗贝尔公式解法也叫行波法,直接求解整条实线上的一维波动方程初值问题。它把解写成向左、向右传播的两个波,因此特别适合看清波速和依赖区域。半无限或有限区间上的一些问题也能结合奇偶延拓、反射等方法处理,但需要把边界条件一并考虑。

弦振动方程的达朗贝尔解法

问题和方程的引入

如果我们所考察的弦线长度很长,而我们需要知道的又只是在较短时间且离开边界较远的一段范围内的振动情况,那么边界条件的影响就可以忽略,所以不妨把所考察弦线的长度视为无限。

我们先写出含外力的一般初值问题,随后利用叠加原理分别处理初值和外力: $$ \begin{cases} u_{tt} = a^2 u_{xx} + f(x, t) & (-\infty < x < +\infty, t > 0), \quad (1) \\ u(x, 0) = \varphi(x), \quad u_t(x, 0) = \psi(x) & (2) \end{cases} $$ 由于方程及初始条件都是线性的,可以把解分成 u = u1 + u2。其中 u1 只由初始位移和速度产生,满足

$$ \begin{cases} (u_1)_{tt}=a^2(u_1)_{xx},&x\in\mathbb R,\ t>0,\quad\text{(3)}\\ u_1(x,0)=\varphi(x),\quad (u_1)_t(x,0)=\psi(x).&\text{(4)} \end{cases} $$

而 u2 只由外力产生,满足

$$ \begin{cases} (u_2)_{tt}=a^2(u_2)_{xx}+f(x,t),&x\in\mathbb R,\ t>0,\quad\text{(5)}\\ u_2(x,0)=0,\quad (u_2)_t(x,0)=0.&\text{(6)} \end{cases} $$

两部分分别求出后相加,便得到原问题的解。

首先考察无外力的齐次方程及给定初值,通过自变量变换求解。

对于齐次的初值问题(方程 utt − a2uxx = 0 配合上述初始条件),开始进行详细的求解

坐标变换与方程化简

为了求解问题 (3)(4),我们引入特征线坐标变换:

$$ \begin{cases} \xi = x - at \\ \eta = x + at \end{cases} \quad \text{(7)} $$

其逆变换为:

$$ \begin{cases} x = \frac{\xi + \eta}{2} \\ t = \frac{\eta - \xi}{2a} \end{cases} \quad \text{(8)} $$

令新的未知函数为 u(ξ, η)。利用复合函数微分法则,我们可以推导出 utt 和 uxx 在新坐标系下的表达式:

首先计算一阶偏导: $$ u_x = \frac{\partial u}{\partial \xi} \frac{\partial \xi}{\partial x} + \frac{\partial u}{\partial \eta} \frac{\partial \eta}{\partial x} = u_\xi \cdot 1 + u_\eta \cdot 1 = u_\xi + u_\eta $$

$$ u_t = \frac{\partial u}{\partial \xi} \frac{\partial \xi}{\partial t} + \frac{\partial u}{\partial \eta} \frac{\partial \eta}{\partial t} = u_\xi \cdot (-a) + u_\eta \cdot a = a(u_\eta - u_\xi) $$

接着计算二阶偏导: $$ u_{xx} = \frac{\partial}{\partial x}(u_\xi + u_\eta) = \frac{\partial}{\partial \xi}(u_\xi + u_\eta) \frac{\partial \xi}{\partial x} + \frac{\partial}{\partial \eta}(u_\xi + u_\eta) \frac{\partial \eta}{\partial x} = (u_{\xi\xi} + u_{\eta\xi}) \cdot 1 + (u_{\xi\eta} + u_{\eta\eta}) \cdot 1 = u_{\xi\xi} + 2u_{\xi\eta} + u_{\eta\eta} \quad \text{(9)} $$

$$ u_{tt} = \frac{\partial}{\partial t}(a(u_\eta - u_\xi)) = \frac{\partial}{\partial \xi}(a(u_\eta - u_\xi)) \frac{\partial \xi}{\partial t} + \frac{\partial}{\partial \eta}(a(u_\eta - u_\xi)) \frac{\partial \eta}{\partial t} = a(u_{\eta\xi} - u_{\xi\xi})(-a) + a(u_{\eta\eta} - u_{\xi\eta})(a) = a^2(u_{\xi\xi} - 2u_{\xi\eta} + u_{\eta\eta}) \quad \text{(10)} $$

将 (9) 和 (10) 代入波动方程 (3) utt − a2uxx = 0: a2(uξξ − 2uξη + uηη) − a2(uξξ + 2uξη + uηη) = 0

化简得到: −4a2uξη = 0 ⇒ uξη = 0  (11)

通解的求解

之后进行通解的求解

方程 (11) 可以通过积分直接求解。先关于 ξ 积分一次,再关于 η 积分一次,得到通解为:

u(ξ, η) = F(η) + G(ξ)  (12) 其中 F 和 G 是具有二阶连续导数的任意函数。

将变量变换 (7) 代回 (12),即可得到方程 (3) 的通解: u(x, t) = F(x + at) + G(x − at)  (13)

利用初始条件确定任意函数

我们利用初始条件 (4) 来确定通解 (13) 中的 F 和 G。

将 t = 0 代入 (13) 和其对 t 的导数: u(x, 0) = F(x) + G(x) = φ(x)  (14)

ut(x, t) = aF′(x + at) − aG′(x − at)

$$ u_t(x,0) = aF'(x) - aG'(x) = \psi(x) \Rightarrow F'(x) - G'(x) = \frac{1}{a}\psi(x) \quad \text{(15)} $$

对 (15) 两边积分,得到: $$ F(x) - G(x) = \frac{1}{a}\int_{x_0}^x \psi(\alpha) d\alpha + C \quad \text{(16)} $$

现在我们有关于 F(x) 和 G(x) 的线性方程组: $$ \begin{cases} F(x) + G(x) = \varphi(x) \\ F(x) - G(x) = \frac{1}{a}\int_{x_0}^x \psi(\alpha) d\alpha + C \end{cases} $$

联立解得: $$ F(x) = \frac{1}{2}\varphi(x) + \frac{1}{2a}\int_{x_0}^x \psi(\alpha) d\alpha + \frac{C}{2} $$

$$ G(x) = \frac{1}{2}\varphi(x) - \frac{1}{2a}\int_{x_0}^x \psi(\alpha) d\alpha - \frac{C}{2} $$

将 x + at 代替 x 到 F(x) 中,将 x − at 代替 x 到 G(x) 中,再代入通解 (13): $$ u(x,t) = \left( \frac{1}{2}\varphi(x+at) + \frac{1}{2a}\int_{x_0}^{x+at} \psi(\alpha) d\alpha + \frac{C}{2} \right) + \left( \frac{1}{2}\varphi(x-at) - \frac{1}{2a}\int_{x_0}^{x-at} \psi(\alpha) d\alpha - \frac{C}{2} \right) $$

$$ u(x,t) = \frac{1}{2}[\varphi(x+at) + \varphi(x-at)] + \frac{1}{2a}\left( \int_{x_0}^{x+at} \psi(\alpha) d\alpha - \int_{x_0}^{x-at} \psi(\alpha) d\alpha \right) $$

$$ u(x,t) = \frac{1}{2}[\varphi(x+at) + \varphi(x-at)] + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \quad \text{(18)} $$

这就是无限长弦自由振动的达朗贝尔公式,或称达朗贝尔解。

这种求解方法就是达朗贝尔解法。

$$ \begin{cases} u_{tt} = a^2 u_{xx} & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(3)} $$

u(x, t) = f(x − at) + g(x + at)  (13)

达朗贝尔解的物理意义

从通解 u(x, t) = F(x + at) + G(x − at) 可见,自由弦振动是两个行波的叠加。这里写成 f(x − at)、g(x + at) 时,f, g 仅表示波形函数,与前面表示外力的 f(x, t) 不是同一个对象。

  • 函数 G(x − at) 描述的振动波形以常速度 a 向右传播,称为右传播波
  • 函数 F(x + at) 描述的振动波形以常速度 a 向左传播,称为左传播波

波动方程 (3) 中出现的常数 a 正好是波的传播速度。 这表明弦上的任意扰动总是以行波形式分别向两个方向传播出去。振动的波动性和传播的有限性是振动的两个基本特征。

考察 u1(x, t) = f(x − at)

证明其为波动方程的解:

这个函数 u1(x, t) = f(x − at) 显然是方程 (3) utt = a2uxx 的解。我们可以通过求偏导来验证这一点:

  • 对 t 求一阶偏导:$\frac{\partial u_1}{\partial t} = f'(x-at) \cdot (-a)$
  • 对 t 求二阶偏导:$\frac{\partial^2 u_1}{\partial t^2} = f''(x-at) \cdot (-a)^2 = a^2 f''(x-at)$
  • 对 x 求一阶偏导:$\frac{\partial u_1}{\partial x} = f'(x-at) \cdot (1)$
  • 对 x 求二阶偏导:$\frac{\partial^2 u_1}{\partial x^2} = f''(x-at) \cdot (1)^2 = f''(x-at)$

将这些代入方程 (3) utt = a2uxx,我们得到 a2f″(x − at) = a2f″(x − at),等式成立。

物理意义:

通过给 t 赋予不同的值,我们可以看出弦在各个时刻相应的振动状态:

  • 在 t = 0 时,我们有 u1(x, 0) = f(x)。这对应于初始时刻的振动状态(相当于弦在初始时刻各点的位移状态),即初始波形。

  • 经过时间 t0 后,在 (x, t0) 处,函数变为 u1(x, t0) = f(x − at0)。在 u − x 平面上,它相当于原来的图形 f(x) 向右平移了一段距离 at0。这意味着波形在空间中整体向右移动。

  • 随着时间的推移,这个波形图将不断地向右移动。这说明当方程 (3) 的解表示成 f(x − at) 的形式时,振动的波形是以常速度 a 向右传播。因此,由函数 f(x − at) 所描述的振动规律,称为右传播波。

考察 u2(x, t) = g(x + at)

证明其为波动方程的解:

与 f(x − at) 类似,u2(x, t) = g(x + at) 也是方程 (3) utt = a2uxx 的解:

  • 对 t 求二阶偏导:$\frac{\partial^2 u_2}{\partial t^2} = g''(x+at) \cdot (a)^2 = a^2 g''(x+at)$
  • 对 x 求二阶偏导:$\frac{\partial^2 u_2}{\partial x^2} = g''(x+at) \cdot (1)^2 = g''(x+at)$

同样满足 utt = a2uxx。

物理意义:

同样,形如 g(x + at) 的函数描述的振动波形是以常速度 a 向左传播。当 t 增加时,为了保持 x + at 不变(即保持波形上的某个特定点),x 必须减小,这表示波形向左移动。它所描述的振动规律,称为左传播波。

总结达朗贝尔解的物理意义

  • 从达朗贝尔公式 u(x, t) = f(x − at) + g(x + at) 可见,弦上的任意扰动总是以行波形式分别向两个方向传播出去。
  • 其传播速度正好是方程 (3) 中出现的常数 a。
  • 振动的波动性和传播的有限性是振动的两个基本特征。
  • 这样,弦振动方程的解为左右传播波的叠加。这体现了波动方程的线性性质,使得其解可以通过叠加不同形式的波来构建。
  • 因此,波动方程的重要特征是传播速度的有限性。

依赖区间、决定区域和影响区域

从达朗贝尔通解公式 u(x, t) = f(x − at) + g(x + at) (公式(13)) 或更具体的达朗贝尔公式: $$ u(x,t) = \frac{1}{2}[\varphi(x-at) + \varphi(x+at)] + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha)d\alpha \quad (18) $$ 我们可以深入理解波动传播的性质。

依赖区间 (Dependence Domain)

问题:初值问题 (3)(4) 的解 u(x, t) 在一点 (x, t) 的数值与初值条件在 x 轴上哪些点的值有关?

结论:从达朗贝尔公式 (18) 可以看到,解在点 (x, t) 的数值仅依赖于 x 轴的区间 [x − at, x + at] 上的初值条件,而与其他点上的初值条件无关。这个区间 [x − at, x + at] 称为点 (x, t) 的依赖区间。 $$ \begin{cases} u_{tt} = a^2 u_{xx} & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(3)} $$

$$ u(x, t) = \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \quad \text{(18)} $$

详细解释:

  • 达朗贝尔公式明确指出,在某个特定的时刻 t 和空间位置 x 处的弦位移 u(x, t) 的值,完全取决于初始位移函数 φ(x) 在 x − at 和 x + at 这两点的值,以及初始速度函数 ψ(x) 在 x − at 到 x + at 的整个区间上的积分值。
  • 这意味着,为了确定 u(x, t) 的值,我们只需要知道 t = 0 时刻 x 轴上落在 [x − at, x + at] 范围内的那些点的初始位移和初始速度。
  • 这个区间就是决定 u(x, t) 值的“过去”的范围,即该点当前状态所“依赖”的初始条件区域。

几何解释(结合下图):

image-20250722145826372
  • 该区间 [x − at, x + at] 是通过点 (x, t) 分别作斜率为 $\frac{1}{a}$ 和 $-\frac{1}{a}$ 的直线(这些直线称为特征线)与 x 轴所交截得的区间。
  • 在 x − t 平面上,从待求点 P(x, t) 向时间轴(t 轴)的反方向(即 t = 0 的 x 轴)引出两条直线:
    • 一条直线连接 P(x, t) 和 x 轴上的点 (x − at, 0),斜率为 $\frac{1}{a}$。
    • 另一条直线连接 P(x, t) 和 x 轴上的点 (x + at, 0),斜率为 $-\frac{1}{a}$。
  • 这两条特征线在 x 轴上截出的区间 [x − at, x + at] 就是点 (x, t) 的依赖区间。
  • 这形象地展示了波动传播的有限速度性质:在 t 时刻 x 处的波形,其信息来源只能是初始时刻距离其 at 范围内的点。这个区间以外的初值不会改变当前的 u(x, t)。

决定区域 (Domain of Determination)

问题:考虑初始轴 t = 0 上的一个区间 [x1, x2],其决定区域是什么?

结论:在 x − t 平面上,过点 (x1, 0) 作斜率为 $\frac{1}{a}$ 的直线 x = x1 + at,过点 (x2, 0) 作斜率为 $-\frac{1}{a}$ 的直线 x = x2 − at,这两条直线在 t > 0 的区域内相交,与 t = 0 时的区间 [x1, x2] 一起构成一个三角形区域。此三角形区域称为区间 [x1, x2] 的决定区域。 $$ \begin{cases} u_{tt} = a^2 u_{xx} & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(3)} $$

$$ u(x, t) = \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \quad \text{(18)} $$

image-20250722145921926

详细解释:

  • 决定区域是“未来”的范围。如果仅仅在初始时刻 t = 0 知道 x 轴上 [x1, x2] 区间内的初始位移和初始速度,那么在这个三角形区域内的所有点 (x, t) 的解 u(x, t) 都可以被完全确定。
  • 这是因为,对于这个三角形区域内的任意一点 (x, t),其依赖区间 [x − at, x + at] 都会完全包含在 [x1, x2] 内部。
  • 这进一步强调了波动传播的有限速度特性:在有限时间内,只有有限区域内的初始信息能够决定有限区域内的解。

影响区域 (Region of Influence)

问题:如果在初始时刻 t = 0 扰动仅在一有限区间 [x1, x2] 上存在,则经过时间 t 后,它所影响到的范围是什么呢?

我们知道,波动是以一定的速度 a 向两个方向传播

$$ \begin{cases} u_{tt} = a^2 u_{xx} & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(3)} $$

$$ u(x, t) = \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \quad \text{(18)} $$

结论:波动是以一定的速度 a 向两个方向传播。因此,经过时间 t 后,它所传到的范围(受初始扰动影响到的范围)由不等式 x1 − at ≤ x ≤ x2 + at 决定,而在此范围之外不受影响。

详细解释:

  • 影响区域是“未来”的扩展范围。与决定区域不同,影响区域关注的是一个初始扰动(例如,一个在 [x1, x2] 区间内非零的初始位移或速度)在未来时刻会影响到哪些空间点。
  • 由于波以速度 a 向左右两边传播,初始在 x1 处的扰动会向左传播,在 t 时刻到达 x1 − at;初始在 x2 处的扰动会向右传播,在 t 时刻到达 x2 + at。
  • 因此,在时刻 t,原始区间 [x1, x2] 内的所有扰动,其影响范围将扩散到更大的区间 [x1 − at, x2 + at]。

几何解释(结合下图):

image-20250722145939508
  • 在 x − t 平面上,不等式 x1 − at ≤ x ≤ x2 + at 所表示的区域称为区间 [x1, x2] 的影响区域。
  • 从 x 轴上的点 (x1, 0) 向上画一条斜率为 $-\frac{1}{a}$ 的直线(方程为 x = x1 − at)。
  • 从 x 轴上的点 (x2, 0) 向上画一条斜率为 $\frac{1}{a}$ 的直线(方程为 x = x2 + at)。
  • 这两条直线向上展开,与初始区间 [x1, x2] 一起围成的区域就是影响区域。
  • 在影响区域内,[x1, x2] 上的初值可能影响 u(x, t);区域外的点则完全不受这段初值影响。
  • 如果初始位移与速度都只在 [x1, x2] 内非零,那么在 [x1 − at, x2 + at] 之外解仍为零。这说的是非零扰动的支集至多以速度 a 向外扩展;区间内部也可能因为左右行波抵消而出现零值。

特殊情况:从一点出发的特征线

特征线示意

当区间缩到 x0 时,两条边界是 x = x0 ± at。对于通常要求连续的初始数据,单独“只改变一点的值”并不构成一个有效的局部扰动;把它理解成越来越窄的扰动区间更合适。若允许 Dirac 脉冲这样的广义初值,初始位移脉冲沿两条特征线传播,而初始速度脉冲可以影响两线之间的区域。由此也能看出,依赖区间和影响区域讲的是“可能影响”,并不保证区间内每一点的解都非零。

齐次化原理

$$ \begin{cases} u_{tt} = a^2 u_{xx} + f(x, t) & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(1)} $$

$$ \begin{cases} u_{tt} = a^2 u_{xx} & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(3)} $$

$$ \begin{cases} u_{tt} = a^2 u_{xx} + f(x, t) & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = 0,\ u_t(x, 0) = 0 \end{cases} \quad \text{(5)} $$

为了求解问题 (5),我们利用齐次化原理,把非齐次方程的求解问题化为相应的齐次方程的情况来处理,从而可以直接利用前面有关齐次方程的结果。

齐次化原理:若 w(x, t; τ) 是初值问题 $$ \begin{cases} w_{tt} = a^2 w_{xx} & (t > \tau), \\ \left. w \right|_{t=\tau} = 0,\ \left. w_t \right|_{t=\tau} = f(x,\tau) \end{cases} \quad \text{(21)} $$ 的解(其中 τ 为参数),则 u(x, t) = ∫0tw(x, t; τ)dτ  (22) 就是初值问题 (5)(6) 的解。

令 t′ = t − τ,并记 $\overline{w}(x,t';\tau) = w(x,t'+\tau;\tau) = w(x,t;\tau)$,则问题 (21) 可化为如下形式: $$ \begin{cases} \partial_{t't'}\overline{w} = a^2 \partial_{xx}\overline{w} & (t' > 0), \\ \left. \overline{w} \right|_{t'=0} = 0,\ \left. \partial_{t'}\overline{w} \right|_{t'=0} = f(x,\tau) \end{cases} \quad \text{(23)} $$

由达朗贝尔公式 (18) $$ u(x,t) = \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \quad \text{(18)} $$ 知问题 (23) 的解为 $$ \overline{w}(x,t';\tau) = \frac{1}{2a}\int_{x-at'}^{x+at'} f(\xi,\tau) d\xi $$

u(x, t) = ∫0tw(x, t; τ)dτ  (22)

$$ \overline{w}(x, t';\tau) = \frac{1}{2a}\int_{x-at'}^{x+at'} f(\xi,\tau) d\xi $$

由 t′ = t − τ,将变量还原得 $$ w(x,t;\tau) = \frac{1}{2a}\int_{x-a(t-\tau)}^{x+a(t-\tau)} f(\xi,\tau) d\xi \quad \text{(24)} $$

再将 (24) 代入公式 (22) 即得初值问题 (5)(6) 的解 $$ u(x,t) = \frac{1}{2a}\int_0^t \int_{x-a(t-\tau)}^{x+a(t-\tau)} f(\xi,\tau) d\xi d\tau \quad \text{(25)} $$

$$ u(x, t) = \frac{1}{2a}\int_0^t \int_{x-a(t-\tau)}^{x+a(t-\tau)} f(\xi,\tau) d\xi d\tau \quad \text{(25)} $$

验证 (25) 是否满足方程,需要用到含参变量积分的求导公式(Leibniz 公式): $$ \frac{d}{dx}\int_{\alpha(x)}^{\beta(x)}F(x,s)\,ds=F(x,\beta(x))\beta'(x)-F(x,\alpha(x))\alpha'(x)+\int_{\alpha(x)}^{\beta(x)}F_x(x,s)\,ds $$

事实上,由 (25) 确定的函数确是问题 (5) 的解。当 f 具有一阶连续导数时,由 (25) 式可得 $$ \begin{aligned} u_t &= \frac{1}{2}\int_0^t\left[f\bigl(x+a(t-\tau),\tau\bigr) +f\bigl(x-a(t-\tau),\tau\bigr)\right]d\tau. \end{aligned} $$

$$ \begin{aligned} u_{tt} &= f(x,t) + \frac{a}{2}\int_0^t [f_x(x+a(t-\tau),\tau) - f_x(x-a(t-\tau),\tau)] d\tau \\ u_x &= \frac{1}{2a}\int_0^t [f(x+a(t-\tau),\tau) - f(x-a(t-\tau),\tau)] d\tau \\ u_{xx} &= \frac{1}{2a}\int_0^t [f_x(x+a(t-\tau),\tau) - f_x(x-a(t-\tau),\tau)] d\tau \end{aligned} $$

于是 utt = a2uxx + f(x, t),即 (25) 满足方程 (5)。

再验证初始条件 (6)。由 (25) 式及 ut 的表达式 $$ u_t = \frac{1}{2}\int_0^t [f(x+a(t-\tau),\tau) + f(x-a(t-\tau),\tau)] d\tau $$ 可得 u|t = 0 = 0, ut|t = 0 = 0

以上证明了由 (25) 确定的函数确是初值问题 (5)(6) 的解。

由叠加原理,可得定解问题 (1)(2)

解可表示为 $$ \begin{aligned} u(x,t) &= \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \\ &\quad + \frac{1}{2a}\int_0^t \int_{x-a(t-\tau)}^{x+a(t-\tau)} f(\xi,\tau) d\xi d\tau \quad \text{(26)} \end{aligned} $$

例:求解下列初值问题 $$ \begin{cases} u_{tt} = u_{xx} + 2x & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \sin x,\ u_t(x, 0) = x \end{cases} $$

解:由公式 (26) 得 $$ \begin{aligned} u(x,t) &= \frac{1}{2}[\sin(x+t) + \sin(x-t)] + \frac{1}{2}\int_{x-t}^{x+t} \alpha d\alpha \\ &\quad + \frac{1}{2}\int_0^t \left( \int_{x-(t-\tau)}^{x+(t-\tau)} 2\xi d\xi \right) d\tau \\ &= \sin x \cos t + xt + xt^2 \end{aligned} $$

达朗贝尔解法解波动方程的总结

无限长弦自由振动问题

定解问题: $$ \begin{cases} u_{tt} = a^2 u_{xx} & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(3)} $$

达朗贝尔解: $$ u(x, t) = \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \quad \text{(18)} $$

通解形式: u(x, t) = f(x − at) + g(x + at)  (13)

方法:行波法或达朗贝尔解法

无限长弦强迫振动问题

定解问题: $$ \begin{cases} u_{tt} = a^2 u_{xx} + f(x, t) & (-\infty < x < +\infty,\ t > 0), \\ u(x, 0) = \varphi(x),\ u_t(x, 0) = \psi(x) \end{cases} \quad \text{(1)} $$

解: $$ \begin{aligned} u(x, t) &= \frac{\varphi(x-at) + \varphi(x+at)}{2} + \frac{1}{2a}\int_{x-at}^{x+at} \psi(\alpha) d\alpha \\ &\quad + \frac{1}{2a}\int_0^t \int_{x-a(t-\tau)}^{x+a(t-\tau)} f(\xi,\tau) d\xi d\tau \quad \text{(26)} \end{aligned} $$

达朗贝尔解的物理意义

image-20250722152129514

初边值问题的分离变量解法

达朗贝尔公式最直接的使用场景是无限长弦。现在换回一根有限长、两端固定的弦:从端点传回来的波会发生反射,单纯把整个实线上的行波公式搬过来,不能自动满足端点条件。分离变量法的想法,是先找出“空间形状不变、只让振幅随时间变化”的基本振动,再把这些基本振动叠加起来。

两端固定弦的定解问题

设弦长为 l,波速为 a > 0,外力暂时取零,初始位移和速度分别为 φ(x)、ψ(x):

$$ \begin{cases} u_{tt}=a^2u_{xx},&0<x<l,\quad t>0,\\ u(0,t)=u(l,t)=0,&t\ge 0,\\ u(x,0)=\varphi(x),\quad u_t(x,0)=\psi(x),&0\le x\le l. \end{cases} $$

问题的关键:方程和两端固定的边界条件都齐次。因此,只要找到一系列同时满足它们的基本解,就可以利用线性叠加去适配一般的初始数据。但“几个基本解相加还是解”不代表“任意无穷级数都可直接求导”;收敛和正则性还要单独检查。

第一步:把空间和时间分开

先尝试非零的乘积解

u(x, t) = X(x)T(t).

代入 utt = a2uxx 得 X(x)T″(t) = a2X″(x)T(t)。在 X, T 不为零的点上除以 a2XT,就有

$$ \frac{T''(t)}{a^2T(t)}=\frac{X''(x)}{X(x)}. $$

左边只依赖 t,右边只依赖 x,而等式对所有允许的 x, t 成立,所以两边必须是同一个常数。按照习惯把它写成 −λ,得到两道常微分方程:

$$ \boxed{X''+\lambda X=0},\qquad \boxed{T''+a^2\lambda T=0}. $$

选择 −λ 而不直接记作 λ,是为了以后出现正特征值时能方便地写出正弦函数。端点固定还要求

X(0)T(t) = X(l)T(t) = 0.

若 T ≡ 0,则 u ≡ 0,对求一般初值没有帮助;因此寻找非平凡解时必须有

X(0) = X(l) = 0.

于是空间部分变成一个边值问题:

$$ \begin{cases} X''(x)+\lambda X(x)=0,&0<x<l,\\ X(0)=X(l)=0. \end{cases} $$

这一步为什么重要?普通的常微分方程可对任意 λ 求解,但再加上两端为零,就只有某些离散的 λ 允许非零解。这些值叫特征值,对应的非零 X 叫特征函数;这就是这里最简单的 Sturm–Liouville 问题。

第二步:求特征值与特征函数

不要一上来就认定 λ > 0,我们把三种可能都检查一遍。

  1. λ < 0:设 $\mu=\sqrt{-\lambda}>0$,则 X″ − μ2X = 0,通解可以写为 X = Asinh (μx) + Bcosh (μx)。由 X(0) = 0 得 B = 0;再由 X(l) = 0 和 sinh (μl) > 0 得 A = 0。只有零解。

  2. λ = 0:这时 X = Ax + B。端点条件先给 B = 0,再给 Al = 0,所以仍只有零解。

  3. λ > 0:设 $k=\sqrt\lambda$,有 X = Acos (kx) + Bsin (kx)。X(0) = 0 迫使 A = 0;要求 X 非零,即 B ≠ 0,便必须满足 sin (kl) = 0。

最后一个条件等价于 kl = nπ,n = 1, 2, 3, …。因此

$$ \boxed{\lambda_n=\left(\frac{n\pi}{l}\right)^2},\qquad \boxed{X_n(x)=\sin\frac{n\pi x}{l}}. $$

n = 0 不能加入:它给出 X ≡ 0,不是一个新的振动模态。可以把 Xn 理解成“端点被固定时允许存在的第 n 种空间形状”。

第三步:求每个空间形状如何随时间振动

把 λn 代入时间方程:

$$ T_n''+\omega_n^2T_n=0, \qquad \omega_n=\frac{n\pi a}{l}. $$

它的通解为

Tn(t) = Cncos (ωnt) + Dnsin (ωnt).

相应的第 n 个基本振动是

$$ u_n(x,t)=\bigl(A_n\cos\omega_nt+B_n\sin\omega_nt\bigr) \sin\frac{n\pi x}{l}. $$

这里把空间函数前的常数吸收到 An, Bn 中。每一项都满足方程,也都自动满足 un(0, t) = un(l, t) = 0;只是单独一项通常不能同时匹配任意给定的 φ 和 ψ。

第四步:叠加各个模态,利用初始条件求系数

自然地考虑级数

$$ \boxed{u(x,t)=\sum_{n=1}^{\infty} \bigl(A_n\cos\omega_nt+B_n\sin\omega_nt\bigr) \sin\frac{n\pi x}{l}}, \qquad \omega_n=\frac{n\pi a}{l}. $$

把 t = 0 代入,首先得到初始位移的正弦级数:

$$ \varphi(x)=\sum_{n=1}^{\infty}A_n\sin\frac{n\pi x}{l}. $$

再对级数逐项求时间导数,并令 t = 0,得到初始速度的正弦级数:

$$ \psi(x)=\sum_{n=1}^{\infty}\omega_n B_n\sin\frac{n\pi x}{l}. $$

正弦函数在 [0, l] 上正交:

$$ \int_0^l\sin\frac{n\pi x}{l}\sin\frac{m\pi x}{l}\,dx =\begin{cases}0,&m\ne n,\\ l/2,&m=n.\end{cases} $$

因此,在相应的 Fourier 展开成立时,系数由下面两个积分给出:

$$ \boxed{A_n=\frac{2}{l}\int_0^l\varphi(x)\sin\frac{n\pi x}{l}\,dx}, $$

$$ \boxed{B_n=\frac{2}{l\omega_n}\int_0^l\psi(x)\sin\frac{n\pi x}{l}\,dx =\frac{2}{n\pi a}\int_0^l\psi(x)\sin\frac{n\pi x}{l}\,dx}. $$

可以从量纲再检查一次:An 是位移的量纲,Bn 也必须是位移的量纲。初始速度积分前应除以角频率 ωn;若漏掉,就无法使 ut(x, 0) = ψ(x)。

关于级数、端点相容与解的意义

上面的级数首先是一个形式解。要把它当成能逐项求二阶导数的古典解,需要保证级数和相关导数级数收敛,并且初值与固定端点相容。最基本的条件是

φ(0) = φ(l) = 0,   ψ(0) = ψ(l) = 0.

前两个来自 u(0, t) = u(l, t) = 0 在 t = 0 的取值;后两个来自对边界条件求一次时间导数。若解足够光滑,还可以对边界条件求两次时间导数,并利用方程得到 φ″(0) = φ″(l) = 0。

例如取 φ ∈ C4[0, l]、ψ ∈ C3[0, l],再要求 φ, φ″, ψ 在两个端点都为零,是一组能够保证上述级数表示古典解的充分条件。它不是唯一或最弱的条件。若初始数据没有这么光滑,仍可在适当的能量或弱解意义下讨论级数;这时不能随意逐项求导到二阶。

所以学习分离变量法时,最好把两个问题分开:怎样找到候选级数,以及什么条件让候选级数成为真正的解。

级数解的物理意义:驻波

再看一个一般项。只要 (An, Bn) ≠ (0, 0),就能把时间部分写成单个正弦函数:

$$ A_n\cos\omega_nt+B_n\sin\omega_nt =R_n\sin(\omega_nt+\theta_n), \qquad R_n=\sqrt{A_n^2+B_n^2}. $$

这里取 sin θn = An/Rn、cos θn = Bn/Rn。使用象限正确的 atan2 (An, Bn) 可确定 θn;仅写 arctan (An/Bn) 会在 Bn = 0 时失效,也可能弄错象限。于是

$$ u_n(x,t)=R_n\sin(\omega_nt+\theta_n)\sin\frac{n\pi x}{l}. $$

固定某个时间 t0,空间图形仍是正弦形,只是前面的整体系数在变化;固定弦上的一点 x0,该点做角频率 ωn 的简谐振动,其振幅与位置 x0 有关。这种空间形状不向左右整体平移的振动,就是驻波。

  • 节点是始终不动的位置。第 n 个模态的节点在 x = ml/n(m = 0, 1, …, n);其中两个端点必然是节点。
  • 腹点是空间振幅最大的地方,在 x = (2m − 1)l/(2n)(m = 1, …, n)。
  • ωn = nπa/l 是角频率;通常说的频率为 νn = ωn/(2π) = na/(2l)。n = 1 是基频,其余是高阶模态。

整体解由不同模态的驻波叠加而成,这也是分离变量法被称为驻波法的原因。它与前面的行波视角并不矛盾:有限弦上的来回反射,可在适当叠加后形成驻波。

例一:两端固定,初始位移是一种模态

取 l = 1,考虑

$$ \begin{cases} u_{tt}=a^2u_{xx},&0<x<1,\ t>0,\\ u(0,t)=u(1,t)=0,\\ u(x,0)=\sin(2\pi x),\quad u_t(x,0)=x(1-x). \end{cases} $$

首先由正交性,位移系数非常简单:A2 = 1,其余 An = 0。所以初始位移本身只激发第 2 个空间模态。但初始速度 x(1 − x) 不是一个单独的正弦模态,必须展开成级数:

$$ B_n=\frac{2}{n\pi a}\int_0^1x(1-x)\sin(n\pi x)\,dx. $$

分部积分两次可得

$$ \int_0^1x(1-x)\sin(n\pi x)\,dx =\frac{2\bigl(1-(-1)^n\bigr)}{n^3\pi^3}, $$

故

$$ B_n=\frac{4\bigl(1-(-1)^n\bigr)}{a n^4\pi^4}. $$

偶数 n 时 Bn = 0;奇数 n 时 Bn = 8/(an4π4)。最后的解可以写成

$$ \boxed{u(x,t)=\cos(2\pi at)\sin(2\pi x) +\sum_{\substack{n\ge1\\n\text{ 为奇数}}} \frac{8}{a n^4\pi^4}\sin(n\pi at)\sin(n\pi x).} $$

这个例子也说明:初始位移只含一种模态,并不代表后续解只有一种模态;初始速度同样会决定哪些模态被激发。对答案做快速检查:t = 0 时级数的时间正弦项消失;对时间求导后,第一项在 t = 0 消失,级数重新给出 x(1 − x) 的正弦展开。

例二:左端固定、右端自由

边界条件一变,空间特征值就会改变。现在令左端 x = 0 固定,右端 x = l 可以沿垂直方向自由移动,因此

u(0, t) = 0,   ux(l, t) = 0.

仍从 u = XT 出发。时间和空间方程与前面相同,但空间边值问题变成

X″ + λX = 0,   X(0) = 0,  X′(l) = 0.

当 λ < 0 时,$X=A\sinh(\sqrt{-\lambda}\,x)$,右端条件要求 $A\sqrt{-\lambda}\cosh(\sqrt{-\lambda}\,l)=0$,所以只剩零解;λ = 0 时 X = Ax,X′(l) = A = 0,也只剩零解。

当 λ > 0 时,$X=B\sin(\sqrt\lambda\,x)$,而 $X'(l)=B\sqrt\lambda\cos(\sqrt\lambda\,l)=0$。要求 B ≠ 0,于是 $\cos(\sqrt\lambda\,l)=0$,得到

$$ k_n=\frac{(2n-1)\pi}{2l},\qquad \lambda_n=k_n^2,\qquad X_n(x)=\sin(k_nx),\qquad \omega_n=ak_n, \quad n=1,2,\ldots. $$

与两端固定的 nπ/l 相比,这里出现的是半整数倍;原因就在于右端要求斜率为零,而不是位移为零。若初值为 φ, ψ,便有

$$ u(x,t)=\sum_{n=1}^{\infty} \bigl[A_n\cos(\omega_nt)+B_n\sin(\omega_nt)\bigr]\sin(k_nx), $$

$$ A_n=\frac{2}{l}\int_0^l\varphi(x)\sin(k_nx)\,dx, \qquad B_n=\frac{2}{l\omega_n}\int_0^l\psi(x)\sin(k_nx)\,dx. $$

这里仍有正交关系 ∫0lsin (knx)sin (kmx)dx = (l/2)δnm。端点相容条件也相应变为 φ(0) = 0、φ′(l) = 0,以及更高光滑性下相应的速度条件。

取一个具体初值来看:

$$ \varphi(x)=x^2-2lx,\qquad \psi(x)=3\sin\frac{3\pi x}{2l}. $$

先算位移系数。由于 knl = (2n − 1)π/2,有 cos (knl) = 0,连续分部积分可得

$$ \begin{aligned} A_n&=\frac{2}{l}\int_0^l(x^2-2lx)\sin(k_nx)\,dx\\ &=-\frac{4}{l k_n^3} =-\frac{32l^2}{(2n-1)^3\pi^3}. \end{aligned} $$

初始速度恰是第 2 个特征函数的 3 倍,所以利用正交性立即得到 Bn = 0(n ≠ 2),而

$$ B_2=\frac{3}{\omega_2}=\frac{2l}{a\pi}. $$

因此

$$ \boxed{\begin{aligned} u(x,t)={}&-\sum_{n=1}^{\infty} \frac{32l^2}{(2n-1)^3\pi^3} \cos\frac{(2n-1)\pi at}{2l} \sin\frac{(2n-1)\pi x}{2l}\\ &+\frac{2l}{a\pi} \sin\frac{3\pi at}{2l}\sin\frac{3\pi x}{2l}. \end{aligned}} $$

虽然右端叫“自由端”,它的位移不一定为零;它满足的是 ux(l, t) = 0。这个区别正是不能把两端固定的系数公式原样套用到此问题的原因。

分离变量法的步骤再梳理

  1. 根据物理条件写清楚方程、初值和边界条件,先确认边界条件是否齐次。
  2. 令 u = XT,把偏微分方程分成空间和时间两个常微分方程。
  3. 先结合边界条件求空间特征值和特征函数;这一步决定可用的模态,不能跳过 λ < 0、λ = 0 的检查。
  4. 把每个特征值代入时间方程,求对应的时间振动。
  5. 利用线性叠加写级数,再用初始条件和特征函数的正交性求系数。
  6. 最后检查端点相容、级数收敛,以及解是否真的能逐项求导到方程要求的阶数。

达朗贝尔解法着重展示波怎样传播;分离变量法则着重展示边界允许哪些频率和空间形状。同一个波动过程可以从这两个角度理解。

高维波动方程的柯西问题

一维弦的达朗贝尔公式已经把“初值怎样沿特征线传播”说得很清楚。不过薄膜振动、空间中的声波等现象不只有一个空间坐标。下面把未知函数写成 u(x, t),其中 x = (x1, …, xn),空间算子是

$$ \Delta u=\sum_{j=1}^{n}\frac{\partial^2u}{\partial x_j^2}. $$

先明确高维问题有哪些定解条件

若研究整个 ℝn,没有有限的物理边界,给出初始位移和初始速度就得到柯西问题:

$$ \begin{cases} u_{tt}-a^2\Delta u=f(\boldsymbol{x},t), &\boldsymbol{x}\in\mathbb R^n,\ t>0,\\ u(\boldsymbol{x},0)=\varphi(\boldsymbol{x}),\quad u_t(\boldsymbol{x},0)=\psi(\boldsymbol{x}), &\boldsymbol{x}\in\mathbb R^n. \end{cases} $$

若只在有界区域 Ω 中研究,则仍要在 Γ = ∂Ω 上附加边界条件。常见的三类分别是 u|Γ = g、∂nu|Γ = g,以及 (∂nu + αu)|Γ = g。此时是初边值问题;无界整空间上的初值问题与它不能直接混用。这里先讨论三维空间、没有外力的柯西问题:

$$ \begin{cases} u_{tt}-a^2\Delta u=0,&\boldsymbol{x}\in\mathbb R^3,\quad t>0,\\ u(\boldsymbol{x},0)=\varphi(\boldsymbol{x}),\quad u_t(\boldsymbol{x},0)=\psi(\boldsymbol{x}). \end{cases} $$

一维时通过 x ± at 把方程化简;三维时我们需要另一种方法,让高维问题也转化为一维波动方程。

球对称的情形:为什么会出现 ru

先假设初始数据只与到原点的距离 r = |x| 有关,于是解也取球对称形式 u(x, t) = U(r, t)。三维拉普拉斯算子作用在径向函数上时,公式是

$$ \Delta U=U_{rr}+\frac{2}{r}U_r =\frac{1}{r^2}\frac{\partial}{\partial r}(r^2U_r), \qquad r>0. $$

因此原方程变成

$$ U_{tt}=a^2\left(U_{rr}+\frac{2}{r}U_r\right). $$

多出来的 2Ur/r 是三维几何带来的,不能直接拿一维的达朗贝尔公式处理。现在设 v(r, t) = rU(r, t),直接求导得到

vrr = 2Ur + rUrr,   vtt = rUtt.

所以

$$ \boxed{v_{tt}=a^2v_{rr}},\qquad r>0. $$

这已经是一维波动方程!代价是只在半直线 r ≥ 0 上讨论,并且若原点处的 U 有界,就要有 v(0, t) = 0。初值也相应地变为

v(r, 0) = rφ(r),   vt(r, 0) = rψ(r).

可以把 rφ(r)、rψ(r) 向负半轴作奇延拓,再对 v 使用一维达朗贝尔公式,最后对 r > 0 除以 r,在 r = 0 处取极限。这个小计算说明了下面球面平均法的核心:把三维的径向部分变成一维波动。

一般初值:引入球面平均

一般的 φ(x) 与 ψ(x) 并不球对称。为此固定一个观察点 x,对以它为球心、半径为 r 的球面取平均。对任意函数 h,定义

$$ M_h(\boldsymbol{x},r) =\frac{1}{4\pi r^2} \int_{|\boldsymbol{y}-\boldsymbol{x}|=r}h(\boldsymbol{y})\,\mathrm dS_{\boldsymbol{y}}, \qquad r>0. $$

把球面写成 y = x + rω,|ω| = 1,同一个平均又能写成

$$ M_h(\boldsymbol{x},r) =\frac{1}{4\pi}\int_{S^2}h(\boldsymbol{x}+r\boldsymbol{\omega}) \,\mathrm d\omega. $$

这里 dSy = r2dω。两个公式分母不同,是因为前者在半径为 r 的球面上积分,后者在单位球面上积分;不能把 4πr2 与 4π 随意互换。

假设 h 足够光滑,则当 r → 0,球面缩到球心,得到

Mh(x, 0) = h(x).

再对单位球面公式求 r 导数:

$$ \frac{\partial M_h}{\partial r} =\frac{1}{4\pi}\int_{S^2} \nabla h(\boldsymbol{x}+r\boldsymbol{\omega})\cdot\boldsymbol{\omega} \,\mathrm d\omega. $$

相反方向的 ω 在 r = 0 时相互抵消,故 ∂rMh(x, 0) = 0。这两条初值性质马上会用到。

球面平均与拉普拉斯算子的关系

为了得到 Mu 所满足的方程,需要把球面平均的 r 导数同空间拉普拉斯联系起来。由上式和高斯公式,记 Br(x) 为半径 r 的球,有

$$ \begin{aligned} \partial_rM_h(\boldsymbol{x},r) &=\frac{1}{4\pi r^2} \int_{|\boldsymbol{y}-\boldsymbol{x}|=r} \partial_nh(\boldsymbol{y})\,\mathrm dS_{\boldsymbol{y}}\\ &=\frac{1}{4\pi r^2} \int_{B_r(\boldsymbol{x})}\Delta h(\boldsymbol{y})\,\mathrm d\boldsymbol{y}. \end{aligned} $$

再对 r 求导。上式的分母产生 −2∂rMh/r,积分上限变化产生球面上 Δh 的平均。因此

$$ \boxed{\partial_{rr}M_h+\frac{2}{r}\partial_rM_h=M_{\Delta h}}. $$

这个恒等式说明:球面平均后的函数,对半径 r 来说,也满足刚才球对称函数的那种“径向拉普拉斯”结构。

现在把 h 换成随时间变化的解 u,并记 M(x, r, t) = Mu(⋅, t)(x, r)。时间微分可以与球面平均交换,所以由 utt = a2Δu 得

$$ M_{tt}=a^2\left(M_{rr}+\frac{2}{r}M_r\right). $$

再设 w(x, r, t) = rM(x, r, t),就像球对称情形那样得到

$$ \boxed{w_{tt}=a^2w_{rr}},\qquad r>0, $$

并且

w(x, 0, t) = 0,   w(x, r, 0) = rMφ(x, r),  wt(x, r, 0) = rMψ(x, r).

原来的三维问题转成了对每个固定球心 x 的一维半直线问题。

用一维公式取回球心处的解

将 rMφ、rMψ 作奇延拓,便可把一维达朗贝尔公式用于 w。因为 M(x, 0, t) = u(x, t),而 w = rM,所以

$$ u(\boldsymbol{x},t)=\left.\frac{\partial w}{\partial r}\right|_{r=0}. $$

对延拓后的一维公式求导,再令 r = 0,可得到非常紧凑的结果:

$$ \boxed{u(\boldsymbol{x},t) =\frac{\partial}{\partial t}\left[tM_{\varphi}(\boldsymbol{x},at)\right] +tM_{\psi}(\boldsymbol{x},at).} $$

把球面平均的定义代回去,就是三维波动方程的Poisson 公式:

$$ \boxed{\begin{aligned} u(\boldsymbol{x},t) ={}&\frac{\partial}{\partial t}\left[ \frac{1}{4\pi a^2t} \int_{|\boldsymbol{y}-\boldsymbol{x}|=at} \varphi(\boldsymbol{y})\,\mathrm dS_{\boldsymbol{y}} \right]\\ &+\frac{1}{4\pi a^2t} \int_{|\boldsymbol{y}-\boldsymbol{x}|=at} \psi(\boldsymbol{y})\,\mathrm dS_{\boldsymbol{y}}, \qquad t>0. \end{aligned}} $$

t = 0 处不要把分母中的 t 直接代成零,而要回到平均形式取极限。用 Mh(x, 0) = h(x)、∂rMh(x, 0) = 0,可以验证 u(x, 0) = φ(x) 且 ut(x, 0) = ψ(x)。在 φ ∈ C3(ℝ3)、ψ ∈ C2(ℝ3) 等足够光滑的条件下,这一表达式给出古典解;唯一性还可以通过后面的能量方法理解。

两个直接代入的例子

例一:常数初始速度。设 φ ≡ 0,ψ ≡ c。常数的球面平均仍为 c,所以 u(x, t) = ct。它满足 utt = 0、Δu = 0,初始速度也正是 c。这个例子可用来检查 Poisson 公式中速度项前面的因子确实是 t。

例二:二次初始位移。设 φ(x) = |x|2,ψ ≡ 0。单位球面上 ω 的平均为零,因而

$$ M_{\varphi}(\boldsymbol{x},r) =\frac1{4\pi}\int_{S^2}|\boldsymbol{x}+r\boldsymbol{\omega}|^2\,\mathrm d\omega =|\boldsymbol{x}|^2+r^2. $$

于是

$$ u(\boldsymbol{x},t) =\frac{\partial}{\partial t} \bigl[t(|\boldsymbol{x}|^2+a^2t^2)\bigr] =|\boldsymbol{x}|^2+3a^2t^2. $$

再直接验算:utt = 6a2,而三维中 Δ|x|2 = 6,所以 utt − a2Δu = 0。通过这种简单数据检验系数,比死记分母 4πa2t 更可靠。

波传播的几何意义及适用范围

三维 Poisson 公式由以 x 为球心、半径为 at 的球面上的初始数据决定。速度项直接取 ψ 的球面平均;位移项还对半径求导,因此也涉及 φ 在该球面处沿径向的变化。若初始数据有紧支集,固定观察点 x 后,足够晚的球面会离开支集,该点的解便回到零。这是三维齐次波动方程鲜明的有限传播特点。

不要把它和一维的情况混同:一维达朗贝尔公式中的初始速度是在整个 [x − at, x + at] 上积分。二维波动方程的传播性质也与三维不同,因此不能直接把上面球面公式的 4π 改成圆周长来得到二维答案。

Poisson 公式在这里处理的是三维、整空间、齐次方程。若存在外力,可沿用前面的齐次化原理:对每个 τ,把 f(⋅, τ) 当作当时施加的初始速度,再积分这些响应。形式上,零初值的受迫解为

uf(x, t) = ∫0t(t − τ)Mf(⋅, τ)(x, a(t − τ)) dτ,

其余初值所产生的解仍由 Poisson 公式给出。若有有限区域的边界,还必须处理反射与边界条件;这时不能只靠整空间的球面平均公式。

能量不等式、波动方程解的唯一性与稳定性

前面的两种显式方法给出了许多解,但对复杂区域或一般初值,我们未必能把解写成整齐的公式。这时还可以先问:假设解存在,能不能从已知数据控制它的大小?同样的数据会不会产生两个不同的解?数据略有变化,解会不会剧烈变化?能量法正是处理这些问题的工具。

这里得到的估计通常叫先验估计:它是在“解已经存在且足够光滑”的前提下推出来的,所以它本身不负责证明存在性,却能为唯一性与稳定性提供依据。

一维两端固定弦:建立能量恒等式

先看 0 < x < l 上的初边值问题:

$$ \begin{cases} u_{tt}-a^2u_{xx}=f(x,t),&0<x<l,\quad 0<t<T,\\ u(0,t)=u(l,t)=0,\\ u(x,0)=\varphi(x),\quad u_t(x,0)=\psi(x). \end{cases} $$

定义(省去一个无关紧要的 1/2)

E(t) = ∫0l(ut(x, t)2 + a2ux(x, t)2) dx.

ut2 对应运动速度的平方,a2ux2 对应弦形变的平方。若还要写成真实的物理能量,需乘上线密度等常数;这里的 E 只作为便于估计的归一化能量。

把方程两边乘以 ut,再对 x ∈ (0, l) 积分:

∫0luttut dx = a2∫0luxxut dx + ∫0lfut dx.

第一项是 $\frac12\frac{d}{dt}\int u_t^2dx$。对含 uxx 的项按 x 分部积分:

$$ \begin{aligned} a^2\int_0^lu_{xx}u_t\,dx &=a^2\bigl[u_xu_t\bigr]_{x=0}^{x=l} -a^2\int_0^lu_xu_{xt}\,dx\\ &=a^2\bigl[u_xu_t\bigr]_{0}^{l} -\frac{a^2}{2}\frac{d}{dt}\int_0^lu_x^2\,dx. \end{aligned} $$

由于端点位移对所有时间恒为零,端点速度 ut(0, t)、ut(l, t) 也为零,边界项消失。于是得到能量恒等式:

$$ \boxed{\frac{dE}{dt}=2\int_0^lf(x,t)u_t(x,t)\,dx.} $$

这一行是后面所有结论的起点。若边界条件不齐次,[uxut]0l 一般不为零,表示边界也能向弦做功;此时不能不加说明地沿用下面的守恒式。

无外力时:能量守恒

当 f ≡ 0,能量恒等式右端为零,所以

E(t) = E(0) = ∫0l[ψ(x)2 + a2φ′(x)2]dx.

解的振动形状虽然在变,但能量并没有凭空产生或消失。固定端点不移动,也就不通过边界做功。如果以后遇到阻尼项或者会运动的端点,能量通常不再守恒,这反过来提示我们去检查方程和边界中有哪些能量输入或损耗。

有外力时:从恒等式走到能量不等式

当 f ≠ 0,外力可能使能量增加。利用 2|pq| ≤ p2 + q2,有

2|∫0lfut dx| ≤ ∫0lf2 dx + ∫0lut2 dx ≤ ∥f(⋅, t)∥L2(0, l)2 + E(t).

于是 E′(t) ≤ E(t)+∥f(⋅, t)∥22。为了看清积分因子,乘上 e−t:

$$ \frac{d}{dt}\bigl(e^{-t}E(t)\bigr) =e^{-t}\bigl(E'(t)-E(t)\bigr) \le e^{-t}\|f(\cdot,t)\|_2^2. $$

从 0 积分到 t 得

$$ \boxed{E(t)\le e^t\left[E(0)+ \int_0^te^{-s}\|f(\cdot,s)\|_2^2\,ds\right].} $$

特别地,对所有 0 ≤ t ≤ T 都有

$$ \boxed{E(t)\le e^T\left[ \int_0^l\bigl(\psi^2+a^2(\varphi')^2\bigr)dx +\int_0^T\int_0^lf(x,s)^2\,dx\,ds\right].} $$

右端只含已知的初始数据和外力;常数 eT 只依赖考察的时间范围。这就是一个典型的能量估计。常数还可以通过其他不等式改进,但“解能被数据控制”的结论才是重点。

也可以从 E′ = 2∫fut 得到另一种直观估计:$E'\le2\|f\|_2\sqrt E$,形式上除以 $2\sqrt E$ 后积分,得

$$ \sqrt{E(t)}\le\sqrt{E(0)}+ \int_0^t\|f(\cdot,s)\|_2\,ds. $$

当 E = 0 的时刻,用 $\sqrt{E+\varepsilon}$ 先推导,再令 ε ↓ 0 即可避免直接除以零。这个形式更能说明:外力在每一小段时间内给能量平方根增加多少。

为什么还要估计 u 本身

E(t) 控制了 ut 和 ux,却没有在定义中直接出现 u。对于零边界条件,可借助 Poincaré 不等式从 ux 控制 u;也可以直接由初始位移和速度来估计。记

H(t) = ∫0lu(x, t)2 dx.

求导并使用 2|uut| ≤ u2 + ut2:

H′(t) = 2∫0luut dx ≤ H(t) + ∫0lut2 dx ≤ H(t) + E(t).

同样乘以 e−t 积分,得到

H(t) ≤ et[H(0) + ∫0te−sE(s) ds].

再代入上一段的能量估计,便能同时控制 ∥u∥L2、∥ut∥L2 和 ∥ux∥L2。其中 H(0) = ∫0lφ2dx;这一项不能因为能量 E(0) 里没有它就随意丢掉。

从能量估计推出唯一性

假设 u1, u2 对应同一组 f, φ, ψ 和零端点条件。令 w = u1 − u2,则 w 满足

$$ \begin{cases} w_{tt}-a^2w_{xx}=0,\\ w(0,t)=w(l,t)=0,\\ w(x,0)=0,\quad w_t(x,0)=0. \end{cases} $$

对 w 使用能量守恒:Ew(t) = Ew(0) = 0。由于被积函数 wt2 + a2wx2 非负,只能有 wt = wx = 0。因此 w 在时空中为常数,再由初始条件 w(x, 0) = 0 知 w ≡ 0,也就是

$$ \boxed{u_1=u_2.} $$

注意这个证明没有给出“解一定存在”;它只说明,如果两个足够规则的解都存在,它们不可能不同。

从差函数理解稳定性

现在允许两组数据稍有不同:ui 对应 φi, ψi, fi(i = 1, 2)。差 w = u1 − u2 仍满足波动方程,只是数据变成差:

wtt − a2wxx = f1 − f2,  w(x, 0) = φ1 − φ2,  wt(x, 0) = ψ1 − ψ2.

把前面的估计应用到 w,便有

$$ \begin{aligned} &\|u_{1,t}(\cdot,t)-u_{2,t}(\cdot,t)\|_2^2 +a^2\|u_{1,x}(\cdot,t)-u_{2,x}(\cdot,t)\|_2^2\\ &\qquad\le e^T\left[ \|\psi_1-\psi_2\|_2^2 +a^2\|\varphi'_1-\varphi'_2\|_2^2 +\int_0^T\|f_1(\cdot,s)-f_2(\cdot,s)\|_2^2\,ds \right]. \end{aligned} $$

再配合对 Hw 的估计,也能控制 ∥u1 − u2∥2。这说明在任意固定的有限时间区间 [0, T] 内,初始位移、初始速度和外力的适当范数趋近时,解也在相应的能量范数中趋近,这就是连续依赖性(稳定性)。

二维有界区域:同一个证明如何推广

取平面上的有界区域 Ω,边界 Γ 足够规则,考察

$$ \begin{cases} u_{tt}-a^2(u_{xx}+u_{yy})=f(x,y,t),&(x,y)\in\Omega,\ t>0,\\ u|_\Gamma=0,\\ u(x,y,0)=\varphi(x,y),\quad u_t(x,y,0)=\psi(x,y). \end{cases} $$

一维的 ux2 变成二维的 |∇u|2 = ux2 + uy2,定义

EΩ(t) = ∫Ω(ut2 + a2|∇u|2) dx dy.

把方程乘以 ut 后在 Ω 上积分,用 Green 公式处理 Δu:

∫Ω(Δu)ut dx dy = ∫Γut∂nu dS − ∫Ω∇u ⋅ ∇ut dx dy.

零 Dirichlet 边界使边界上的 ut = 0,于是

EΩ′(t) = 2∫Ωfut dx dy.

和一维完全同样地使用 Cauchy–Schwarz、不等式 2pq ≤ p2 + q2 与积分因子,就有

$$ \boxed{E_\Omega(t)\le e^T\left[ \int_\Omega\bigl(\psi^2+a^2|\nabla\varphi|^2\bigr)\,dx\,dy +\int_0^T\int_\Omega f^2\,dx\,dy\,ds \right],\quad 0\le t\le T.} $$

若 f = 0,上式可加强为 EΩ(t) = EΩ(0)。对两个解的差应用该估计,立即得到古典解的唯一性;允许数据不同,则得到 H1 位移、L2 速度意义下的稳定性。这里的 H1 简单理解为函数和它的一阶空间导数都可用平方积分控制。

边界条件的提醒:上面的边界项消失依靠 u|Γ = 0。齐次 Neumann 条件 ∂nu = 0 也会让这个边界项消失,但 E 对常数位移不敏感;证明唯一性时还需用初始位移确定该常数。非齐次边界条件应保留边界做功项,或先做变量代换把边界齐次化。

整个空间中的柯西问题:收缩的影响锥

有界区域上的能量证明使用了固定边界。整空间没有固定边界,但达朗贝尔公式告诉我们传播速度是有限的。要在二维柯西问题中同样利用这一点,可围绕目标点 (x0, t0) 画一个向过去张开的特征锥。

设 x = (x, y) ∈ ℝ2,考虑齐次方程 utt − a2Δu = 0。对 0 ≤ s ≤ t0,取时刻 s 的圆盘

Ωs = {x : |x − x0| ≤ a(t0 − s)}.

s = 0 时它的半径是 at0;随着 s 接近 t0,圆盘收缩到观察点。定义圆盘内的能量

$$ \mathcal E(s)=\frac12\int_{\Omega_s} \bigl(u_t^2+a^2|\nabla u|^2\bigr)\,dx\,dy. $$

求导时有两个变化:被积函数本身随时间变化,积分区域的边界又以速度 a 向内移动。对第一部分使用方程和 Green 公式,对第二部分计算移动边界带走的能量,可得

$$ \begin{aligned} \mathcal E'(s) ={}&\int_{\Omega_s}u_t(u_{tt}-a^2\Delta u)\,dx\,dy\\ &+\int_{\partial\Omega_s} \left[a^2u_t\partial_nu -\frac a2\bigl(u_t^2+a^2|\nabla u|^2\bigr)\right]dS. \end{aligned} $$

内部积分为零。对边界积分,先用 |∂nu| ≤ |∇u|,再用 2pq ≤ p2 + q2,便有

$$ a^2u_t\partial_nu \le\frac a2u_t^2+\frac{a^3}{2}|\nabla u|^2. $$

因此 ℰ′(s) ≤ 0,即

$$ \boxed{\mathcal E(s)\le\mathcal E(0),\qquad 0\le s<t_0.} $$

这是特征锥能量不等式:向目标点收缩时,锥内能量不能凭空增加。如果两组解在底面 Ω0 的初始位移和速度完全相同,其差 w 的底面能量为零,则锥内 wt = ∇w = 0;再由底面上 w = 0 可知整锥内 w = 0。特别地,两组解在目标点的值相同。由此既得到局部唯一性,也从另一个角度证明了有限传播速度:底面圆盘以外的初始扰动无法在 t0 前影响圆锥顶点。

这里的锥底半径必须是 at0,锥在时刻 s 的半径必须是 a(t0 − s)。把半径写成 as 会变成向未来扩大的圆盘,就不再是同一个“过去依赖区域”。

柯西问题的稳定性和最后的区分

若两组初始数据只是接近而不完全相同,仍对差 w = u1 − u2 使用同一个收缩锥。由 ℰw(s) ≤ ℰw(0),锥内的速度和空间梯度的平方积分范数被底面上 ψ1 − ψ2、∇(φ1 − φ2) 控制。再利用

w(x, s) = w(x, 0) + ∫0swt(x, τ) dτ,

也能把位移差控制住。于是只要在依赖圆盘内的初始位移、初始速度及适当导数足够接近,有限时间锥内的两组解也会接近。

把这一章和前面的显式方法放在一起看,可以分清三件事:

  1. 构造解:达朗贝尔公式、分离变量级数和三维 Poisson 公式回答“解能怎样写出来”。
  2. 唯一性:对两个解的差做零数据能量估计,回答“会不会有另一个不同的解”。
  3. 稳定性:对两个解的差做非零但很小的数据估计,回答“数据改变一点,解会不会失控”。

能量估计通常不能代替构造解的存在性证明;而单纯写出一个形式级数,也不能代替收敛、边界相容和唯一性检查。二者结合,才构成对波动方程定解问题比较完整的理解。