1. 引言
不可混溶两相流(Immiscible Two-phase Flow)是指由固、液、气等相态中任意两相互不混溶的物质所组成的流动体系。由于相界面不断变形和相变的发生,两相流问题在模型建立和数值求解方面面临巨大挑战。扩散界面模型通过在宏观动量方程中引入相场耦合项,无需显式追踪界面,已广泛应用于金属凝固、枝晶长大、油藏热采中固–液、油–蒸汽等界面问题。
关于可压缩与不可压缩Navier-Stokes/Allen-Cahn (简称NSAC)系统的数学理论,已有多项重要成果。Feireisl等[1]在三维可压缩情形下证明了全局弱解的存在性,Chen等[2]研究了稳态问题并给出了弱解存在;Zhao [3]在小扰动条件下获得了全局强解及其衰减速率;Yan等[4]针对一维变粘度可压缩模型建立了全局强解理论;Chen等[5]构造了可压缩三维系统的全局重正化弱解,并在后续工作中(Chen等[6])给出了爆破准则及强/经典解的全局存在唯一性;Chen等[7]进一步证明了小扰动下解的指数型衰减。与此同时,Frigeri与Grasselli [8]针对非局域Cahn-Hilliard-Navier-Stokes系统建立了弱解存在性,而Ganesan等[9]则提出了用于多孔介质的扩散界面两相流模型。
在不可压缩扩散界面模型方面,Abels等[10]以及Yang等[11]分别研究了含表面活性剂与曲率相关移动率的两相不可压缩流并证明了弱解存在性;Abels等[12]利用匹配渐近展开严格证明了不可压缩NSAC系统向经典尖界面(Sharp-Interface)模型的收敛;Abels等[13]则从更一般的框架完善了界面极限分析。这些理论成果共同构成了NSAC系统在可压缩/不可压缩、多物理耦合及界面极限等方向上的核心基础,为我们开展稳态边值问题的研究及后续PINN数值验证提供了重要支撑。
2. 研究背景
我们研究一类非混相两相流扩散界面模型的一维情形,即一维等熵可压缩Navier-Stokes/Allen-Cahn (简称NSAC)系统的界面极限分析与模拟计算。考虑一维情况
,该系统在欧拉坐标系下有如下形式:
(1)
其中
表示混合流体的总密度,
表示混合流体的平均流速,
为两组分的浓度差值,
为化学势,
表示粘性系数而参数
表示其混合界面区域的厚度。问题(1)的稳态方程如下:
(2)
由
可得恒定质量通量
,于是
将
表达式代入式(2),可得到关于
和
的二阶常微分方程组:
(3)
其中
为积分常数。我们考虑稳态一维NSAC系统在区间
上如下边值条件:
稳态连续方程给出恒定质量通量
,代入动量方程后,可将稳态系统化为
其中
通过左端点动量平衡唯一确定为
。为避免直接处理端点条件并保证数值稳定性,我们对解进行如下参数化:
在
处:
,故修正函数
自动消失;
在
时逼近
,自然满足相场边值要求。因此
的边界条件全部被硬性嵌入,无需额外施加。系统中的第二个方程
是标准的Allen-Cahn稳态方程。令
,
,代入得到无量纲方程
其唯一满足
的解为
。这给出稳态相场在
下的精确界面结构。我们将相场近似视为已知场,将
的表达式代入第一方程,可写成:
该方程右端在
上平滑,且边界条件已经自动满足,故只需求解修正函数
的常微分方程。接下来将解析相场结构与速度方程的解结合,可得稳态解的如下形式:
其中
项来自
的求解。这些表达式提供了一个全域上一致精度为
的近似稳态解,同时满足了所有边界条件。
3. 解的存在性与稳定性分析
我们从方程组(3)出发,令
,则(3)可等价写为如下的一阶常微分方程组:
给定如下定解条件:
其中
是待定参数。将上述方程组改写为
(4)
其中
定理1 (局部存在唯一性与速度正性)设常数
,且边界速度
。取初始点
及初始值
,则存在一个含
的开区间
,使得由方程组(4)与初值
所构成的初值问题在
上存在且仅存在一条
解。进一步,若解在某开区间内满足
,则该解可延拓到整个
,且在整个区间上始终有
。
证明:选取矩形域
其中
,
,
足够大,使得初始值
落入
。在
上计算偏导数并证明有界性。对第一分量
:
对第二分量
,有
由于在
上
,
,
,因此这些偏导数在
上有上界(依赖于
等常数)。对第三分量
,有
由于在
上
且
,上述偏导数在
上也有界。由
在包含初值的闭有界域
上满足Lipschitz条件,可直接应用Picard-Lindelöf定理,即存在
使得初值问题在区间
上存在且唯一的
解,并且该解保持在
内。这证明了解的局部存在唯一性。
进一步,证明解的速度全局正性。令
可直接计算
由于
当
,故存在
使得
另一方面,系统的第一积分常数为
。因此当
时,
。稳态速度方程可写为
若假设存在
使得
,则在靠近
的区间内有
,由上式得
,从而
但由于
,且
在趋向零点之前为严格递减,与初值矛盾。故
在整个
上恒正。
4. 稳态问题的PINN方法
本节重点介绍基于物理信息神经网络(PINN)的方法及其在数值模拟中的应用,我们采用基于PyTorch搭建的PINN网络。
4.1. PINN方法概览
我们采用如下PINN架构与训练流程来求解稳态ODE边值问题(3)。令
(5)
其中
由含三层50单元的全连接网络输出,网络参数记为
。定义PDE残差为
构造加权残差损失
以及边界损失
。
最终总损失为
(6)
在
上等距采样
个collocation点,采用Adam优化器迭代5000轮,直至
收敛。图1展示了PINN求解流程示意图:
Figure 1. A PINN-based solution procedure for steady-state ODE boundary value problems
图1. 基于PINN的稳态ODE边值问题求解流程
4.2. 特殊情形:
为验证算法正确性,首先考察最特殊的情况
:此时质量通量
,方程(3)中对
的方程简化为
边界条件
。该简化模型对应纯相场平衡态,与经典Allen-Cahn方程的稳态解一致。图2和图3分别展示了
时PINN求解的
剖面与训练损失收敛曲线:
Figure 2. Profile of
obtained by the PINN under the case of
图2.
情形下PINN求得的
剖面
Figure 3. Convergence of the PINN training Loss
with respect to the number of iterations in the case
图3.
情形下PINN训练损失
随迭代轮数的收敛情况
在该最简模型中,PINN以高精度捕捉到
基解,验证了算法与代码实现的正确性。后续小节将进一步展示
时的完整耦合解及界面极限分析结果。
4.3. 一般情形:耦合解及界面极限分析
在验证了最简模型
的正确性之后,我们进一步求解完整的稳态方程组(3),考察
、不同界面宽度
下的耦合剖面与锐界面极限行为。取流量
、粘性系数
、常数
、
,边界速度
、
,界面宽度
,其余PINN超参与
情形一致。图4给出了不同
下的相场
剖面:
Figure 4. Profiles of
for different interface widths
图4. 不同界面宽度
的
剖面
由图可见,当
较大时(如0.5),界面平滑、过渡较宽;随着
,界面层迅速收缩至
宽度,并逼近
型分布。图5与图6分别给出对应的速度
和密度
分布:
Figure 5. Distributions of
for different interface widths
图5. 不同界面宽度
的
分布
Figure 6. Distributions of
for different interface widths
图6. 不同界面宽度
的
分布
可观察到,当
较大时,
平滑过渡于
至
,
相应平滑衔接;当
时,速度剖面在界面处呈现越来越尖锐的折线结构,密度在界面处出现尖峰后快速恢复至远场常值,吻合界面模型预期。为定量验证
的锐界面极限,定义相场误差为
,其中
。图7展示了
随
的变化:
Figure 7. Convergence behavior of the
error
with respect to the interface width
图7.
误差
随界面宽度
的收敛行为
结果表明误差呈
量级收敛,验证了PINN在锐界面极限下的数值稳定性。
至此,我们通过PINN方法成功获得了一维稳态NSAC系统在多种
情形下的数值解,并验证了相场与动量分布在锐界面极限的正确收敛性与物理一致性。后续将在第5节中进一步对比经典BVP解算器结果,并分析数值效率与精度。
5. 数值实验与结果分析
本节从多角度评估PINN方法在求解一维稳态NSAC常微分方程组边值问题上的性能,包括训练损失收敛、相场和速度解的误差指标、以及不同超参对结果的影响。首先考察总损失
随迭代轮数的变化。图8展示了
时,前5000轮迭代内
、
及
的收敛情况:
Figure 8. Convergence curve of the PINN training loss with respect to the number of iterations
图8. PINN训练损失随迭代轮数收敛曲线
由图8可见,总损失在迭代2000轮后趋于稳定,PDE残差损失降至10−5量级,说明PINN对稳态方程组的拟合精度满足要求。为量化PINN解的精度,定义
、
(其中
由经典BVP解算器如scipy.solve_bvp获得)。表1列出不同
情形下的
和
数值:
Table 1. Global
errors for different interface widths
表1. 不同界面宽度下的全域
误差
|
|
|
0.5 |
3.382 × 10−1 |
5.418 × 10−1 |
0.01 |
2.686 × 10−4 |
5.744 × 10−1 |
0.001 |
6.133 × 10−5 |
1.036 |
注:
时
突增,源于窄界面层下速度场梯度骤增,当前全域等距采样对高梯度区域覆盖不足,后续可通过界面层局部加密采样优化。面向PINN的采样策略至关重要。固定
,令collocation点数
,考察误差
和训练时间。图9给出
随
的变化:
Figure 9. Effect of the number of sampling points on
图9. 采样点数对
的影响
由图9可知,当
从1000增至10,000时,
从3 × 10−4降至6 × 10−5,且训练时间未呈线性增长(采样点数翻倍时,训练时间增幅约40%),体现了PINN的优势。为验证单调性与光滑性惩罚项的重要性,设置组合
,记录对应的
和
,结果如表2所示:
Table 2. Results of the penalty-weight ablation study
表2. 惩罚权重消融测试结果
|
|
|
|
0 |
0 |
6.983 × 10−4 |
4.653 × 10−1 |
0.1 |
0 |
3.873 × 10−4 |
6.584 × 10−1 |
0 |
0.001 |
6.019 × 10−4 |
6.626 × 10−1 |
0.1 |
0.001 |
1.024 × 10−3 |
4.424 × 10−1 |
在表2中,当
、
时,
显著高于
、
时的6.983 × 10−4,这一现象需结合物理约束与网络拟合特性解释:1) 平滑性惩罚项
过度“平滑”相场界面梯度,削弱了tanh型过渡特征,导致与参考解偏差增大;2) 单调性惩罚项
虽保证速度场单调性,但叠加平滑性惩罚后,网络优化目标偏向“整体平滑”而非“局部界面精度”。
综上,权重选择需权衡物理需求:优先保证界面精度可采用
,
;需避免速度反向波动则保留
,并将
降至0.0001以下。
至此,我们从多个角度验证了PINN方法在求解一维稳态NSAC常微分方程组边值问题上的收敛性与鲁棒性,并为后续与经典方法的比较奠定了量化基础。
需要指出的是,对于本文研究的一维稳态Navier-Stokes/Allen-Cahn系统边值问题,在纯正向求解的计算效率方面,PINN相比于传统的BVP求解器(如基于打靶法或配置法的经典算子)并无明显优势。传统求解器利用高度优化的数值线性代数库,通常能在秒级时间内获得收敛解,而PINN的训练则需要数千轮迭代。因此,PINN的核心贡献和应用潜力往往体现在非传统数值方法所长的领域,例如处理反问题(参数识别)等。特别地,反问题在实际工程应用中,系统的物理参数(如界面宽度
或粘性系数
)往往难以预先精确确定。传统的BVP方法在处理此类反问题时,通常需要复杂的灵敏度分析和嵌套优化循环。而PINN的框架可以无缝集成部分观测数据,通过在损失函数中添加数据残差项,将待识别参数设为网络的可学习变量,在求解场函数的同时直接反演物理参数。
为此,我们设计一个研究通过物理信息神经网络(PINN)开展参数反演的实验,旨在从稀疏且含噪声的观测数据中准确识别出难以直接测量的相场模型参数
。实验以稳态Allen-Cahn方程
为物理约束,在区间
内随机选取50个含1%高斯噪声的观测点,并将
设为网络的可学习参数进行优化。
Table 3. Results of PINN parameter inversion under different noise conditions
表3. 不同噪声条件下PINN参数反演结果
实验参数 |
真实值 |
初始值 |
PINN识别值 |
相对误差 |
无噪声数据 |
0.0100 |
0.0500 |
0.01008 |
0.80% |
高斯噪声数据 |
0.0100 |
0.0500 |
0.01052 |
5.20% |
实验结果如表3所示,在无噪声的理想条件下,PINN能够从严重偏离真实值400%的初始猜测(0.05)出发,成功反演出接近真实的参数值,相对误差控制在1%以内。然而,当观测数据中包含1%的高斯噪声时,反演精度出现显著下降,相对误差增至5.20%。这印证了噪声对反问题求解带来的固有挑战。尽管如此,该方法仍能从噪声数据中稳定地识别出具有明确物理意义的参数值,避免了传统梯度求解伴随方程的复杂性,展现了其作为工程反问题求解工具的潜力与鲁棒性。
6. 非稳态问题的推广
在上一节中,我们针对一维稳态NSAC常微分方程组边值问题验证了PINN方法的有效性。为了处理原始的非稳态PDE系统,本节介绍如何将上述框架扩展到含时间维度的全时域PINN。
非稳态PDE系统定义在时空区间
,采用随机采样与分段重点采样相结合的方式,内部点在
区域内随机采样
个collocation点
;边界点在两侧边界
上随机采样
个时间点
;初始点在初始面
上随机采样
个空间点
。
使用多输入多输出的全连接网络
,网络架构为:
(7)
其中
分别对应密度、速度、相场、化学势。总损失由四部分构成:
。在内部点
上,计算四个方程的残差:
(8)
并设置
。在点
上强制Dirichlet边界
、
,定义
在初始面
上匹配给定初值
、
、
,定义
为提升解的物理一致性,添加单调性和平滑性项:
。进行训练,对网络参数
初始化后,用零梯度在所有点上做一次应用,创建优化器内置slot;循环在所有三类点上计算
;用GradientTape自动求参,对
进行Adam更新;每隔若干轮计算损失与
误差,并在验证集上绘制曲线。图10和图11分别展示了
和
在
上的3D曲面效果图:
Figure 10. Spatiotemporal evolution of
in the unsteady NSAC system
图10. 非稳态NSAC系统中
时空演化
Figure 11. Spatiotemporal evolution of
in the unsteady NSAC system
图11. 非稳态NSAC系统中
时空演化
以上即为将PINN方法从稳态ODE推广到非稳态PDE的全过程,展示了采样策略、网络结构、损失设计和训练流程,并以典型实例验证了方法的可行性与准确性。
为了进一步定量评价全时空PINN模型在处理非稳态Navier-Stokes/Allen-Cahn系统时的精度,我们采用基于五阶WENO格式的空间离散和三阶Runge-Kutta时间步进的高分辨率有限体积法(FVM)作为基准解。计算网格数设为
,以确保能够完全解析
量级下的界面细微结构。
Figure 12. Comparison of relative
errors and interface trajectories between PINN and FVM ground truth for the unsteady NSAC system
图12. 非稳态NSAC系统下PINN与FVM基准解的
相对误差及界面轨迹对比
我们计算了相场变量
在整个模拟时域内的
相对误差。结果如图12所示,在界面演化的剧烈阶段,PINN的最大相对误差保持在10−3量级以下,证明了物理信息约束能够有效引导网络捕捉复杂的动力学过程。界面位置
定义为
的零水平集。图12中还展示了PINN预测的界面移动轨迹与FVM基准解的对比。两条曲线高度重合,能够准确捕捉到由于速度场
驱动引起的界面迁移以及最终向稳态位置靠拢的过程。这有力地证明了本文提出的PINN框架在解决全时空两相流问题上的有效性。
7. 结论
本研究围绕一维稳态Navier-Stokes-Allen-Cahn (NSAC)系统的界面极限问题展开,系统地分析了稳态扩散界面模型的数学性质与数值实现方法。通过推导稳态条件下的常微分方程组,并结合匹配渐近展开方法,获得了界面层内的tanh型基解,揭示了解在界面宽度
时的渐近行为。基于此,构建了精确的边界值问题(BVP)解算框架,并引入物理信息神经网络(PINN)方法,实现了对界面层结构的高效、稳定求解。数值结果验证了解在锐界面极限下的收敛性与稳定性,从数学上证明了稳态解的唯一性及其渐近稳定性,同时也探讨了广义Dirichlet边界条件对系统行为的影响。
总体而言,本研究在自由边界背景下深化了对NSAC系统稳态结构及其界面特征的理解,为复杂多相流体系中界面演化的建模与计算提供了理论依据与方法支撑。研究结果不仅在数学分析上具有理论意义,也为油气流动、水利工程、大气海洋模拟以及化工与材料加工等工程领域的多相界面计算提供了潜在应用价值。展望未来的工作,研究可以进一步扩展至更高维度的情形,探讨曲率效应与非平衡动力学下的界面行为,也可以结合时间依赖的NSAC系统,研究界面随流体运动的耦合机制与演化规律。同时,物理信息神经网络方法仍有进一步发展的空间,特别是在多尺度界面结构的高精度建模、自适应采样策略以及训练稳定性方面的改进,均值得深入探索。未来的研究还可在非均匀介质与复杂边界条件下,结合理论分析与数值模拟,揭示相场模型在更广泛实际问题中的适用性与普适规律,从而为多相流动的数学建模与智能计算开辟新的路径。
NOTES
*通讯作者。