相场理论基础

Foundation of Phase Field Modeling
- 基本原理
- - 相场参量（序参量）
  - 扩散界面模型
  - 自由能泛函
  - 相场法的优势
- 体自由能密度
- - 1. 二元混合物中的相分离
  - 2. 相变的朗道自由能泛函
  - 3. 凝固
- 平衡
- 相场演化动力学
- - Allen-Cahn方程
  - Cahn-Hilliard 方程
  - Landau-Lifshitz-Gilbert 方程
  - 随机涨落
  - 界面特性
- 相场模型的计算
- - 有限差分算法
  - - 应用
    - - 一维瞬态热传导
      - 二元合金的斯皮诺达分解（Spinodal Decomposition）
      - 晶粒生长
      - 固相烧结
      - 枝晶凝固
      - 多细胞系统的结构演化
  - 傅里叶谱法
  - - 应用
    - - 多元合金的析出行为
      - 二元合金相分离行为
      - 合金中的晶格缺陷
  - 有限元方法

Foundation of Phase Field Modeling

你好！这是你第一次使用 Markdown编辑器 所展示的欢迎页。如果你想学习如何使用Markdown编辑器, 可以仔细阅读这篇文章，了解一下Markdown的基本语法知识。

基本原理

对于特定物理体系的状态或过程，平衡状态可以视为有约束（或无约束）条件下的全局能量最小化的结果。 基于广义力或通量的运动学模型是用来表示微结构演化的经典模型。

相场参量（序参量）

在相场方法中，根据相场参量的取值而区分不同的相。 对于特定的物理问题，需先选择一个变量作为相场参量，它可以是模型中的具体物理量，也可以是一组变量的函数。通过定义好的相场参量，可以反应体系的物理性质，比如结晶度（degree of crystallinity）、原子排列（atomic ordering）、晶体对称性（crystal symmetry）、状态（state）、磁化强度（magnetization）、电极化量（electric polarization）、浓度（concentration）、组成（composition）、损伤程度（damage level）等。

因为相场法研究的更多是微观/介观尺度，而非原子/电子尺度，所以相场参量是相的一种特定热力学性质的平均值。相场参量也被称为序参量。

扩散界面模型

在传统相变模型中，常因为追踪不同相之间界面的变化而导致大量繁琐的计算。相场法假定相之间的界面是扩散性（diffusive）的，介于相邻的相之间的扩散界面由用 ϕ\phiϕ 表示的相场参量的连续插值描述。

在扩散界面模型下，微观结构类似于相场变量的一个分布函数 ϕ(x)\phi(x)ϕ(x)，其中 xxx是位置向量。这种模型提供了用公式表示体积内能量密度的可能。在传统的尖锐界面下，界面对能量的贡献与界面能密度和界面拓扑结构有关，而界面形态（interface morphology）轨迹，特别是三维情况下，可能会十分复杂。而在扩散界面模型下，界面能够通过相场参量的梯度进行识别。因此界面形态和界面能量都能够很自然地基于梯度项定义。

自由能泛函

每个相的自由能密度都依赖于相的热力学状态，因此每个相的体自由能密度（bulk free energy density）都能够由相场参量ϕ\phiϕ 用公式描述。界面上的自由能密度是不同于相内部的，所以界面位置能够用相场参量的梯度来描述。系统 A 的自由能 FFF 能够由以下相场自由能泛函给定：
F(ϕ,∇ϕ)=∫A[f(ϕ)+12κ∣∇ϕ∣2]dVF\left( {\phi ,\nabla \phi } \right) = \int_A {\left[ {f\left( \phi \right) + {1 \over 2}\kappa {{\left| {\nabla \phi } \right|}^2}} \right]dV} F(ϕ,∇ϕ)=∫A[f(ϕ)+21κ∣∇ϕ∣2]dV其中，f(ϕ)f\left( \phi \right)f(ϕ) 是体自由能密度， 12κ∣∇ϕ∣2{{1 \over 2}\kappa {{\left| {\nabla \phi } \right|}^2}}21κ∣∇ϕ∣2是包含参量κ\kappaκ的梯度能密度。

体积自由能密度f(ϕ)f\left( \phi \right)f(ϕ) 依赖于局域的相场参量。例如，一个多相系统的体积自由能密度图像通常呈多势阱的形状。梯度项12κ∣∇ϕ∣2{{1 \over 2}\kappa {{\left| {\nabla \phi } \right|}^2}}21κ∣∇ϕ∣2是非局域的，并且会在界面处引起能量损耗。需要注意的是，仅当κ→0\kappa \to 0κ→0 时，尖锐界面才是可能的。

