肿瘤增长数学模型

曲面偏微分方程的一个实际应用——肿瘤增长

简介

基本想法是利用表面有限元数值方法求解演化曲面上的反应扩散方程。不断增长的生物表面上的模式形成，是通过求解曲面上的反应扩散方程得到的，可以参考图灵的经典文章 ¹。这种方程表现出空间均匀结构的扩散驱动不稳定性，导致空间不均匀的图案。这项工作的进一步的发展，有了更多的数学模型，其中的一些应用包括九头蛇的图案形成、脊椎动物的连续图案形成、动物的皮肤表面的颜色标记、蝴蝶翅膀上的色素沉着图案等等。可以参考综述文章 ²。

图灵的原始工作没有考虑几何的增长和变化。然而，这些因素在自然界中所观察到的模式的发展中很重要。增长的区域上的模式形成已被广泛研究。例如 Crampin 等人 ³ 考虑一维区域的区域增长，并表明它可能是一种机制，用于增加模式形成时，相对于初始条件的稳健性。在了解在一种名叫 Pomacanthus 的鱼中观察到的条纹进化的背景下⁴，二维增长平面区域被考虑⁵。也可以看看⁶，其中图灵扩散驱动的不稳定性条件被推广到具有缓慢、各向同性增长的反应扩散系统。本文中提出的数值方法是一种在演化表面上探索此类反应动力学的工具。

许多有趣的工作都是关于增长平面区域的。然而，在生物的应用中，很自然地考虑在不断发展的空间三维区域的表面上形成图案，例如我们参考⁷，和⁸等，其中研究了曲率
对固定球体和生长锥体的图案形成的影响。另一项初步研究是在这篇文章 ⁹ 中，它提出了一种数值方法，用于解决球形表面图案形成问题，并应用于实体肿瘤生长。他们的数值方法基于参数化球体上的表面，然后使用线方法。

下面提到的方法基于表面有限元方法，它是有限元方法⁶的自然延伸，并且能够处理复杂的几何形状和形状^101112131415。这个想法是对表面进行三角剖分，并使用基于三角剖分的分段线性表面有限元空间来逼近偏微分方程组。在这种情况下，三角剖分的顶点以规定的或受某些进化规律支配的速度移动。为了做到这一点，我们需要在表面给出一个适当的守恒定律。当以适当的变分形式给出时，演化表面有限元方法利用了该守恒定律的特殊特征。

接下来，我们要介绍的内容陈述如下。

首先，我们给出演化表面上的反应扩散方程。在这个框架内，我们给出了表面梯度的符号表达，这是演化表面上的第一原理反应扩散系统推导出来的关键。我们考虑经过充分研究的活化剂耗尽衬底模型¹⁶¹⁷¹⁸，它也被称为 Brusselator 模型。同样，在此框架中可以轻松处理任何其他合理的反应动力学。我们还讨论了表面生长模型的建模。

其次，我们看看一个演化表面有限元方法，该方法应用于改变形状的演化表面上的反应扩散系统。变分公式的一个特点是不会出现某些几何量，例如曲面的平均曲率和法线。然后通过有限元方法以自然的方式利用这一点。在这种方法中，在导出反应扩散系统或利用表面有限元方法时不需要转换或参数化。在空间中离散会产生一个非自治常微分方程组，该方程组在时间上是离散的。特别是因为这个应用涉及半线性抛物方程，我们线性化非线性动力系统[^马兹瓦缪斯]，以获得线性代数系统。然后使用 GMRes 算法¹⁹求解。

最后，在固定和进化的表面上呈现出图灵模式图案。我们展示了表面有限元方法的普遍适用性，该方法具有处理生物模型可能导致的各向异性生长的能力。结果证实，结合区域增长增强了模式选择过程。该方法也适用于将表面演化与表面上的反应扩散系统耦合的模型。这在再下一个章节中有说明，它应用于模拟实体瘤的生长。

部分结论性意见在最后一个章节。其他也考虑表面偏微分方程数值解的文章的例子包括^2021222324。这篇文章²⁵ 对计算演化表面的不同方法的回顾。

演化表面上反应扩散方程的推导

曲面梯度

