Song's Blog

外观

边界积分表示公式与 Green 函数推导

约 1958 字大约 7 分钟

2026-04-11

做 Poisson 或 Laplace 方程时，一个很自然的问题是：能不能不在整个区域里求解，而只在边界上做文章？

Green 函数和边界积分方法的核心思路就是这件事：先用基本解把体内信息“拉”出来，再把问题尽可能压缩到边界上。

这篇文章想说明两件事：

  • 为什么解可以写成“体积分 + 边界积分”；
  • 为什么当观测点逼近边界时，会自然冒出一个 12\frac12 的系数。

先从 Poisson 方程开始：

−Δu=fin Ω- \Delta u = f \quad \text{in } \Omega

目标是把解 u(x)u(x) 写成

源项的体积分 + 边界数据的边界积分

Green 函数

Green 函数 G(x,y)G(x,y) 定义为：

−ΔyG(x,y)=δ(x−y)- \Delta_y G(x,y) = \delta(x - y)

其中：

  • δ(x−y)\delta(x-y)：Dirac delta 函数
  • xx：观测点
  • yy：源点

Green 函数表示：

在点 yy 放置一个单位点源时，在 xx 处的响应

自由空间 Green 函数（全空间无边界条件）依维度不同而不同：

二维（x,y∈R2x, y \in \mathbb{R}^2）：

G(x,y)=−12πln⁡∣x−y∣G(x,y) = -\frac{1}{2\pi} \ln |x-y|

三维（x,y∈R3x, y \in \mathbb{R}^3）：

G(x,y)=14π∣x−y∣G(x,y) = \frac{1}{4\pi |x-y|}

两者均满足 −ΔyG(x,y)=δ(x−y)-\Delta_y G(x,y) = \delta(x-y)，只是奇异性的衰减行为不同：二维为对数奇异，三维为代数奇异。

推导：Green 第二恒等式

关键工具是 Green 第二恒等式。对足够光滑的 u,vu,v，有

∫Ω(uΔv−vΔu) dy=∫∂Ω(u∂v∂n−v∂u∂n) ds\int_\Omega (u \Delta v - v \Delta u)\,dy = \int_{\partial \Omega} \left( u \frac{\partial v}{\partial n} - v \frac{\partial u}{\partial n} \right)\,ds

现在把 vv 取成 Green 函数 G(x,y)G(x,y)。由于

  • ΔyG(x,y)=−δ(x−y)\Delta_y G(x,y) = -\delta(x-y)
  • Δu=−f\Delta u = -f

左边就变成

∫Ω(−u(y)δ(x−y)+G(x,y)f(y)) dy=−u(x)+∫ΩG(x,y)f(y) dy\int_\Omega \left( - u(y)\delta(x-y) + G(x,y) f(y) \right)\,dy = -u(x) + \int_\Omega G(x,y) f(y)\,dy

右边则是标准边界项

∫∂Ω(u∂G∂n−G∂u∂n) ds\int_{\partial \Omega} \left( u \frac{\partial G}{\partial n} - G \frac{\partial u}{\partial n} \right)\,ds

把两边整理一下，就得到最常见的表示公式：

u(x)=∫ΩG(x,y)f(y) dy+∫∂Ω(G(x,y)∂u∂n(y)−u(y)∂G∂ny(x,y)) dsy\boxed{ u(x) = \int_\Omega G(x,y) f(y)\,dy + \int_{\partial \Omega} \left( G(x,y) \frac{\partial u}{\partial n}(y) - u(y) \frac{\partial G}{\partial n_y}(x,y) \right)\,ds_y }

它告诉我们，u(x)u(x) 可以拆成两部分来看：

  • 体积分：来自源项 ff
  • 边界积分：来自边界上的 Dirichlet/Neumann 数据

特别是当 f=0f=0，也就是 Laplace 方程时，体积分消失，解完全由边界贡献决定。

边界积分方程（BIE）

真正有意思的地方在下一步：如果把观测点 xx 一路推到边界上，会发生什么？

下面为了把符号写清楚，专门以三维 Laplace 方程（f=0f=0）为例。这时自由空间 Green 函数是

G(x,y)=14π∣x−y∣G(x,y)=\frac{1}{4\pi |x-y|}

表示公式变成

u(x)=∫∂ΩG(x,y)∂u∂n(y) dsy⏟单层势−∫∂Ωu(y)∂G∂ny(x,y) dsy⏟双层势u(x) = \underbrace{\int_{\partial\Omega} G(x,y)\frac{\partial u}{\partial n}(y)\,ds_y}_{\text{单层势}} - \underbrace{\int_{\partial\Omega} u(y)\frac{\partial G}{\partial n_y}(x,y)\,ds_y}_{\text{双层势}}

这两个积分看起来都在边界上，但它们在 x→x0∈∂Ωx\to x_0\in\partial\Omega 时的行为并不一样。

单层势比较温和。核函数

G(x,y)∼1∣x−y∣G(x,y)\sim \frac{1}{|x-y|}

在二维边界曲面上是弱奇异的，所以极限可以直接取到边界上，不会产生跳跃项。

双层势才是关键。它的核是

∂G∂ny(x,y)=−14π(y−x)⋅ny∣x−y∣3\frac{\partial G}{\partial n_y}(x,y) = -\frac{1}{4\pi}\frac{(y-x)\cdot n_y}{|x-y|^3}

如果只盯着这个式子看，会觉得它比单层势更危险，因为分母是三次方。真正正确的理解是：

  • 它不是简单的“不可积爆炸”；
  • 它需要按主值积分的方式理解；
  • 当 xx 从区域内部逼近边界时，会额外产生一个跳跃项。