另外，在某些问题中，还需要考虑梯度的高阶项，例如∣∇ϕ∣4|\nabla\phi|^4∣∇ϕ∣4。由于界面能的空间对称性，更高阶项应该是偶数阶的。

相场法的优势

在相场理论中，系统的状态和动力学由能量最小原理决定，系统处于平衡状态时，序参量ϕeq(x){\phi _{eq}}\left( x \right)ϕeq(x)应当使得自由能 FFF 取得最小值，即自由能的一阶变分导数为零：δF(ϕ,∇ϕ)δϕ=0{{\delta F\left( {\phi ,\nabla \phi } \right)} \over {\delta \phi }} = 0δϕδF(ϕ,∇ϕ)=0

对于不受约束的热力学量，可以假定微结构演化过程是以整体能量或势最快降低的方向进行的。在相场理论中，这种关于微结构变的整体能量变化即由这个一阶变分导数给定。因此相场参量下的动力学模型能够直接由该导数用公式明确表示。注意到该一阶变分导数能够由下述变分得到：

δF=δ∫A[f(ϕ)+12κ∣∇ϕ∣2]dV=∫A[δf(ϕ)+12κδ∣∇ϕ∣2]dV=∫A[∂f(ϕ)∂ϕδϕ+κ∂iϕ∂iδϕ]dV=∫A{∂f(ϕ)∂ϕδϕ+κ[∂i(∂iϕδϕ)−∂i∂iϕδϕ]}dV=∫A[∂f(ϕ)∂ϕ−κ∂iiϕ]δϕdV+∫Aκ∂i(∂iϕδϕ)dV=∫A[∂f(ϕ)∂ϕ−κ∂iiϕ]δϕdV+∫∂Aκ∂iϕniδϕdS=∫A[∂f(ϕ)∂ϕ−κΔϕ]δϕdV+∫∂Aκ∇ϕ⋅nδϕdS\begin{aligned} \delta F &= \delta \int_A {\left[ {f\left( \phi \right) + {1 \over 2}\kappa {{\left| {\nabla \phi } \right|}^2}} \right]{\rm{d}}V} \\ & = \int_A {\left[ {\delta f\left( \phi \right) + {1 \over 2}\kappa \delta {{\left| {\nabla \phi } \right|}^2}} \right]{\rm{d}}V} \\ & = \int_A {\left[ {{{\partial f\left( \phi \right)} \over {\partial \phi }}\delta \phi + \kappa {\partial _i}\phi {\partial _i}\delta \phi } \right]{\rm{d}}V} \\ & = \int_A {\left\{ {{{\partial f\left( \phi \right)} \over {\partial \phi }}\delta \phi + \kappa \left[ {{\partial _i}\left( {{\partial _i}\phi \delta \phi } \right) - {\partial _i}{\partial _i}\phi \delta \phi } \right]} \right\}{\rm{d}}V} \\ & = \int_A {\left[ {{{\partial f\left( \phi \right)} \over {\partial \phi }} - \kappa {\partial _{ii}}\phi } \right]\delta \phi {\rm{d}}V} + \int_A {\kappa {\partial _i}\left( {{\partial _i}\phi \delta \phi } \right){\rm{d}}V} \\ & = \int_A {\left[ {{{\partial f\left( \phi \right)} \over {\partial \phi }} - \kappa {\partial _{ii}}\phi } \right]\delta \phi {\rm{d}}V} + \int_{\partial A} {\kappa {\partial _i}\phi {n_i}\delta \phi {\rm{d}}S} \\ & = \int_A {\left[ {{{\partial f\left( \phi \right)} \over {\partial \phi }} - \kappa \Delta \phi } \right]\delta \phi {\rm{d}}V} + \int_{\partial A} {\kappa \nabla \phi \cdot {\bm{n}}\delta \phi {\rm{d}}S} \end{aligned}δF=δ∫A[f(ϕ)+21κ∣∇ϕ∣2]dV=∫A[δf(ϕ)+21κδ∣∇ϕ∣2]dV=∫A[∂ϕ∂f(ϕ)δϕ+κ∂iϕ∂iδϕ]dV=∫A{∂ϕ∂f(ϕ)δϕ+κ[∂i(∂iϕδϕ)−∂i∂iϕδϕ]}dV=∫A[∂ϕ∂f(ϕ)−κ∂iiϕ]δϕdV+∫Aκ∂i(∂iϕδϕ)dV=∫A[∂ϕ∂f(ϕ)−κ∂iiϕ]δϕdV+∫∂Aκ∂iϕniδϕdS=∫A[∂ϕ∂f(ϕ)−κΔϕ]δϕdV+∫∂Aκ∇ϕ⋅nδϕdS