我们假设 $\Gamma$ 可以全局表示为函数 $d$ 的零水平集，即存在开集$\mathcal{U}\subset {\mathbb{R}}^3$和函数$d\in C^2(\mathcal{U})$，使得
U ∩ Γ = { x ∈ U : d ( x ) = 0 } ,  and  ∇ d ( x ) ≠ 0 , ∀ x ∈ U ∩ Γ \mathcal{U} \cap \Gamma=\{x \in \mathcal{U}: d(x)=0\}, \quad \text { and } \quad \nabla d(x) \neq 0, \quad \forall x \in \mathcal{U} \cap \Gamma U∩Γ={x∈U:d(x)=0}, and ∇d(x)​=0,∀x∈U∩Γ
我们定义单位法向量场为，
$$\mathbf{v}(x)=\frac{\nabla d(x)}{|\nabla d(x)|}$$
我们定义切向梯度为，
$${\nabla }_{\Gamma }\eta (x)=\nabla \eta (x)-\nabla \eta (x)\cdot \mathbf{v}(x)\mathbf{v}(x),x\in \Gamma$$
它只依赖于 $\eta$ 在曲面上的值。它具有三个分量，表示为：
$${\nabla }_{\Gamma }\eta (x)=\left({\underline{D}}_1\eta (x),{\underline{D}}_2\eta (x),{\underline{D}}_3\eta (x)\right)$$
Laplace-Beltrami 算子就定义为切向梯度的切向散度：
$${\Delta }_{\Gamma }\eta (x)={\nabla }_{\Gamma }\cdot {\nabla }_{\Gamma }\eta (x)=\sum _{i=1}^3{\underline{D}}_i{\underline{D}}_i\eta (x)$$

演化表面上的反应扩散系统

让 $\Gamma (t)$ 是嵌入到三维空间 $\Omega (t)$ 的二维曲面。曲面上每一点的速度表示为，
$$\mathbf{v}:=V\mathbf{v}+{\mathbf{v}}_T$$

u = { u i } i = 1 m \boldsymbol{u}=\left\{u_{i}\right\}_{i=1}^{m} u={ui​}i=1m​ 是个矢量场，$\mathcal{R}(t)$ 是曲面的某一部分，由物质守恒定律，

$$\frac{d}{dt}{\int }_{\mathcal{R}(t)}{u}_i=-{\int }_{\partial \mathcal{R}(t)}{\mathbf{q}}_i\cdot \mathbf{\mu }+{\int }_{\mathcal{R}(t)}{f}_i(\mathbf{u})$$

这里边的 ${\mathbf{q}}_i$ 表示曲面通量，$f_i(\mathbf{u})$ 表示曲面内的净生产率。

$$\frac{d}{dt}{\int }_{\mathcal{R}(t)}{u}_i=-{\int }_{\partial \mathcal{R}(t)}{\mathbf{q}}_i\cdot \mathbf{\mu }+{\int }_{\mathcal{R}(t)}{f}_i(\mathbf{u})$$

通过分部积分公式，
$${\left({\nabla }_{\Gamma }\cdot \mathbf{w},v\right)}_{\Gamma }={\left({\mathbf{n}}_{\partial \Gamma }\cdot \mathbf{w},v\right)}_{\partial \Gamma }-{\left(\mathbf{w},{\nabla }_{\Gamma }v\right)}_{\Gamma }+{\left(\mathbf{w},H{\mathbf{n}}_{\Gamma }v\right)}_{\Gamma }$$

我们得到，
$${\int }_{\partial \mathcal{R}(t)}\mathbf{q}\cdot \mathbf{v}={\int }_{\mathcal{R}(t)}{\nabla }_{\Gamma }\cdot \mathbf{q}+{\int }_{\mathcal{R}(t)}\mathbf{q}\cdot \mathbf{v}H={\int }_{\mathcal{R}(t)}{\nabla }_{\Gamma }\cdot \mathbf{q}$$

同样，由经典的雷诺输运公式，
$$\frac{\mathrm{d}}{\mathrm{d}t}{\int }_{\Omega (t)}f\mathrm{d}\mathrm{x}={\int }_{\Omega (t)}\dot{f}+f\mathrm{div}\mathbf{v}\mathrm{d}\mathrm{x}$$

我们有，
$$\frac{d}{dt}{\int }_{\mathcal{R}(t)}\eta ={\int }_{\mathcal{R}(t)}{\partial }^{\bullet }\eta +\eta {\nabla }_{\Gamma }\cdot \mathbf{v}$$

