Keller-Segel 趋化性模型保结构有限元方法的超收敛性分析

《Communications in Nonlinear Science and Numerical Simulation》：Superconvergence analysis of a structure-preserving finite element method for the Keller-Segel chemotaxis model

　　

**摘要**

本文提出了一种针对Keller-Segel (KS)趋化模型的结构保持有限元(FE)方法。首先，通过对数变换导出一个等价系统，以确保细胞密度的正性。基于该重构，提出了一种线性化的Crank-Nicolson (CN)有限元格式。随后，利用误差拆分技术和数学归纳法建立了数值解的有界性。进一步严格推导了l∞(L2)和l∞(H1)范数下的无条件最优误差估计以及超接近性结果。最后，通过插值后处理技术获得了超收敛结果。数值算例验证了理论发现并展示了所提格式的有效性。

**引言**

趋化是细胞、细菌及其他生物体响应环境中化学浓度梯度而进行的定向运动。这一基本机制主导着伤口愈合、免疫应答和神经元迁移等关键生物学过程。KS模型作为最简单的趋化模型，由Evelyn F. Keller和Lee A. Segel于1970年提出，其形式如下：

$$\begin{cases} \rho_t - \Delta\rho + \nabla \cdot (\rho \nabla v) = 0, & \text{在} \ \Omega \times (0, T], \\ v_t - \Delta v + v = \rho, & \text{在} \ \Omega \times (0, T], \\ \frac{\partial v}{\partial n} = 0, \quad \frac{\partial \rho}{\partial n} - \rho \frac{\partial v}{\partial n} = 0, & \text{在} \ \partial\Omega \times (0, T], \\ \rho|_{t=0} = \rho_0, \quad v|_{t=0} = v_0, & \text{在} \ \Omega. \end{cases}$$

在该系统中，$\Omega \subset \mathbb{R}^2$ 是具有Lipschitz边界 $\partial\Omega$ 的有界矩形区域，$\mathbf{n}$ 是 $\partial\Omega$ 上的外法向单位向量。变量 $\rho(x,t)$ 和 $v(x,t)$ 分别表示细胞密度和细胞分泌的化学物质的浓度，而输运项 $\nabla \cdot (\rho \nabla v)$ 表示 $\rho$ 对 $v$ 的趋化敏感性。

经典KS模型及其变体的系统性综述见文献[1]，涵盖了整体存在性和有限时间爆破等基本理论问题。与此同时，文献[2]发展了一个统一的数学框架，将生物学基础与解行为的严格分析联系起来。除这些一般理论之外，在具体模型特征的理解方面也取得了相当大的进展：文献[3]建立了三维Patlak-Keller-Segel系统的整体适定性；文献[4]为含旋转项的二维模型识别了区分爆破与有界性的尖锐临界条件；文献[5]证明了具有奇异敏感性的两物种系统的整体存在性和物种持续性；文献[6]为血管生成的趋化-对流模型推导了先验有界性理论。这些研究从多种分析和建模视角显著深化和拓展了趋化动力学的理论图景。

KS趋化模型具有三个基本物理原理：正性保持、质量守恒和自由能耗散。由于这些性质对于正确捕捉解的定性行为至关重要，因此人们投入了大量努力设计结构保持数值方法，包括有限差分方法[7]-[13]、有限元方法[14]-[21]和有限体积方法[22]-[28]。为了保持非负性和有界性，已引入了质量凝聚[14,15]和对数变换[16]等技术。为了保证能量稳定性，发展了标量辅助变量(SAV)方法[17]和不变能量二次化(IEQ)方法[18]。除了保持这些内在结构外，准确捕捉有限时间爆破现象等奇异行为对于趋化模型的数值模拟也至关重要。在这方面，文献[19]为KS模型发展了一种内部罚不连续Galerkin方法，并严格建立了完全离散格式的稳定性和收敛性。

除结构保持外，KS模型有限元格式的严格误差分析也引起了相当多的关注。文献[29]提出了一种后向Euler (BE)有限元格式，严格建立了正性保持、质量守恒和误差估计。在文献[30]中，作者首先对连续问题进行了详细分析，然后发展了一个时间半离散格式以及对化学方程的混合有限元逼近。基于该混合有限元框架，文献[31]进一步给出了趋化-排斥模型两种离散化的详细误差估计和比较研究，即标准BE协调有限元格式和混合有限元格式。此外，针对KS模型及其与流体方程耦合系统的完全离散有限元格式的误差分析已在[32]-[39]中得到研究。为了实现更高阶精度，超收敛技术也在KS型系统中得到了探索。在文献[40]中，为KS模型提出了一种线性化Crank-Nicolson有限元方法，建立了超接近性性质，并通过插值后处理技术进一步获得了全局超收敛性。然而，该工作仅关注KS模型的原形式，该形式无法保持正性。