给定κ∇ϕ⋅n=0\kappa \nabla \phi \cdot {\bm{n}} = 0κ∇ϕ⋅n=0，则由上式可得：
δFδϕ=∂f(ϕ)∂ϕ−κΔϕ{{\delta F} \over {\delta \phi }} = {{\partial f\left( \phi \right)} \over {\partial \phi }} - \kappa \Delta \phi δϕδF=∂ϕ∂f(ϕ)−κΔϕ 其中，Δ\DeltaΔ是拉普拉斯算子。

体自由能密度

体自由能密度f(ϕ)f\left( \phi \right)f(ϕ)依赖于相场参量的意义。下面举例一些常用的体积自由能。

1. 二元混合物中的相分离

在这个问题中，我们将其中一种组分的密度作为相场参量ϕ\phiϕ ，基于二元混合物的平均场自由能理论，体积自由能采用以下形式：
f(ϕ)=εν2ϕ(1−ϕ)+kBT[ϕln⁡ϕ+(1−ϕ)ln⁡(1−ϕ)]+bϕ−νεAA2f\left( \phi \right) = {{\varepsilon \nu } \over 2}\phi \left( {1 - \phi } \right) + {k_B}T\left[ {\phi \ln \phi + \left( {1 - \phi } \right)\ln \left( {1 - \phi } \right)} \right] + b\phi - {{\nu {\varepsilon _{AA}}} \over 2}f(ϕ)=2ενϕ(1−ϕ)+kBT[ϕlnϕ+(1−ϕ)ln(1−ϕ)]+bϕ−2νεAA 其中， kB{k_B}kB 是玻尔兹曼常量，TTT是温度， ε,εAA,ν{\varepsilon}, {\varepsilon_{AA}}, {\nu}ε,εAA,ν 是材料参数。

如上图左可见，随着温度 TTT 升高并超过临界温度 TC=εν4kB{T_C} = {{\varepsilon \nu } \over {4{k_B}}}TC=4kBεν ，体积自由能也从双势阱变为了单势阱。低于临界温度时，存在两种平衡状态，分别对应两个势阱。当 b=0b=0b=0时，这两个平衡状态实际上是局部最小值状态。而对于b≠0b\ne 0b=0，这两种平衡状态能够被余切定理控制，如上图右所示。当T→TCT \to {T_C}T→TC时，这两个状态连续地汇聚为一个。

2. 相变的朗道自由能泛函

在朗道理论中，假定自由能在相 ϕ=0\phi=0ϕ=0 处的级数展开：
f(ϕ)=f(T,ϕ=0)+∑n=2Man(T)nϕnf\left( \phi \right) = f\left( {T,\phi = 0} \right) + \sum\limits_{n = 2}^M {{{{a_n}\left( T \right)} \over n}} {\phi ^n}f(ϕ)=f(T,ϕ=0)+n=2∑Mnan(T)ϕn
对于对称相变，可以使用
f(ϕ)=a(T)+a2(T)2ϕ2+a4(T)4ϕ4+O(ϕ)f\left( \phi \right) = a\left( T \right) + {{{a_2}\left( T \right)} \over 2}{\phi ^2} + {{{a_4}\left( T \right)} \over 4}{\phi ^4} + O\left( \phi \right)f(ϕ)=a(T)+2a2(T)ϕ2+4a4(T)ϕ4+O(ϕ)
其中，系数an(T){a_n}\left( T \right)an(T)也能够被级数展开为线性项。例如
a2(T)=a20(T−TC)a4(T)=a40+b40(T−TC)\begin{aligned} & {a_2}\left( T \right) = a_2^0\left( {T - {T_C}} \right) \\ & {a_4}\left( T \right) = a_4^0 + b_4^0\left( {T - {T_C}} \right) \\ \end{aligned} a2(T)=a20(T−TC)a4(T)=a40+b40(T−TC) 高于临界温度时，a2(T),a4(T)>0{a_2}\left( T \right),{a_4}\left( T \right) > 0a2(T),a4(T)>0。

当温度高于临界温度时，只存在一个稳态 ϕ=0\phi = 0ϕ=0，而当温度低于临界温度时，存在两个对称的稳态 ϕ≈±a20a40(TC−T)\phi \approx \pm \sqrt {{{a_2^0} \over {a_4^0}}\left( {{T_C} - T} \right)} ϕ≈±a40a20(TC−T)
当T→TCT \to {T_C}T→TC 时，两个稳态向 ϕ=0\phi = 0ϕ=0 汇聚。上述代表了一个典型的二阶相变。