联立上述结果，我们有，
$${\int }_{\mathcal{R}(t)}{\partial }^{\bullet }{u}_i+u_i{\nabla }_{\Gamma }\cdot \mathbf{v}+{\nabla }_{\Gamma }\cdot {\mathbf{q}}_i={\int }_{\mathcal{R}(t)}{f}_i(\mathbf{u})$$

由 $\mathcal{R}(t)$ 的任意性，我们得到，
$${\partial }^{\bullet }{u}_i+u_i{\nabla }_{\Gamma }\cdot \mathbf{v}+{\nabla }_{\Gamma }\cdot {\mathbf{q}}_i=f_i(\mathbf{u}),\text{ on }\Gamma (t)$$

我们假设 $\mathcal{D}$ 是一个对角的扩散张量矩阵，即${\mathcal{D}}_{ij}=d_i{\delta }_{ij}$，并且假设化学物质根据菲克定律扩散，即
$${\mathbf{q}}_i=-{\mathcal{D}}_{ij}{\nabla }_{\Gamma }{u}_j$$
那么，方程就变成了，
$${\partial }^{\bullet }{u}_i+u_i{\nabla }_{\Gamma }\cdot \mathbf{v}={\nabla }_{\Gamma }\left({\mathcal{D}}_{ij}{\nabla }_{\Gamma }{u}_j\right)+f_i(\mathbf{u}),\text{ on }\Gamma (t)$$
写成向量的形式，演化曲面上的反应扩散方程系统可以写为，
$${\partial }^{\bullet }\mathbf{u}+\mathbf{u}{\nabla }_{\Gamma }\cdot \mathbf{v}=D{\Delta }_{\Gamma }\mathbf{u}+\mathbf{f}(\mathbf{u})$$
采用的初值条件是，
$$\mathbf{u}(\cdot ,0)={\mathbf{u}}_0(\cdot )\text{ on }\Gamma (0)$$

注意到，这里的速度场可能是依赖于曲面本身，也可能是依赖于反应扩散方程组的解。

模型说明

为了简单起见，我们考虑一个非常简单的模型，叫做 Brusselator 模型（活化剂耗尽的衬底模型）。假设 $\mathbf{u}=(u,w)$ 是二维的，表示两种化学物质的浓度。

$$\mathbf{f}(\mathbf{u})={\left(f_1,f_2\right)}^T={\left(\gamma \left(a-u+u^{2}w\right),\gamma \left(b-u^{2}w\right)\right)}^T$$

对应于 $\mathbf{f}(\mathbf{u})=0$，有唯一的解，
$$u=a+b,w=b/(a+b{)}^2$$
这刚好对应了不动曲面的唯一稳态。对于给定的模型，可能存在一组参数，其值属于图灵空间，对于这些参数，固定取域上的模型由于空间均匀结构的扩散驱动不稳定性而展现出空间模式的形成。区域增长增加了生物形态发生对的范围，与固定区域相比，它们更有可能诱发图灵模式。

专注于曲面演化的三种增长模型：

• 各向同性增长
  $$X(t)=\rho (t)X_0$$
• 各向异性增长，我们用水平集函数的零水平集来刻画曲面增长
  $$\varphi (x,t)=x_1^2+x_2^2+a(t{)}^{2}G\left(x_3^2/L(t{)}^2\right)-a(t{)}^2$$
  通过规定不同的$G$、$a$ 和 $L$ 能获得不同的曲面增长。
• 化学驱动的表面演化
  我们设想，把曲面的生长过程与物质的在表面浓度结合起来，可以得到化学物质驱动的曲面演化，
  $$V=g(\mathbf{u},\mathbf{v},H)$$
  这里的 $H$ 是平均曲率，约定法向指向外时，平均曲率为正。下面只考虑速度沿着法向运动。

曲面有限元方法

下面考虑在演化曲面上，使用曲面有限元方法求解反应扩散系统。

变分形式

