English

理解扩散项（一）：从二阶导数定义到离散化
以二阶导数定义作为切入点

在应用数学和物理学中最引人入胜的方程之一是 反应-扩散方程。它支配了从 化学反应 到 生物模式形成（如动物毛皮图案）再到 热传递 等广泛现象，还记得我刚开始接触的时候是真得感概自己的才疏学浅，真没想到，连动物身上的花纹模式也是可以算出来的！

刚开始看到这个方程确实非常令人迷惑，我在此记录一下我的思考过程，从 菲克第二定律、扩散项、二阶导数的定义以及如何离散化它们 开始。 描述图像的替代文本

1. 菲克第二定律与扩散项

菲克第二定律描述了扩散如何使物质浓度随时间演化。它的表达式为：

\[\frac{\partial u}{\partial t} = D \nabla^2 u\]

其中：

这个方程表明，浓度随时间的变化率 与 浓度分布的空间曲率 成正比。简单来说，扩散会将高浓度区域分散，并填补低浓度区域。

理解拉普拉斯算子

乍一看，拉普拉斯算子 可能显得令人生畏，但它仅仅是 所有空间维度上的二阶导数之和。在一维情况下，它简化为：

\[\nabla^2 u = \frac{d^2u}{dx^2}\]

在二维或三维中，它考虑了所有空间方向的变化：

\[\nabla^2 u = \frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial z^2}\]

这个算子本质上量化了一个函数与其邻近值的偏差。值得注意的是，如果某点的函数具有 负拉普拉斯值（类似峰值），其值倾向于 减小；而 正拉普拉斯值（类似谷底）则导致值 增加。这一特性使拉普拉斯算子成为建模扩散的理想数学工具，因为物质自然从高浓度区域扩散到低浓度区域。可能这里说的还是有点抽象，下面会从导数定义的角度去理解这个特性。

2. 二阶导数与扩散

要理解为什么 拉普拉斯算子（二阶导数）出现在扩散项中，让我们先看看 一阶导数 的定义：

\[\frac{du}{dx} = \lim_{\Delta x \to 0} \frac{u(x+\Delta x) - u(x)}{\Delta x}\]

二阶导数 测量了一阶导数本身的变化：

\[\frac{d^2u}{dx^2} = \lim_{\Delta x \to 0} \frac{\frac{du}{dx} (x+\Delta x) - \frac{du}{dx} (x)}{\Delta x}\]

展开一阶导数的定义，我们还可以用函数值来表达二阶导数：

\[\frac{d^2u}{dx^2} = \lim_{\Delta x \to 0} \frac{\left(\frac{u(x+\Delta x) - u(x)}{\Delta x}\right) - \left(\frac{u(x) - u(x-\Delta x)}{\Delta x}\right)}{\Delta x}\]

整理得：

\[\frac{d^2u}{dx^2} = \lim_{\Delta x \to 0} \frac{u(x+\Delta x) - 2u(x) + u(x-\Delta x)}{(\Delta x)^2}\]

如何将其与浓度的变化联系起来？ 描述图像的替代文本 我们以拉普拉斯值 大于 0 的情况为例。假设 $u(x)$ 表示特定时间 $ t $ 的浓度。在下一个时间 $t+1$，我们期望浓度从 较高浓度区域向较低浓度区域 扩散。

在上图中，我们观察到点 c 的浓度高于 b，而 b 高于 a。这意味着物质自然从 b 流向 a，从 c 流向 b。然而，由于从 c 到 b 的物质流动多于从 b 到 a 的流动,即图中的$ \Delta u_1 > \Delta u_2$，点 b 的净积累增加。

因此，在下一时间 $ t+1 $，点 b 的浓度应高于 $t$ 时的浓度。换句话说，当以下条件成立时，浓度会增加：

\[\frac{d^2u}{dx^2} > 0\]

这表明 正拉普拉斯值 为何会导致某点浓度的增加。而当$\frac{d^2u}{dx^2} < 0$时，也可以类似推出物质在该点浓度会降低。同时也可以体会到，一阶导数在浓度扩散中并没有给到足够的信息。

3. 离散化扩散项

为了数值求解扩散方程，我们使用 有限差分法 对二阶导数进行 离散化。对于 $ \frac{d^2u}{dx^2} $ 根据(2)中的推导，中心差分近似为：

\[\frac{d^2u}{dx^2} \approx \frac{u(x+\Delta x) - 2u(x) + u(x-\Delta x)}{(\Delta x)^2}\]

这个表达式表示 差分的差分，捕捉了浓度相对于其邻近点的变化。

在 数值模拟 中，我们将空间表示为网格，并使用以下公式更新值：

\[\frac{\partial u}{\partial t} \approx \frac{u^{n+1}_i - u^n_i}{\Delta t} = D \frac{u^n_{i+1} - 2u^n_i + u^n_{i-1}}{(\Delta x)^2}\]

其中 $ u^n_i $ 表示时间步 $ n $ 时网格点 $ i $ 的浓度。这使我们能够随时间计算扩散。

离散化的伪代码

initialize u[x]  # 在空间网格上定义浓度分布

for each time step:
    for each spatial point i (excluding boundaries):
        # du/dt = D * (u(i+1) - 2u - u(i-1)) / dx^2
        # （ u_new[i] - u[i] ）/ dt = D * (u[i+1] - 2*u[i] + u[i-1]) / dx^2
        u_new[i] = u[i] + D * dt / dx^2 * (u[i+1] - 2*u[i] + u[i-1])

    update u with u_new  # 进入下一时间步

4. 添加反应项

为了模拟化学反应与扩散的耦合，我们引入反应项 $ R(u) $：

\[\frac{\partial u}{\partial t} = D \frac{\partial^2 u}{\partial x^2} + R(u)\]

反应项 $ R(u) $ 代表以下过程：

这使得反应-扩散方程成为多个科学领域的强大工具。

5. 完整的反应-扩散方程及其应用

在复杂系统中，多种物质相互作用，导致耦合方程的出现。完整的反应-扩散系统通常如下：

\[\frac{\partial u}{\partial t} = D_u \nabla^2 u + f(u, v)\] \[\frac{\partial v}{\partial t} = D_v \nabla^2 v + g(u, v)\]

其中：

该系统可解释模式形成现象，例如：

Key takeaways

*****

☕ Happy learning journey~ 🛠️