通过增加奇数阶项，可以获得非对称相图的相变体积自由能。此外，通过调节an(T){a_n}\left( T \right)an(T) 的符号也能获得一阶相变的多项式。

3. 凝固

对纯材料单晶的凝固过程，体自由能密度可以设为以下形式：
f(ϕ,T)=f0−SL(T−Tm)+Hϕ2(1−ϕ)2+[SL−LTm(3−2ϕ)ϕ2]ϕ2(T−Tm)f\left( {\phi ,T} \right) = {f_0} - {S_L}\left( {T - {T_m}} \right) + H{\phi ^2}{\left( {1 - \phi } \right)^2} + \left[ {{S_L} - {L \over {{T_m}}}\left( {3 - 2\phi } \right){\phi ^2}} \right]{\phi ^2}\left( {T - {T_m}} \right)f(ϕ,T)=f0−SL(T−Tm)+Hϕ2(1−ϕ)2+[SL−TmL(3−2ϕ)ϕ2]ϕ2(T−Tm) 其中，SLS_LSL 是液相的体积熵密度； TmT_mTm是熔点；HHH 是双势阱高度； LLL 是熔化潜热。
右式第三项形成了一个双势阱，而最后一项破坏了关于熔点的对称性，如下图所示：

实际上，使用多项式来形成双势阱体自由能泛函在相场模型中是非常常见的。

平衡

平衡微结构的系统全局能量或势最小化过程依赖于约束条件。一般来说，在平衡状态下，相场参量的驱动力或通量会变为。而在相场模型下，驱动力或通量可以表示为自由能泛函关于相场参量的一阶变分导数，所以平衡状态微结构应该满足以下关系式：
δFδϕ=∂f(ϕ)∂ϕ−κΔϕ=0{{\delta F} \over {\delta \phi }} = {{\partial f\left( \phi \right)} \over {\partial \phi }} - \kappa \Delta \phi = 0δϕδF=∂ϕ∂f(ϕ)−κΔϕ=0
以一个简单的体积自由能的多项式双势阱函数f(ϕ)=4Hϕ2(1−ϕ)2f\left( \phi \right) = 4H{\phi ^2}{\left( {1 - \phi } \right)^2}f(ϕ)=4Hϕ2(1−ϕ)2 为例。它有两个稳定的相ϕ=0\phi = 0ϕ=0，ϕ=1\phi = 1ϕ=1，并且能量势垒为HHH。一个边界条件为ϕ∣x=−∞=0{\left. \phi \right|_{x = - \infty }} = 0ϕ∣x=−∞=0，ϕ∣x=+∞=1{\left. \phi \right|_{x = + \infty }} = 1ϕ∣x=+∞=1的一维问题，对上一个等式的解为：
ϕeq(x)=12[1+tanh⁡(x2l)]{\phi _{eq}}\left( x \right) = {1 \over 2}\left[ {1 + \tanh \left( {{x \over {2l}}} \right)} \right]ϕeq(x)=21[1+tanh(2lx)]
式中，l=κHl = \sqrt {{\kappa \over H}} l=Hκ

体自由能和序参量的图像如下
该模型中包括两个序参量：HHH代表稳态之间的势垒，$\kappa 为梯度能项的系数。另外，模型中包含两个界面相关的物理参量：为梯度能项的系数。另外，模型中包含两个界面相关的物理参量：为梯度能项的系数。另外，模型中包含两个界面相关的物理参量：l为界面厚度，为界面厚度，为界面厚度，G$为界面能密度。

对于给定的双势阱函数f(ϕ)=4Hϕ2(1−ϕ)2f\left( \phi \right) = 4H{\phi ^2}{\left( {1 - \phi } \right)^2}f(ϕ)=4Hϕ2(1−ϕ)2，有l=κHl=\sqrt {\frac{\kappa}{H}}l=Hκ。