对于，
$${\partial }^{\bullet }{u}_i+u_i{\nabla }_{\Gamma }\cdot \mathbf{v}=d_i{\Delta }_{\Gamma }{u}_i+f_i(\mathbf{u})$$
使用标准的有限元变分，可以得到，
$$\begin{array}{rl}{\int }_{\Gamma (t)}{f}_i(\mathbf{u})\phi & ={\int }_{\Gamma (t)}{\partial }^{}{u}_i\phi +u_i\phi {\nabla }_{\Gamma }\cdot v-d_i{\int }_{\Gamma (t)}\phi {\Delta }_{\Gamma }{u}_i\\ & ={\int }_{\Gamma (t)}{\partial }^{\bullet }{u}_i\phi +u_i\phi {\nabla }_{\Gamma }\cdot v+d_i{\int }_{\Gamma (t)}{\nabla }_{\Gamma }{u}_i\cdot {\nabla }_{\Gamma }\phi -{\int }_{\partial \Gamma (t)}\phi {\nabla }_{\Gamma }{u}_i\cdot \mu \end{array}$$
不妨假设曲面是没有边界的，那么上面的最后一项就可以丢掉了，进一步利用输运定理，
$$\begin{array}{rl}{\int }_{\Gamma (t)}{f}_i(\mathbf{u})\phi & ={\int }_{\Gamma (t)}{\partial }^{\bullet }{u}_i\phi +u_i\phi {\nabla }_{\Gamma }\cdot \mathbf{v}+d_i{\int }_{\Gamma (t)}{\nabla }_{\Gamma }{u}_i\cdot {\nabla }_{\Gamma }\phi \\ & ={\int }_{\Gamma (t)}{\partial }^{\bullet }\left(u_i\phi \right)-u_i{\partial }^{\bullet }\phi +u_i\phi {\nabla }_{\Gamma }\cdot \mathbf{v}+d_i{\int }_{\Gamma (t)}{\nabla }_{\Gamma }{u}_i\cdot {\nabla }_{\Gamma }\phi \\ & =\frac{d}{dt}{\int }_{\Gamma (t)}{u}_i\phi -{\int }_{\Gamma (t)}{u}_i{\partial }^{\bullet }\phi +d_i{\int }_{\Gamma (t)}{\nabla }_{\Gamma }\phi \cdot {\nabla }_{\Gamma }{u}_i\end{array}$$
变分形式就变成了，寻找 $u_i\in H^1(\Gamma (t))$，使得
$$\begin{array}{rl}& \frac{d}{dt}{\int }_{\Gamma (t)}{u}_i\phi -{\int }_{\Gamma (t)}{u}_i{\partial }^{\bullet }\phi +d_i{\int }_{\Gamma (t)}{\nabla }_{\Gamma }\phi \cdot {\nabla }_{\Gamma }{u}_i\\ & ={\int }_{\Gamma (t)}{f}_i(\mathbf{u})\phi ,\forall \phi \in H^1(\Gamma (t))\end{array}$$

演化曲面有限元方法

对于参数化的有限元方法，我们选择三角逼近面片的顶点的运动来刻画曲面的运动，
$${\dot{X}}_j(t)=v\left(X_j(t),t\right)$$
有限元空间表示如下，
$$S_h(t)=\left\{\varphi \in C^0\left({\Gamma }_h(t)\right):{\varphi |}_{T_k}\text{ is linear affine for each }{T}_k\in {\mathcal{T}}_h(t)\right\}$$
我们用 { χ j ( ⋅ , t ) } j = 1 N \left\{\chi_{j}(\cdot, t)\right\}_{j=1}^{N} {χj​(⋅,t)}j=1N​ 来表示移动的节点基函数。

我们定义离散的物质速度，将其表示为节点速度的线性插值，即
$${\mathbf{v}}_h=\sum _{j=1}^N{\dot{X}}_j(t){\chi }_j$$
并且定义离散的物质导数，
$${\partial }_h^{\bullet }\varphi ={\varphi }_t+v_h\cdot \nabla \varphi$$
使用标准代数离散，我们可以得到半离散的代数系统，
$$\frac{d}{dt}\left(\mathcal{M}(t){\mathbf{\alpha }}_i\right)+d_i\mathcal{S}(t){\mathbf{\alpha }}_i={\mathbf{F}}_i(t)$$
其中，
$$\mathcal{M}(t{)}_{jk}={\int }_{{\Gamma }_h(t)}{\chi }_j{\chi }_k$$
$$\mathcal{S}(t{)}_{jk}={\int }_{{\Gamma }_h(t)}{\nabla }_{{\Gamma }_h}{\chi }_j\cdot {\nabla }_{{\Gamma }_h}{\chi }_k$$
$${\mathbf{F}}_{ij}={\int }_{{\Gamma }_h(t)}{f}_i(\mathbf{U}){\chi }_j$$

