1. 引言
粘弹性流体 [1] 存在于日常生活的方方面面,比如,食物中的番茄酱,蛋黄酱等,熔化的塑料,生物流体等。在工业上,向水中加入水溶性高分子聚合物可以得到粘弹性流体。粘弹性流体不仅具有粘性属性,又具有弹性行为,是一种介于粘弹性流体与弹性固体之间的特殊流体。粘弹性流体属于非牛顿流体,有一些特殊的物理现象,比如weissenberg爬杆效应,剪切稀变效应等。
不可压扩散peterlin粘弹性流体模型 [2] [3] [4] [5] [6] 是由动量守恒的NS方程,质量守恒的不可压缩条件
和利用peterlin 逼近描述构象张量d的本构方程耦合而成的,多变量、强非线性的方程组,如下
所示. 给定有界凸多边形区域
,有限时间T,求解速度u:
,压力p:
,对称构象张量d:
。
(1)
其中
表示张量d的迹,I是单位矩阵。物理参数
和
分别表示流体粘性系数和弹性张量系
数。为了简单起见,速度u和构象张量d赋于Dirichlet齐次边界条件,即
,
,并给定相
应的初始条件
,
。
文献 [2] [3] 讨论了该流体模型解的存在唯一性以及正则性。关于该流体模型的数值研究,可看文献 [4] [5] [6] ,文献 [5] 讨论了问题(1)时间离散一阶欧拉线性化压力校正投影格式,并进行了稳定化和先验误差分析。本文将继续讨论问题(1)的时间离散二阶BDF2压力校正投影格式。BDF2是线性多步法中的3步方法 [7] ,具有2阶收敛精度。压力修正投影格式,首先是由Chorin在1960年提出的 [8] ,主要是为了避开不可压缩流体的不可压条件。关于该方法针对不可压缩流体的综述文章见文献 [9] 。
本文结构如下。第二节给出该流体模型相应的预备知识和二阶压力校正投影格式。第三节对二阶BDF2压力校正投影格式进行了稳定性分析。第四节给出数值算例验证格式的稳定性。最后给出了小结。
2. 预备知识和数值格式
设
表示区域
上的勒贝格空间,其范数为
。符号
和
分别表示
和
空间的范数。
表示索伯列夫空间
,其范数为
。我们为了描述问题(1)的压力校正投影方
法,首先定义下列函数空间:
为了给出问题(1)的时间离散BDF2压力修正投影方法,首先对时间区间[0, T]进行等距剖分。设
是区间[0, T]的等距剖分,时间步长为
。记
和
。
格式2.1:(时间离散的BDF2压力校正投影)
第一步:通过下面的方程组求解
(2.1)
(2.2)
第二步:通过下面的方程组求解
(2.3)
注1. 由于BDF2是三步格式,需要知道前两步的函数值
。 我们可以用其它格式,例如Backward欧拉格式,得到
。
注2. 由于格式2.1是非线性全隐格式,因此在后面的算例计算时,需要进行线性化迭代,这里我们采用Newton 迭代。
注3. 为了使压力与速度解耦计算,本文用了BDF2压力校正投影方法。第一步是关于中间速度
和构象张量d的非线性椭圆型问题。第二步是关于速度和压力的线性方程组。
注4. 如果对第二步的第一个方程求散度
,可以得到
,即
。这是一个关于压力的Possion方程。在求得压力
后,可以通过下式得到
最终速度
于是时间离散的BDF2压力校正投影格式2.1,可以写成三步形式。
格式2.2:(时间离散的BDF2压力校正投影三步方法)
第一步:与格式2.1的第一步相同
第二步:通过下面的方程求解压力
第三步:通过下面的方程求解速度
3. 二阶BDF2压力校正投影格式2.1的稳定性分析
定理3.1. 设
是数值格式2.1的解。格式2.1几乎无条件稳定,即当时间步长
时,
有下面的稳定性结论
(3.1)
证明. 方程(2.1)两端乘以
并在区域
上积分,利用Green公式和三线性项
,我们得到
(3.2)
同理,方程(2.2)两端乘以
并在区域
上积分, 利用Green公式和三线性项
,我们得到
(3.3)
方程(3.2)和(3.3)相加并且利用关系式
可得
(3.4)
我们首先对式子(3.4)的等号右端项进行估计。利用关系式
可得
(3.5)
利用hölder和young’s不等式,我们得到
(3.6)
(3.7)
将不等式估计(3.5),(3.6)和(3.7)代入(3.4)可得:
(3.8)
从格式2.1的第二步(2.3),我们可以得到:
对上式两边分别取内积,可得:
(3.9)
应用关系式
(3.10)
和
(3.11)
我们将关系式(3.9),(3.10),(3.11)代入(3.8)得到:
(3.12)
我们对不等式(3.12)关于n从
到
求和可以得到
当时间步长
时,利用Gronwall不等式,可以得到定理3.1. 证明完毕
4. 数值算例
顶盖驱动方腔流(Lid Driven Cavity Flow)模型 [10] 常常用来检验代码和评估算法的优劣。计算区域
为单位正方形。本文给出如下的边界和初值条件。构象张量d的边界条件是
;在方腔的底端,以及左右边界上,水平速度u1和垂直速度u2都等于0;而在顶盖y = 1上,水平速度
,垂直速度u2 = 0。初值条件是u0 = 0,
。本文主要计算自由能(free energy)来验证算法的稳定性。
自由能定义如下:
分别取时间步长
时,参数
。时间方向离散采用有限差分,空间方向离散采用有限元。分别用(P1b, P1, P1)协调有限元来逼近速度,压力,构象张量。空间网格h = 1/80。
图1给出了三个不同的时间步长
,自由能随着时间t的动态演化过程。从图1容易看出,自由能曲线逐渐趋于稳定状态。
![](//html.hanspub.org/file/1-2620954x94_hanspub.png)
Figure 1.Dynamical evolution process of the free energy functional with t for three different time steps
图1. 对于三个不同时间步长,自由能函数随着时间t的动态演化过程
当时间步长
,和空间网格h = 1/80,图2和图3分别给出了最终时刻T = 5的构象张量d的各分量d11, d12, d22;压力p,速度u的分量u1, u2的轮廓图。这些结果与已有文献 [4] [5] [6] 完全吻合。说明了本文算法的数值稳定性。
![](//html.hanspub.org/file/1-2620954x96_hanspub.png)
Figure 2.From left to right, there are the contours d11, d12, d22 of the conformation tensor d, respectively
图2. 从左到右,依此是构象张量d的分量d11, d12, d22轮廓图
![](//html.hanspub.org/file/1-2620954x97_hanspub.png)
Figure 3. From left to right,there are the contours of pressure p, velocity u1, u2, respectively
图3. 从左到右,依次是压力p,速度u1, u2的轮廓图
5. 小结
本文对扩散Peterlin粘弹性流体提出了时间离散的BDF2压力校正投影方法,避开了速度与压力的不可压缩约束条件。本文证明了该格式的几乎(almost)无条件稳定性。最后通过顶盖驱动方腔流模型检验了该算法的稳定性。下一步我们将对该方法进行误差分析。
基金项目
河南科技大学大学生研究训练计划(SRTP)经费资助:2018188。