另外，因为两个稳态处的体自由能密度为000，所以体自由能密度f(ϕeq)+12κ∣∇ϕeq(x)∣2f\left( {{\phi _{eq}}} \right) + {1 \over 2}\kappa {\left| {\nabla {\phi _{eq}}\left( x \right)} \right|^2}f(ϕeq)+21κ∣∇ϕeq(x)∣2 在整个域上的积分即为界面上的界面能：G=∫−∞+∞[f(ϕeq(x))+12κ∣∇ϕeq(x)∣2]dx=13κHG = \int_{ - \infty }^{ + \infty } {\left[ {f\left( {{\phi _{eq}}\left( x \right)} \right) + {1 \over 2}\kappa {{\left| {\nabla {\phi _{eq}}\left( x \right)} \right|}^2}} \right]} {\rm{d}}x = {1 \over 3}\sqrt {\kappa H} G=∫−∞+∞[f(ϕeq(x))+21κ∣∇ϕeq(x)∣2]dx=31κH
对于一般形式的体积自由能，界面参量依赖于模型参量：
lκH，lκHl~\sqrt {\frac{\kappa}{H}}，l~\sqrt {\kappa H}l Hκ，l κH

用界面参量GGG，lll 来表达的相场模型有更加明晰的物理解释，比如裂纹扩展的相场模型经常用如下形式表示：F=∫AGl(g(ϕ)+Gl∣∇ϕeq(x)∣2)dVF = \int_A {{G \over l}\left( {g\left( \phi \right) + Gl{{\left| {\nabla {\phi _{eq}}\left( x \right)} \right|}^2}} \right)} {\rm{d}}VF=∫AlG(g(ϕ)+Gl∣∇ϕeq(x)∣2)dV

如果需要同时考虑超过一种物理特性，可以直接在体积自由能中增加相应的能量项，如果需要的话也可以添加耦合项，所以，相场模型能够被扩展应用到解决多物理特性问题上。

相场演化动力学

微观结构的演化能够通过定义了与时间相关的相场变化率的运动学模型进行模拟。根据序参量所受到的约束不同，对应的运动学模型也不同。

Allen-Cahn方程

根据 Allen-Cahn 模型（也称 Langevin-type 模型），对于无约束的序参量，微结构ϕ(x){\phi}\left( x \right)ϕ(x) 以最陡的泛函梯度方向使全局自由能最小化（或全局系统熵最大化），即
∂ϕ∂t=−M∂F∂ϕ=−M(δf(ϕ)δϕ−κΔϕ){{\partial \phi } \over {\partial t}} = - M{{\partial F} \over {\partial \phi }} = - M\left( {{{\delta f\left( \phi \right)} \over {\delta \phi }} - \kappa \Delta \phi } \right)∂t∂ϕ=−M∂ϕ∂F=−M(δϕδf(ϕ)−κΔϕ) 其中，MMM 是迁移率参数。上式是一个二阶偏微分方程。

值得注意的是，该模型只适用于纯弛豫过程的序参量。例如，非保守序参量的磁化强度、损伤变量和凝固。

Cahn-Hilliard 方程

如果序参量遵循以下形式的通量守恒法则：
∂ϕ∂t=−∇⋅j，j=−D∇μ{{\partial \phi } \over {\partial t}} = - \bm{\nabla} \cdot \bm{j}，\bm{j} = - D\nabla \mu ∂t∂ϕ=−∇⋅j，j=−D∇μ
则称序参量是保守的。

其中，j\bm{j}j是通量向量，它依赖于化学势 μ\muμ 的梯度。在这种情况下，假定化学势采取以下一阶变分导数形式：
μ=δFδϕ\mu = {{\delta F} \over {\delta \phi }}μ=δϕδF
因此
∂ϕ∂t=∇⋅M∇∂F∂ϕ{{\partial \phi } \over {\partial t}} = \nabla \cdot M\nabla {{\partial F} \over {\partial \phi }}∂t∂ϕ=∇⋅M∇∂ϕ∂F 这是一个四阶微分方程。

值得说明的是，浓度是典型的保守序参量。

Landau-Lifshitz-Gilbert 方程

在一些弛豫过程中，相场向量MMM 的模恒定。特别地，对于温度远低于居里点TCT_CTC 的磁化过程，磁化强度几乎不变。因此该磁化仅仅在微结构演化中旋转。对应的演化方程由 L-L-G 方程给定：
∂M∂t=−γ1M×δFδM−γ2M×(M×δFδM){{\partial M} \over {\partial t}} = - {\gamma _1}M \times {{\delta F} \over {\delta M}} - {\gamma _2}M \times \left( {M \times {{\delta F} \over {\delta M}}} \right)∂t∂M=−γ1M×δMδF−γ2M×(M×δMδF)
叉乘保证了磁化向量的驱动力与该向量保持垂直。

随机涨落

如果演化动力学方程中包含了热扰动或化学干扰，则需要考虑序参量的涨落，此时可以引入一个随机噪声项。例如，对于 A-C 型方程：
∂ϕ∂t=−M∂F∂ϕ+η(x,t){{\partial \phi } \over {\partial t}} = - M{{\partial F} \over {\partial \phi }} + \eta \left( {x,t} \right)∂t∂ϕ=−M∂ϕ∂F+η(x,t)
根据扰动产生的原因，噪声项应该遵循特定的时空随机分布。