最近，基于对数变换，文献[41]提出了一种无条件正性保持、质量守恒和能量稳定的完全离散BE有限元格式，并严格建立了L2范数下的最优误差估计。该框架随后在文献[42]中扩展至一般的k阶后向差分(BDFk)公式。然而，仍有几个问题需要进一步研究。首先，[41]中未建立最优的H1范数误差估计，而后者对于描述数值解的梯度行为至关重要。其次，[42]中的高阶收敛分析主要涉及时间精度，空间和时间两方面完整的高阶最优误差估计仍然缺乏。此外，CN方法是经典的二阶时间步进格式，但将其纳入[41,42]的正性保持框架仍探索不足。特别是，[41,42]中开发的误差估计能否推广至CN格式值得进一步研究。

为填补这些空白，本文提出并分析了一类基于[41,42]正性保持框架的完全离散CN有限元格式。所提格式无条件正性保持、质量守恒且能量稳定。此外，我们为它们建立了完整的收敛理论。本文的主要贡献有两方面：

**(a) L2和H1范数下的最优误差估计。** 我们为所提格式推导了严格的先验误差估计，并证明了l∞(L2)范数下的最优收敛性。更重要的是，将分析扩展至l∞(H1)范数，为这类正性保持格式得出了无条件的最优H1误差估计。

**(b) 通过后处理实现超收敛。** 我们进一步建立了超收敛结果，并利用插值后处理技术将数值解的收敛阶提升一阶。这一改进无需增加额外的计算复杂度或存储需求。

在全文中，一般正常数$C$独立于网格尺寸$h$和时间步长$\tau$，且在不同出现处可能取不同的值。

**章节片段：KS趋化模型的结构保持有限元格式**

设$W^{m,p}$、$H^m$和$L^2$为标准Sobolev空间，$W^{m,p}$和$H^m$上的范数分别记为$\|\cdot\|_{m,p}$和$\|\cdot\|_m$，$L^2$上的内积记为$(\cdot, \cdot)$。对于任何Banach空间$X$，我们在时间区间$[0,T]$上定义相应的$L^p(X)$和$l^p(X)$空间，详细定义见[43]。将$\rho = \exp(u)$代入方程(1)得：

$$\begin{cases} \rho = \exp(u), & \text{在} \ \Omega \times (0, T], \\ v_t - \Delta v + v = \rho, & \text{在} \ \Omega \times (0, T], \\ u_t - \Delta u = \nabla u \cdot \nabla u - \nabla u \cdot \nabla v - \Delta v, & \text{在} \ \Omega \times (0, T], \\ \frac{\partial v}{\partial n} = \frac{\partial u}{\partial n} = 0, & \text{在} \ \partial\Omega \times (0, T], \\ v|_{t=0} = v_0, \quad u|_{t=0} = u_0 = \ln(\rho_0), & \text{在} \ \Omega. \end{cases}$$

设$M = H^1(\Omega)$

**收敛性与超收敛分析**

本节给出数值解在l∞(L2)和l∞(H1)范数下的误差估计。令

$$\sigma \varepsilon_k = \varepsilon_k - G_k, \quad \tilde{\sigma\rho}_k = \rho_k - \tilde{P}_k,$$

$$e^{\varepsilon h}_k = G_k - G^h_k = G_k - R_h G_k + R_h G_k - G^h_k := \eta^{\varepsilon h}_k + \theta^{\varepsilon h}_k,$$

$$e^{\tilde{\rho h}}_k = \tilde{P}_k - \tilde{P}^h_k = \tilde{P}_k - R_h \tilde{P}_k + R_h \tilde{P}_k - \tilde{P}^h_k := \eta^{\tilde{\rho h}}_k + \theta^{\tilde{\rho h}}_k,$$

其中$k = 1, 2, \cdots, N$，$\varepsilon = u, v, \rho$，$G = U, V, P$。$U_k$和$V_k$分别是$u_k$和$v_k$的时间半离散近似，$(P_k, \tilde{P}_k)$是$\rho_k$的时间半离散近似。

为后续误差分析作如下假设：

**假设1** 假设……

**作者贡献声明（CRediT）**

Kesen Zhao：撰写—初稿、软件、数据整理。Minghao Li：撰写—审阅与编辑、方法论、概念化。Liuchao Xiao：指导、资助获取。

**利益冲突声明**

作者声明不存在已知的可能影响本文工作的竞争性经济利益或个人关系。

Kesen Zhao | Minghao Li | Liuchao Xiao

河南工业大学 数学与统计学院，郑州，450001，中国