这个跳跃项就是后面出现的 12u(x0)\frac12 u(x_0)。

挖去奇异点

为了看清这个 12\frac12 从哪来，可以在边界点 x0x_0 附近挖去一个半径为 ε\varepsilon 的小球，记

Ωε=Ω∖Bε(x0)\Omega_\varepsilon=\Omega\setminus B_\varepsilon(x_0)

这样新的边界由两部分组成：

∂Ωε=(∂Ω∖Bε(x0))∪Σε\partial\Omega_\varepsilon = (\partial\Omega\setminus B_\varepsilon(x_0))\cup \Sigma_\varepsilon

其中 Σε\Sigma_\varepsilon 是新出现的那片小半球帽。对光滑边界来说，当 ε\varepsilon 很小时，Σε\Sigma_\varepsilon 可以近似看成半个球面。

关键点在这里：12\frac12 不是来自原边界上被挖掉的那一小块，而是来自这片新出现的半球帽边界。

在完整球面上，基本解法向导数的积分给出

∫∂Bε(x0)∂G∂n ds=−1\int_{\partial B_\varepsilon(x_0)} \frac{\partial G}{\partial n}\,ds = -1

而现在只剩半个球面，所以对应贡献就是一半：

∫Σε∂G∂n ds→ε→0−12\int_{\Sigma_\varepsilon} \frac{\partial G}{\partial n}\,ds \xrightarrow{\varepsilon\to 0} -\frac12

这个 12\frac12 本质上来自局部几何中的“半个立体角”。

双层势的跳跃关系

把上面的极限过程认真做完，就得到双层势的跳跃关系。若从区域内部逼近边界点 x0x_0，则

lim⁡x→x0−∫∂Ωu(y)∂G∂ny(x,y) dsy=−12u(x0)+P.V.⁡∫∂Ωu(y)∂G∂ny(x0,y) dsy\lim_{x \to x_0^-} \int_{\partial\Omega} u(y)\frac{\partial G}{\partial n_y}(x,y)\,ds_y = -\frac{1}{2}u(x_0) + \operatorname{P.V.}\int_{\partial\Omega} u(y)\frac{\partial G}{\partial n_y}(x_0,y)\,ds_y

这里 P.V. 表示 Cauchy 主值积分。

这个公式里最容易记错的地方就是符号：对于本文采用的外法向约定和

u=S ⁣(∂u∂n)−Duu=S\!\left(\frac{\partial u}{\partial n}\right)-Du

这样的写法，内部极限对应的是 -\frac12 u，不是 +\frac12 u。

把这个极限代回表示公式，再把项移到左边，就得到边界积分方程

12u(x0)+P.V.⁡∫∂Ωu(y)∂G∂ny(x0,y) dsy=∫∂ΩG(x0,y)∂u∂n(y) dsy\boxed{ \frac{1}{2}u(x_0) + \operatorname{P.V.}\int_{\partial\Omega} u(y)\frac{\partial G}{\partial n_y}(x_0,y)\,ds_y = \int_{\partial\Omega} G(x_0,y)\frac{\partial u}{\partial n}(y)\,ds_y }

这就是最经典的一类 BIE：区域里的 PDE 被压缩成了边界上的积分方程。

从计算角度看，它的含义很直接：

  • 如果边界条件给的是 uu，那未知量更偏向于 ∂nu\partial_n u
  • 如果边界条件给的是 ∂nu\partial_n u，那未知量更偏向于 uu

总之，未知函数被搬到了边界上。

若边界在 x0x_0 处不光滑（如角点），立体角不再是 2π2\pi，系数 12\frac{1}{2} 变为依赖局部几何的 c(x0)∈(0,1)c(x_0) \in (0,1)。

二维与三维的对比

二维（R2\mathbb{R}^2）三维（R3\mathbb{R}^3）
Green 函数−12πln⁡∣x−y∣-\dfrac{1}{2\pi}\ln|x-y|14π∣x−y∣\dfrac{1}{4\pi|x-y|}
奇异性类型对数奇异代数奇异（1/r1/r）
边界 ∂Ω\partial\Omega曲线（1D）曲面（2D）
边界积分线积分面积分
BIE 未知量标量函数（定义在曲线上）标量函数（定义在曲面上）

两种情况下的总体框架是一致的：都是 Green 恒等式、主值积分、边界极限和跳跃关系。不同之处主要在于：

  • Green 函数的具体形式不同；
  • 边界的维数不同；
  • 奇异核的局部行为写出来不一样。

与有限元方法的关系

两者都是求解 PDE 的数值方法，但出发点和适用场景不同：

FEMBIE / BEM
离散对象整个区域 Ω\Omega（体网格）仅边界 ∂Ω\partial\Omega（面/线网格）
基函数局部多项式（hat function 等）Green 函数（全局、奇异）
矩阵稀疏（局部耦合）稠密（全局耦合）
适用方程几乎任意 PDE须有已知 Green 函数（线性、常系数）
非均匀介质天然支持困难
无界域需要截断/PML天然满足远场条件

FEM 更适合： 非线性问题、非均匀材料、复杂内部结构。

BEM 更适合： 无界域（声学散射、电磁）、只关心边界量（裂纹、应力集中）。

两者也可以耦合使用（FEM-BEM coupling）：近场非线性区域用 FEM，无界远场用 BEM，在人工边界上匹配。例如地震波模拟中，地下非均匀介质用 FEM，无限远传播用 BEM。

FEM 是“分而治之”——把区域切碎，局部求解；BEM 是“化体为边”——用 Green 函数把体问题压缩到边界。

版权归属：Guisong Wu

许可证：署名 4.0 国际 (CC-BY-4.0)