界面特性

界面的各向异性和迁移率对演化的影响十分重要。在特定的情况下，梯度项系数和迁移率参数可以由晶格取向，甚至由相场参量梯度方向表示的界面取向决定。例如：M→M(θ)，κ→κ(θ)，θ=arctan⁡(∇ϕ)M \to M\left( \theta \right)，\kappa \to \kappa \left( \theta \right)，\theta {\rm{ = }}\arctan \left( {\nabla \phi } \right)M→M(θ)，κ→κ(θ)，θ=arctan(∇ϕ)

相场模型的计算

有限差分算法

有限差分算法为偏微分方程的数值解提供了比任何其他方法都更直接的方法。有限差分算法编码简单、计算经济、易于并行，适用于分布式计算环境。然而，它们在精度和施加复杂边界条件方面也有缺点。

有限差分算法的核心思想是将每个导数替换为一个差商。对于一维网格上连续且相对光滑的函数uuu，网格点间距为hhh，有限差分算法采用如下公式：

后向差分：(Δu)i−=ui−ui−1h\left( {\Delta u} \right)_i^ - = {{{u_i} - {u_{i - 1}}} \over h}(Δu)i−=hui−ui−1
前向差分：(Δu)i+=ui+1−uih\left( {\Delta u} \right)_i^ + = {{{u_{i + 1}} - {u_i}} \over h}(Δu)i+=hui+1−ui
中心差分：(Δu)i±=ui+1−ui−12h\left( {\Delta u} \right)_i^ \pm = {{{u_{i + 1}} - {u_{i - 1}}} \over {2h}}(Δu)i±=2hui+1−ui−1
二阶中心差分：(Δ2u)i=ui+1−2ui+ui−1h2\left( {\Delta^2 u} \right)_i = {{{u_{i + 1}} -2{u_i}+ {u_{i - 1}}} \over {h^2}}(Δ2u)i=h2ui+1−2ui+ui−1

经过泰勒展开后，各差分的误差项可表示为如下形式：
(Δu)i−−u′(xi)=O(h)(Δu)i+−u′(xi)=O(h)(Δu)i±−u′(xi)=O(h2)(Δ2u)i−u′′(xi)=O(h2)\begin{aligned} & \left( {\Delta u} \right)_i^ - - u'\left( {{x_i}} \right) = O\left( h \right) \\ & \left( {\Delta u} \right)_i^ + - u'\left( {{x_i}} \right) = O\left( h \right) \\ & \left( {\Delta u} \right)_i^ \pm - u'\left( {{x_i}} \right) = O\left( {h^2} \right) \\ & \left( {{\Delta^2} u} \right)_i - u''\left( {{x_i}} \right) = O\left( {h^2} \right) \end{aligned} (Δu)i−−u′(xi)=O(h)(Δu)i+−u′(xi)=O(h)(Δu)i±−u′(xi)=O(h2)(Δ2u)i−u′′(xi)=O(h2)

这些差分公式，特别是近似拉普拉斯算子的中心二阶差分，将用于相场模型的数值解。
对于二维网格，在网格内点(xi,yj)\left( {{x_i},{y_j}} \right)(xi,yj)处，拉普拉斯算子为：
(∇2u)i,j=(Δxx2u)i,j+(Δyy2u)i,j=ui+1,j−2ui,j+ui−1,jhx2+ui,j+1−2ui,j+ui,j−1hy2\begin{aligned} {\left( {{\nabla ^2}u} \right)_{i,j}} & = {\left( {\Delta _{xx}^2u} \right)_{i,j}} + {\left( {\Delta _{yy}^2u} \right)_{i,j}} \\ & = {{{u_{i + 1,j}} - 2{u_{i,j}} + {u_{i - 1,j}}} \over {h_x^2}} + {{{u_{i,j + 1}} - 2{u_{i,j}} + {u_{i,j - 1}}} \over {h_y^2}} \end{aligned} (∇2u)i,j=(Δxx2u)i,j+(Δyy2u)i,j=hx2ui+1,j−2ui,j+ui−1,j+hy2ui,j+1−2ui,j+ui,j−1

这种近似称为五点法，只涉及五个格点，这可以从方程和图1中推断出来。