求解这样一个代数系统，本质上就得到了原问题的离散逼近。

时间离散

考虑最简单的两种物质浓度，$\mathbf{u}=(u,w)$，它的时间离散表示为$\left(U^n,W^n\right)$，令$d_1=1,d_2=d$，则时间离散可以写为，
$$\{\begin{array}{l}\frac{1}{\tau }{\int }_{{\Gamma }_h^{n+1}}{U}^{n+1}{\chi }_j^{n+1}+{\int }_{{\Gamma }_h^{n+1}}{\nabla }_{{\Gamma }_h^n}{U}^{n+1}\cdot {\nabla }_{{\Gamma }_h^n}{\chi }_j^{n+1}\\ =\frac{1}{\tau }{\int }_{{\Gamma }^n}{U}^n{\chi }_j^n+{\int }_{{\Gamma }_h^{n+1}}^n{f}_1\left(U^{n+1},W^{n+1}\right){\chi }_j^{n+1}\\ \frac{1}{\tau }{\int }_{{\Gamma }_h^{n+1}}{W}^{n+1}{\chi }_j^{n+1}+d{\int }_{{\Gamma }_h^{n+1}}{\nabla }_{{\Gamma }_h^n}{W}^{n+1}\cdot {\nabla }_{{\Gamma }_h^n}{\chi }_j^{n+1}\\ =\frac{1}{\tau }{\int }_{{\Gamma }^n}{W}^n{\chi }_j^n+{\int }_{{\Gamma }_h^{n+1}}{f}_2\left(U^{n+1},W^{n+1}\right){\chi }_j^{n+1}\end{array}$$
将 $f$ 的表达式代入，再利用简单的线性化技巧，
$${\left(U^{n+1}\right)}^2\approx U^n{U}^{n+1}$$
最后得到一个线性化系统，
$$\{\begin{array}{c}\left(\frac{1}{\tau }+\gamma \right){\int }_{{\Gamma }_h^{n+1}}{U}^{n+1}{\chi }_j^{n+1}+{\int }_{{\Gamma }_h^{n+1}}{\nabla }_{{\Gamma }_h^n}{U}^{n+1}\cdot {\nabla }_{{\Gamma }_h^n}{\chi }_j^{n+1}\\ -\gamma {\int }_{{\Gamma }_h^{n+1}}{U}^n{W}^n{U}^{n+1}{\chi }_j^{n+1}=\frac{1}{\tau }{\int }_{{\Gamma }^n}{U}^n{\chi }_j^n+\gamma a{\int }_{{\Gamma }_h^n}{\chi }_j^n,\\ \begin{array}{c}\frac{1}{\tau }{\int }_{{\Gamma }_h^{n+1}}{W}^{n+1}{\chi }_j^{n+1}+d{\int }_{\Gamma }^{n+1}{\nabla }_{{\Gamma }_h^n}{W}^{n+1}\cdot {\nabla }_{{\Gamma }_h^n}{\chi }_j^{n+1}\\ +{\gamma }_{\Gamma }^{n+1}{\left(U^{n+1}\right)}^2{W}^{n+1}{\chi }_j^{n+1}=\frac{1}{\tau }{\int }_{{\Gamma }^n}{W}^n{\chi }_j^n+\gamma b{\int }_{{\Gamma }_h^n}{\chi }_j^n\end{array}\end{array}$$
代数系统的形式为，
$$\{\begin{array}{l}\left(\left(\frac{1}{\tau }+\gamma \right){\mathcal{M}}^{n+1}+{\mathcal{S}}^{n+1}-\gamma {\mathcal{M}}_1^{n+1}\right){\mathbf{U}}^{n+1}=\frac{1}{\tau }{\mathcal{M}}^n{\mathbf{U}}^n+{\mathbf{F}}_1^n\\ \left(\frac{1}{\tau }{\mathcal{M}}^{n+1}+d{\mathcal{S}}^{n+1}+\gamma {\mathcal{M}}_2^{n+1}\right){\mathbf{W}}^{n+1}=\frac{1}{\tau }{\mathcal{M}}^n{\mathbf{W}}^n+{\mathbf{F}}_2^n\end{array}$$
其中，
$$\begin{array}{rl}{\mathcal{M}}_{ij}^n& ={\int }_{{\Gamma }_h^n}{\chi }_i^n{\chi }_j^n,{\mathcal{M}}_{1ij}^n={\int }_{{\Gamma }_h^n}{U}^n{W}^n{\chi }_i^n{\chi }_j^n,{\mathcal{M}}_{2ij}^n={\int }_{{\Gamma }_h^n}{\left(U^{n+1}\right)}^2{\chi }_i^n{\chi }_j^n\\ {\mathcal{S}}_{ij}^n& ={\int }_{{\Gamma }_h^n}{\nabla }_{{\Gamma }_h^n}{\chi }_i^n\cdot {\nabla }_{{\Gamma }_h^n}{\chi }_j^n,{\mathbf{F}}_{1i}^n=\gamma a{\int }_{{\Gamma }_h^n}{\chi }_i^n,{\mathbf{F}}_{2i}^n=\gamma b{\int }_{{\Gamma }_h^n}{\chi }_i^n\end{array}$$

该方法需要对初始曲面进行三角剖分。例如，球体按以下方式进行三角剖分。使用均匀放置在球体上的 6 个点手动获得非常粗略的三角剖分。然后通过将每个三角形细分。然后将新顶点投影到球体表面。重复这种细化，直到获得可接受的最大网格尺寸。

肿瘤增长

考虑肿瘤生长模型来说明 ESFEM 在耦合曲面演化与曲面反应扩散的系统中的应用。

模型建立

该模型的曲面演化是化学驱动型的。这里 $u$ 表示生长促进因子，$w$ 表示生长抑制因子。此外，我们认为表面演化受的法向速度控制，表达如下，
$$V=\delta u-ϵH$$
其中 $\delta$ 是生长速率，$ϵ$ 可以被认为是模拟肿瘤表面张力的参数。我们也可以将 $ϵ$ 视为一个小的数学上的参数，它规范了进化规律。这很有用，因为很容易看出，即使在法线方向上表面速度恒定的简单情况下，光滑表面也可能在短时间内形成拐角。ESFEM 的一个优点是我们只需要三角形顶点的位置来进行计算，从而导致该方法具有吸引力的灵活性。

数值方法

数值方法分成两步。
第一步，把曲面的线性逼近写成基函数的线性组合，组合系数为节点的位置，然后，使用标准的有限元方法，根据当前逼近曲面顶点的位置，得到下一步逼近曲面的顶点的位置。
$$\left(\frac{1}{\tau }{\mathcal{M}}^n+ϵ{\mathcal{S}}^n\right){\mathbf{X}}^{n+1}=\frac{1}{\tau }\mathcal{M}{\mathbf{X}}^n+{\mathbf{G}}^n$$
$${\mathbf{G}}_i^n={\int }_{{\Gamma }^n}\delta {\left({\mathbf{U}}_h^n\right)}_1{\mathbf{v}}^n{\chi }_i$$
其中，${\mathbf{X}}^n$ 表示节点坐标。
第二步，根据得到的下一步的分片曲面的顶点位置，我们可以得到新的基函数，利用基函数计算刚度矩阵和质量矩阵，进而执行之前提到的离散格式：

数值结果

选定参数 $(u(\cdot ,0),w(\cdot ,0))=(1.0,0.9)$，及

d	$\gamma$	a	b	$\overline{(}t)$	$\delta$	$ϵ$	Timestep	Solver	Solver Tol.
10	30	0.1	0.9	5.	0.1	0.01	0.001	GMRes	10^(-9)

网格自适应和重新划分网格

原则上一段时间后，由于三角形的顶点根据速度规律演化，演化表面的三角剖分可能不适用于有限元模拟。例如，它们可以出现非常小的角度从而出现的某种意义上退化。另一方面，三角剖分可能无法解析浓度快速变化的子区域中的解的特征。由于夹断和自相交，表面的拓扑结构也可能发生变化。

上面提出的模拟中，我们要保证初始三角剖分足够细足够均匀足以进行计算。可以通过以自适应方式细化和粗化三角剖分来执行更有效的计算，以反映曲面偏微分方程的解的行为并避免三角剖分的退化。

---

1. Turing A (1952) The chemical basis of morphogenesis. Phil Trans R Soc Lond B 237:37–72 ↩︎

2. Murray JD (2002) Mathematical biology I and II, 3rd edn. Springer, Berlin ↩︎

3. Crampin EJ, Gaffney EA, Maini PK (1999) Reaction and diffusion on growing domains: scenarios for robust pattern formation. Bull Math Biol 61:1093–1120 ↩︎

4. Kondo S, Asai R (1995) A reaction-diffusion wave on the skin of the marine anglefish, Pomacanthus. Nature 376:765–768 ↩︎

5. Maini PK, Painter KJ, Chau HNP (1997) Spatial pattern formation in chemical and biological systems. Faraday Trans 93:3601–3610 ↩︎

6. Madzvamuse A, Maini PK, Wathen AJ (2003) A moving grid finite element method applied to a model biological pattern generator. J Comp Phys 190:478–500 ↩︎ ↩︎

7. Aragón JL, Barrio RA, Varea C (1999) Turing patterns on a sphere. Phys Rev E 60:4588–4592 ↩︎

8. Barrio RA, Maini PK, Padilla P, Plaza RG, Sánchez-Garduno F (2004) The effect of growth and curvature on pattern formation. J Dyn Diff Equ 4:1093–1121 ↩︎

9. Chaplain MAJ, Ganesh M, Graham IG (2001) Spatio-temporal pattern formation on spherical surfaces: numerical simulation and application to solid tumour growth. J Math Biol 42:387–423 ↩︎

10. Dziuk G, Elliott CM (2007) Finite elements on evolving surfaces. IMA J Num Anal 27:262–292 ↩︎

11. Dziuk G, Elliott CM (2007) Surface finite elements for parabolic equations. J Comp Math 25:430–439 ↩︎

12. Eilks C, Elliott CM (2008) Numerical simulation of dealloying by surface dissolution via the evolving surface finite element method. J Comp Phys 227:9727–9741 ↩︎

13. Barrett JW, Garcke H, Nürnberg (2008) Parametric approximation of Willmore flow and related geometric evolution equations. SIAM J Sci Comput 31:225–253 ↩︎

14. Barreira R (2009) Numerical solution of non-linear partial differential equations on triangulated surfaces, D.Phil Thesis, University of Sussex ↩︎

15. Elliott CM, Stinner B (2010) Modeling and computation of two phase geometric biomembranes using surface finite elements. J Comp Phys 229:6585–6612 ↩︎

16. Gierer A, Meinhardt H (1972) A theory of biological pattern formation. Kybernetik 12:30–39 ↩︎

17. Prigogine I, Lefever R (1968) Symmetry breaking instabilities in dissipative systems. II. J Chem Phys 48:1695–1700 ↩︎

18. Schnakenberg J (1979) Simple chemical reaction systems with limit cycle behaviour. J Theor Biol 81:389–400 ↩︎

19. Golub GH, Van Loan CF (1996) Matrix Computations. JHU Press, Baltimore ↩︎

20. Greer J, Bertozzi AL, Sapiro G (2006) Fourth order partial differential equations on general geometries. J Comput Phys 216:216–246 ↩︎

21. Calhoun DA, Helzel C (2009) A finite volume method for solving parabolic equations on logically cartesian curved surface meshes. SIAM J Sci Comput 6:4066–4099 ↩︎

22. Dziuk G, Elliott CM (2010) An Eulerian approach to transport and diffusion on evolving surfaces. Comput Vis Sci 13:17–28 ↩︎

23. Lefevre J, Mangin J-F (2010) A reaction-diffusion model of the human brain development. PLoS Comput Biol 6:e1000749 ↩︎

24. Elliott CM, Stinner B, Styles VM (2010) Numerical computation of advection and diffusion on evolv-ing diffuse interfaces. IMA J Num Anal. Advance Access published on 11 May 2010. ↩︎

25. Deckelnick KP, Dziuk G, Elliott CM (2005) Computation of Geometric PDEs and Mean Curvature Flow. Acta Numerica 14:139–232 ↩︎