九点法可以离散化误差降低到O(h4)O\left( {{h^4}} \right)O(h4)，因此可以使用更大的hhh，而误差小于五点法。然而，这种情况会导致计算成本的增加。因此，在相场模型求解时，通常选用五点法来进行空间离散化。

如果在空间中 xxx 和 yyy 方向上的格点间距彼此相等，h=hx=hyh=h_x=h_yh=hx=hy，则前式可以继续化简为：
(∇2u)i,j=ui+1,j+ui−1,j+ui,j+1+ui,j−1−4ui,jh2\left( {\nabla^2 u} \right)_{i,j} = {{{u_{i + 1,j}} + {u_{i - 1,j}} + {u_{i,j+1}} + {u_{i ,j-1}}-4{u_{i,j}}} \over {h^2}}(∇2u)i,j=h2ui+1,j+ui−1,j+ui,j+1+ui,j−1−4ui,j

对于图1所示等间距的9个格点，在不考虑任何边界条件的情况下，可以使用矩阵向量积的形式将上式重写为：

(∇2u)i,j=1h2(−4101000001−4101000001−4001000100−4101000101−4101000101−4001000100−4100000101−4100000101−4)(u1,1u2,1u3,1u1,2u2,2u3,2u1,3u2,3u3,3)(1){\left( {{\nabla ^2}u} \right)_{i,j}} = {1 \over {{h^2}}} \left( \begin{matrix} { - 4} & 1 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 1 & { - 4} & 1 & 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 1 & { - 4} & 0 & 0 & 1 & 0 & 0 & 0 \\ 1 & 0 & 0 & { - 4} & 1 & 0 & 1 & 0 & 0 \\ 0 & 1 & 0 & 1 & { - 4} & 1 & 0 & 1 & 0 \\ 0 & 0 & 1 & 0 & 1 & { - 4} & 0 & 0 & 1 \\ 0 & 0 & 0 & 1 & 0 & 0 & { - 4} & 1 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 1 & { - 4} & 1 \\ 0 & 0 & 0 & 0 & 0 & 1 & 0 & 1 & { - 4} \end{matrix} \right) \left( \begin{matrix} {{u_{1,1}}} \\ {{u_{2,1}}} \\ {{u_{3,1}}} \\ {{u_{1,2}}} \\ {{u_{2,2}}} \\ {{u_{3,2}}} \\ {{u_{1,3}}} \\ {{u_{2,3}}} \\ {{u_{3,3}}} \end{matrix} \right) \tag{1}(∇2u)i,j=h21−4101000001−4101000001−4001000100−4101000101−4101000101−4001000100−4100000101−4100000101−4u1,1u2,1u3,1u1,2u2,2u3,2u1,3u2,3u3,3(1)

通常，该公式可以重写为：
(∇2u)i,j=1h2MUui,j{\left( {{\nabla ^2}u} \right)_{i,j}} = {1 \over {{h^2}}} M U_{u_{i,j}} (∇2u)i,j=h21MUui,j
式中，
M=(AIOOOIAI⋯OOOIAOO⋮⋱⋮OOO⋯IA)M = \left( \begin{matrix} A & I & O & & O & O \\ I & A & I & \cdots & O & O \\ O & I & A & & O & O \\ & \vdots & & \ddots & \vdots \\ O & O & O & \cdots & I & A \\ \end{matrix} \right)M=AIOOIAI⋮OOIAO⋯⋱⋯OOO⋮IOOOA 其中，III是一个3×33 \times 33×3 的单位阵，OOO是一个3×33 \times 33×3 的零矩阵，而A为：
A=(−4101−4101−4)A = \left( \begin{matrix} -4 & 1 & 0 \\ 1 & -4 &1 \\ 0& 1 & -4 \end{matrix} \right)A=−4101−4101−4
值得注意的是，MMM具有一些独特的性质：

MMM是一个稀疏矩阵，每行最多只有5个非零元素；
它是一个块三对角矩阵（block tri-diagonal），有三对角线和对角线块；
对称；
它是一个对角占优矩阵（diagonally dominant matrix）：矩阵中每个主对角元素的模都大于与它同行的其他元素的模的总和；
MMM的对角元素为负，而其他元素为正；
它是一个负定矩阵：其特征值均为负。

应用

一维瞬态热传导

ρcρ∂T∂t=∂∂x(λ∂T∂x)\rho {c_\rho }{{\partial T} \over {\partial t}} = {\partial \over {\partial x}}\left( {\lambda {{\partial T} \over {\partial x}}} \right)ρcρ∂t∂T=∂x∂(λ∂x∂T) Tin+1=Tin+μΔt(Ti+1n−2Tin+Ti−1nΔx2)T_i^{n + 1} = T_i^n + \mu \Delta t\left( {{{T_{i + 1}^n - 2T_i^n + T_{i - 1}^n} \over {\Delta {x^2}}}} \right)Tin+1=Tin+μΔt(Δx2Ti+1n−2Tin+Ti−1n)

二元合金的斯皮诺达分解（Spinodal Decomposition）

使用有限差分算法在相场模型中求解保守Cahn-Hilliard方程。

晶粒生长

使用有限差分算法在相场模型中求解多分量非守恒Allen-Cahn方程。

固相烧结

使用有限差分法求解守恒Cahn-Hillard方程和多个非守恒Allen-Cahn方程的耦合解。
该模型基于两类场变量对微观结构的表示。一个是守恒密度场ρ\rhoρ，它在固相处为1，在孔隙处为0，在固-孔界面处迅速变化。另一个是非守恒阶参数ηi\eta_{i}ηi，用于区分微观结构中的不同颗粒。

枝晶凝固

处理界面能和晶界迁移率的各向异性。

多细胞系统的结构演化

通过修改自由能泛函，使相场模型中用非守恒序参数描述的每个单元在保留其原始体积的同时，能够发生大变形和形状变化。公式仍遵循总能量最小化的原则，但是当细胞聚集在一起时，它们会发生集体变形，而不是收缩、生长或相互吸收。

傅里叶谱法

应用

多元合金的析出行为

在固态微观结构演变过程中，多元合金系统会同时发生相变和相分离。相变采用非守恒的Allen - Cahn公式，相分离过程中合金成分守恒采用Cahn - hilliard公式。

二元合金相分离行为

在微观组织演化的过程中，演化相的力学性能和施加的应力场的不匹配将导致新的动态力学平衡出现。

合金中的晶格缺陷

在考虑弹性不均匀性的基础上，还将考虑晶格缺陷：位错和晶界。

有限元方法

相场理论基础-Foundation of Phase Field Modeling相关推荐

matlab 相场法,相场法数值模拟及答案.ppt
相场法数值模拟phase-field modeling 内容相场法数值模拟介绍 (Introduction) 相场变量(Phase-field variables) 热力学势函数(th ...
Godot 4.0中的基于有向距离场SDF（Signed Distance Field）的实时全局光照技术
原文地址 Godot官网 <Godot 4.0 gets SDF based real-time global illumination> 作者 Juan Linietsky 6月28日 ...
matlab 相场法,晶体相场法模拟.pdf
晶体相场法模拟物理学报 ActaPhys．Sin． Vo1．60,No．8 (2011) 088104 多晶凝固及后续调幅分解过程的晶体相场法模拟张琪王锦程张亚丛杨根仓 (西北 ...
COMSOL——相场模拟
相场模拟基础相场方程经典Kobayashi模型模型概要 COMSOL实现 COMSOL中追踪移动界面的方法有变形几何DG.ALE.水平集.相场.不过此篇博客要说的不是其自带的用于CFD中的相场接 ...
COMSOL两相流（相场法）
COMSOL是一款功能非常强大的多物理场仿真软件,拥有可视化和后处理工具,可方便有效的展示您的仿真效果,具有灵活.高效.友好的特点,软件能够为用户提供各类工程问题的解决方案,从而帮助广大研究和工程人员 ...
【物理应用】基于相场法模拟金属镍的晶粒的长大过程附matlab完整代码
1 内容介绍金属镍是支撑信息技术产业和光伏产业发展的基础材料,它的性能直接影响和决定着相关高技术产业的发展.随着直拉硅单晶的直径不断增大以及集成电路特性线宽持续减小,金属镍的缺陷控制及其质量保证尤为 ...
定向凝固各向异性枝晶生长相场模拟
定向凝固各向异性枝晶生长相场模拟 1.依据Kobayashi的经典模型,实现定向凝固各向异性枝晶生长(可定量修改相关参数) 2.matlab手写代码,利用快速求解方法求解方程,代码注释详细 3.可利 ...
comsol相场枝晶模拟，初学者必不可少的模型案例
comsol相场枝晶模拟,初学者必不可少的模型案例 ID:6959657002308946[雨过@天晴]
RTS核心技术:流场寻路详解(Flow Field Pathfinding)
RTS里面经常会有很多角色,群体一起寻路到目的地附近,这种寻路是如何实现的,今天给大家详细的讲解基于流场寻路的算法.在本教程中,我将解释向量场寻路及其相对于Dijkstra等传统寻路算法的优势.对Di ...

相场理论基础-Foundation of Phase Field Modeling