数学物理学报, 2026, 46(6): 2264-2276

非定常微极流体方程的一阶算子分裂方法

阿妮柯孜·奥斯曼,*

喀什大学 新疆喀什 844000

First-Order Operator Splitting Method for Unsteady Micropolar Fluid Equations

Anikezi Aosiman,*

Kashi University, Xinjiang Kashi 844000

通讯作者: 阿妮柯孜·奥斯曼, E-mail: 2720187148@qq.com

收稿日期: 2025-08-20   修回日期: 2026-09-15  

基金资助: 微极流体方的有限元方研究((2024) 2920)

Received: 2025-08-20   Revised: 2026-09-15  

Fund supported: Research on Finite Element Methods for Micropolar Fluid Equations((2024) 2920)

摘要

该文讨论 2D/3D 非定常微极流体方程的一阶算子分裂方法. 该方法在时间方向用一阶欧拉格式离散, 空间方向用协调有限元方法离散. 该文给出了一阶半离散算子分裂方法的无条件稳定性和误差估计, 并提出了相应的全离散格式. 证明了一阶全离散格式的无条件稳定性结论. 最后, 通过数值算例验证了该方法的准确性和有效性.

关键词： 微极流体方程; 算子分裂法; 有限元法; 稳定性分析; 误差估计

Abstract

This paper discusses the first-order operator splitting method for 2D/3D unsteady micropolar fluid equations. The method is discretized by the first-order Euler scheme in the time direction and the coordinated finite element method in the space direction. The unconditional stability and error estimation of the first-order semi-discrete operator splitting method are given, and the corresponding fully discrete scheme is proposed. The unconditional stability conclusion of the first-order fully discrete scheme is proved. Finally, the accuracy and effectiveness of the method are verified by numerical examples.

Keywords： microscale fluid equations; operator splitting method; finite element method; stability analysis; error estimation

PDF (808KB) 元数据 多维度评价 相关文章 导出 EndNote| Ris| Bibtex  收藏本文

本文引用格式

阿妮柯孜·奥斯曼. 非定常微极流体方程的一阶算子分裂方法[J]. 数学物理学报, 2026, 46(6): 2264-2276

Anikezi Aosiman. First-Order Operator Splitting Method for Unsteady Micropolar Fluid Equations[J]. Acta Mathematica Scientia, 2026, 46(6): 2264-2276

1 引言

微极流体方程是经典 Navier-Stokes 方程的推广, 当考虑各向异性或可极化流体时用角动量守恒来描述此类流体的流动. 例如, 用微极流体方程来描述动物血液, 液晶[1,2]和轴承润滑的流动[3,4]. 本文考虑了基于一阶算子分裂方法的非定常微极流体方程, 为了确保数值求解的收敛性和简化模型的复杂度, 本文选择区域 $\Omega\subset\mathbb{R}^{d} (d=$2 或 3) 是具有光滑边界 $\partial\Omega$ 的凸连通区域, 则非定常微极流体方程的控制方程为:

$\begin{equation} \left \{ \begin{array}{ll} u_t-\nu_{0}\Delta u +(u\cdot\nabla)u+\nabla p=f+2\nu_{r}\nabla\times w, \qquad & \quad\Omega\times(0,T],\\ \nabla \cdot u=0, \qquad & \quad\Omega\times(0,T],\\ j w_t-c_{1} \Delta w+j(u\cdot\nabla)w-c_{2}\nabla\nabla\cdot w+4\nu_{r} w=g+2\nu_{r}\nabla\times u, \qquad & \quad\Omega\times(0,T],\\ u=0,\quad w=0, \qquad & \quad\partial \Omega\times(0,T],\\ u=u_{0},\quad w=w_{0}, \qquad & \quad\Omega\times\{0\}. \end{array} \right. \end{equation}$

方程中未知函数是速度 $u$, 角速度 $w$, 压力 $p$. $\nu>0$ 是通常的粘性系数, $j$, $\nu_r$, $c_0$, $c_a$, $c_d$ 是与非对称应力张量有关的运动粘度系数, $c_a$, $c_d$ 和 $c_0$ 是满足不等式 $c_0 + c_d -c_a >0$. 方程中 $\nu_0=\nu+\nu_r$, $c_1=c_a+c_d$, $c_2 = c_0+ c_d-c_a$. 当空间维数 $d=2$ 时, 速度 $u=(u_1,u_2,0)$, 角速度 $w=(0,0,w_3)$. 当 $d=3$ 时, $u=(u_1,u_2,u_3)$, $w=(w_1,w_2,w_3)$[5,6]. 初边值条件如下:

$\begin{equation*} \left \{ \begin{array}{ll} u=u^{0},\quad w=w^{0}, \qquad & \quad\Omega\times\{0\},\\ u=0,\quad w=0, \qquad & \quad\partial \Omega\times(0,T], \end{array} \right. \end{equation*}$

微极流体方程是一个强耦合的非线性系统, 方程 (1.1) 中的第一个方程是线性动量守恒, 第二个方程是质量守恒, 第三个方程是角动量守恒方程. 方程 (1.1) 有强耦合性, 非线性性和不可压缩性等一些难点, 所以此方程的分析和求解是一个难题. 近年来, 微极流体方程的研究受到广泛关注. Sizopa 和 Qing 等人证明了二维微极流体方程弱解的存在唯一性[7,8]. Salgado[9] 提出了考虑微极流体方程的全离散分数时间步进法. Ortega-Torres 和 Rojas-Medar[10] 提出了全离散惩罚有限元法讨论微极流体方程并证明了最优误差估计. Nochetto 和 salago[11] 引入了一阶半隐式全离散有限元方法来研究微极流体方程. Bo 和 Jiang[12] 提供了微极流体方程的完全解耦方案和解耦惩罚投影方法. Maimaiti 和 Liu[13] 提出了时间相关微极流体方程的压力修正投影方法. Yang[14] 讨论了微极流体方程的一些投影方法. 投影边界条件与原始问题边界条件之间的不相容涉及压力上方的数值边界层[14]. 沈[15,16]讨论了大型耦合系统是个鞍点问题, 求解鞍点问题的成本高, 计算规模大等一些缺陷. 因此, 为微极流体方程设计一种高效、准确和小规模的数值解耦算法是有必要的.

算子分裂法应用于偏微分方程的数值求解. 特别是当算子分裂法用于求解不可压缩流动问题时, 将非线性和不可压缩分解为子步骤.

在每个子步骤中, 问题的原始边界条件被强制执行, 这能保证该方法的中间速度和步骤结束速度都收敛到连续解[17,18]. Dai[19] 提出了求解 Navier-Stokes 方程的二阶算子分裂方法, 并提供了沿时间方向的二阶误差估计. 在本文中, 采用类似的思想, 算子分裂方法推广来求解非定长微极性流体系统 (1.1).

2 预备知识

这里引入一般的标量 Sobolev 空间 $H^m(\Omega)$, $m=0,1,2$ 和向量 Sobolev 空间 $H^m(\Omega)^d$, 不引起混淆的情形下记它们的范数为 $\|\cdot\|_{m}$. 特别地, 当 $m=0$ 时, 令 $L^{2}(\Omega)=H^0(\Omega)$ 和 $L^{2}(\Omega)^d=H^0(\Omega)^d$, 具有对应的内积和范数为 $(\cdot,\cdot)$, $\|\cdot\|$, 记:

$\begin{equation*} \begin{array}{l} X={H}_{0}^1(\Omega)^d=\{v\in {H}^1(\Omega)^d,\quad v\mid_{\partial\Omega}=0\},\\ M=L_0^2(\Omega)=\{q\in L^2(\Omega),\quad \displaystyle{\int_\Omega q \text {d}x=0\}},\\ W={H}_{0}^1(\Omega)^d,\quad V=\{v\in X, \quad \nabla\cdot v=0\},\\ H=\{v\in L^{2}(\Omega)^d,\quad \nabla\cdot v=0,\quad v\cdot n|_{\partial\Omega}=0\}.\\ \end{array} \end{equation*}$

对于时间相关的序列 $\{\phi_i\}_{i=0}^n$, 引入以下一阶和二阶后向差分格式:

$\begin{array}{l} \delta_{t} \phi_{i+1}=\phi_{i+1}-\phi_{i}, \\ \delta_{t}^{2} \phi_{i+1}=\phi_{i+1}-2 \phi_{i}+\phi_{i-1}, \end{array}$

以及如下关系:

$\begin{equation*} \begin{array}{l} (\nabla\times w, u)=(w, \nabla\times u), \quad \forall u,w\in X,\\ \|\nabla u\|^{2}=\|\nabla\times u\|^{2}+\|\nabla\cdot u\|^{2}, \quad \forall u\in X,\\ \|\nabla\times u\|^{2}\leq \|\nabla u\|^{2}, \quad \forall u\in X. \end{array} \end{equation*}$

求 $(u,p,w)\in X\times M\times W$, 并 $(v,q,s)\in X \times M\times W$, 方程 (1) 的变分问题为:

$\left\{\begin{array}{l} \left(u_{t}, v\right)+B((u, p) ;(v, q))+N(u ; u, v)=(f, v)+d(w, v) \\ \left(w_{t}, s\right)+\tilde{a}(w, s)+j \tilde{N}(u ; w, s)+c(w, s)+4 \nu_{r}(w, s)=(g, s)+\tilde{d}(u, s) \end{array}\right.$

相应的双线性和三线性形式定义为

$\begin{equation*} \begin{array}{l} B((u,p);(v,q))=a(u,v)-b(v,p)+b(u,q),\quad \tilde{d}(u,s)=2\nu_r(\nabla\times u,s)\\ a(u,v)=\nu_{0}(\nabla u,\nabla v),\quad b(v,p)=(\nabla\cdot v,q),\quad\tilde{a}(w,s)=c_1(\nabla w,\nabla s)\\ c(w,s)= c_2(\nabla\cdot w,\nabla\cdot s),\quad d(w,v)=2\nu_r(\nabla\times w,v),\\ N(u;u,v)=(u\cdot\nabla)u,v),\quad \tilde{N}(u;w,s)=(j(u\cdot\nabla)w,s). \end{array} \end{equation*}$

关于三线性项有以下结论[22]

$\left\{\begin{array}{ll} N(u ; w, v)=-N(u ; v, w), & \forall u \in H, v, w \in X, \\ N(u ; v, v)=0, & \forall u \in H, v, w \in X, \end{array}\right.$
$\begin{equation*} |N(u;w,v)|\leq \left \{ \begin{array}{ll} \|\nabla u\|\|\nabla v\|\|\nabla w\|,&\quad \quad u, v, w\in X,\\ \|u\|_2\|\nabla v\|\|\nabla w\|, &\quad \quad u\in H^2(\Omega)^d, v\in X, w\in L^2(\Omega)^d,\\ \|u\|_2\|\nabla v\|\|w\|, &\quad \quad u\in L^2(\Omega)^d, v\in X, w\in L^2(\Omega)^d,\\ \|\nabla u\|\|v_2\|\|w\|, &\quad \quad u\in L^2(\Omega)^d, v\in H^2(\Omega)^d, w\in X.\\ \end{array} \right. \end{equation*}$

为了问题的讨论, 设初值和右端项满足如下正则性:

$\begin{equation*} \begin{array}{l} u^0\in H^2(\Omega)^d\cap V,\quad w^0\in H^2(\Omega)^d\cap X,\quad f,g\in L^\infty(0,T;L^2(\Omega)^d)\cap L^2(0,T;H^1(\Omega)^d). \end{array} \end{equation*}$

这里空间 $Z$ 是 Banach 空间, 我们使用缩写形式 $L^p(Z)=L^p(0,T;Z)$, 引入离散空间 $l^p(Z)$; 设 $w=\{w^n:n=1,\ldots,N\}$, 当 $p=\infty$ 时, 范数 $\|\cdot\|_{l^{(p)}(Z)}=\max_{0\leq n\leq N-1}\|w^{n+1}\|_Z$.

在上述假设下, 方程 (1.1) 的弱解满足如下正则性[20,21]:

引理 2.1

$\begin{equation} \left \{ \begin{array}{ll} \mathop{\sup}\limits_{t \in[T]}\|u(\cdot, t)\|_2+{\|w(\cdot, t)\|_2+\|u_t(\cdot, t)\|+\|w_t(\cdot, t)\|+\|\nabla p(\cdot, t)\|}\leq c,\\ \displaystyle{\int_0^T \|\nabla u_t(\cdot, t)\|^2 \text {d}t+\int_0^T\|\nabla w_t(\cdot, t)\|^2\| \text {d}t}\leq c.\\ \end{array} \right. \end{equation}$

这里 $c$ 是一个正常数, 它依赖于 $u^0, w^0, f, g$.

为了对时间相关序列进行分析, 引入以下离散 Gronwall 引理.

引理 2.2 (离散 Gronwall 引理) 设 $y^n$, $h^n$, $g^n$ 和 $f^n$ 是非负序列, 对 $m=0,1,\cdots$, $M=[T/\Delta t]$, 成立:

$\begin{equation*} \begin{array}{l} y^{m}+\Delta t\sum\limits _{n=0}^{m} h^{n} \leq B+\Delta t\sum\limits_{n=0}^{m}(g^{n}y^{n}+f^{n}),\qquad \Delta t\sum\limits_{n=0}^{M} g^{n}\leq K,\\ \end{array} \end{equation*}$

进一步假设 $\Delta t g^{n}<1$, 令 $\rho =\max\limits _{0\leq n \leq T/k}(1-\Delta t g^{n})^{-1}$. 则成立:

$\begin{equation*} \begin{array}{l} y^{m}+\Delta t\sum\limits_{n=0}^{m} h^{n}\leq {\rm e} ^{\rho K}(B+\Delta t\sum\limits_{n=0}^{m} f^{n}),\qquad 0\leq m\leq M. \end{array} \end{equation*}$

3 算子分裂方法

本节介绍了微极性流体方程的一阶算子分裂方法, 同时对一阶时间半离散方法的稳定性分析和误差估计进行了论证, 同时给出了一阶全离散算子分裂方法.

3.1 一阶时间半离散算子分裂方法:

算法 3.1 如 $({u}^n, p^{n}, w^{n})$, $(n=1,2\cdots$) 已知, 由以下三步求 $({u}^{n+1}, p^{n+1}, w^{n+1})$.

步骤 1 : 求一个中间速度 $ u^{n+\frac{1}{2}}\in X$:

$\begin{equation} \frac{u^{n+\frac{1}{2}}-u^n}{\Delta t}-\nu_0\Delta u^{n+\frac{1}{2}}+(u^n\cdot\nabla)u^{n+\frac{1}{2}}=f(t_{n+1})+2\nu_{r}\nabla\times w^{n}. \end{equation}$

步骤 2 : 求 $u^{n+1}\in X $ 和 $p^{n+1}\in M$:

$\begin{equation} \left \{ \begin{array}{ll} \frac{u^{n+1}-u^{n+\frac{1}{2}}}{\Delta t}-\nu_{0} \Delta(u^{n+1}-u^{n+\frac{1}{2}})+\nabla p^{n+1}=0,\\ \nabla\cdot u^{n+1}=0.\\ \end{array} \right. \end{equation}$

步骤 3 : 求 $w^{n+1}\in W$:

$\begin{equation} \begin{array}{ll} ~~~j\frac{w^{n+1}-{w}^{n}}{\Delta t}-c_{1} \Delta w^{n+1}+j(u^{n+1}\cdot\nabla)w^{n+1}-c_2\nabla\nabla\cdot w^{n+1}+4\nu_{r} w^{n+1}\\ =g(t_{n+1})+2\nu_{r}\nabla\times {u}^{n+1}.\\ \end{array} \end{equation}$

该算法在步骤 1 和步骤 3 中求解线性椭圆问题, 在步骤 2 中求解鞍点问题, 在步骤 1 和步骤 3 中分别求解动量守恒方程和角动量守恒方程.

定理 3.1 (稳定性) 假设 $0\leq N< M$, 算法 3.1 有一下无条件稳定性

$\begin{equation*} \begin{array}{l} ~~~|u^{N+\frac{1}{2}}\|^2+(j+4\nu_r\Delta t)\|w^N\|^2+\sum\limits_{n=0}^{N-1}(\|u^{n+1}-u^n\|^{2}+j\|w^{n+1}-w^{n}\|^{2})\\ ~~+\Delta t \sum\limits_{n=0}^{N-1}(\nu_0\|\nabla u^{n+\frac{1}{2}}\|^2+c_1\|\nabla w^{n+1}\|^2+c_2\|\nabla\cdot w^{n+1}\|^2)\\ \leq \|u^0\|^2+(j+4\nu_{r}\Delta t)\| w^0\|^2+\frac{ C }{\nu}\|f\|^{2}_{l^\infty(L^2)}+\frac{ C }{\nu}\|g\|^{2}_{l^\infty(L^2)}. \end{array} \end{equation*}$

证 用 $2\Delta t u^{n+\frac{1}{2}}$ 与 (3.1) 做内积, 用 $2\Delta t w^{n+1}$ 与 (3.3) 做内积, 有

$\begin{equation} \begin{array}{l} ~~~|u^{n+\frac{1}{2}}\|^{2}-\|u^n\|^{2}+\|u^{n+\frac{1}{2}}-u^n\|^{2}+2\nu_0\Delta t\|\nabla u^{n+\frac{1}{2}}\|^2\\ ~~+ ((u^n\cdot\nabla)u^{n+\frac{1}{2}},2\Delta t u^{n+\frac{1}{2}})\\ =2\Delta t(f(t_{n+1}),u^{n+\frac{1}{2}})+4\nu_{r}\Delta t(\nabla\times w^{n},u^{n+\frac{1}{2}}). \end{array} \end{equation}$
$\begin{equation} \begin{array}{l} ~~~j\|w^{n+1}\|^{2}-j\|w^{n}\|^{2}+j\|w^{n+1}-w^{n}\|^{2}+2c_{1}\Delta t\|\nabla w^{n+1}\|^{2}\\ ~~+2c_{2}\Delta t\|\nabla\cdot w^{n+1}\|^{2}+8\nu_{r}\Delta t\|w^{n+1}\|^{2}\\ =4\nu_{r} \Delta t(\nabla\times u^{n+1},w^{n+1})+2\Delta t(g(t_{n+1}),w^{n+1}). \end{array} \end{equation}$

(3.4) 和 (3.5) 相加并用非线性项的性质, 得

$\begin{equation} \begin{array}{l} ~~~|u^{n+\frac{1}{2}}\|^{2}-\|u^n\|^{2}+\|u^{n+\frac{1}{2}}-u^n\|^{2}+2\nu_0\Delta t\|\nabla u^{n+\frac{1}{2}}\|^2\\ ~~+j\|w^{n+1}\|^{2}-j\|w^{n}\|^{2}+j\|w^{n+1}-w^{n}\|^{2}\\ ~~+2c_{1}\Delta t\|\nabla w^{n+1}\|^{2}+2c_{2}\Delta t\|\nabla\cdot w^{n+1}\|^{2}+8\nu_{r}\Delta t\|w^{n+1}\|^{2}\\ =4\nu_{r}\Delta t(\nabla\times w^{n},u^{n+\frac{1}{2}})+2\Delta t(f(t_{n+1}),u^{n+\frac{1}{2}})\\ \quad+4\nu_{r} \Delta t(\nabla\times u^{n+1},w^{n+1})+2\Delta t(g(t_{n+1}),w^{n+1}). \end{array} \end{equation}$

估计 (3.6) 右侧项, 有

$ \begin{align*} 4\nu_{r} \Delta t(\nabla\times w^{n}, u^{n+\frac{1}{2}})&=4\nu_{r} \Delta t(w^n,\nabla\times u^{n+\frac{1}{2}})\\ &\leq 4\nu_{r} \Delta t\|w^{n}\| \|\nabla\times u^{n+\frac{1}{2}}\|\\ &\leq 4\nu_{r} \Delta t\|w^{n}\|^{2}+\nu_{r} \Delta t \|\nabla u^{n+\frac{1}{2}}\|^{2},\\ 2\Delta t(f(t_{n+1}),u^{n+\frac{1}{2}})&\leq 2\Delta t\|f(t_{n+1})\|\|\nabla u^{n+\frac{1}{2}}\|\\ &\leq \frac{\Delta t C}{\nu}\|f(t_{n+1})\|^2+\nu\Delta t\|\nabla u^{n+\frac{1}{2}}\|^2,\\ 4\nu_{r} \Delta t(\nabla\times u^{n+1},w^{n+1})&\leq 4\nu_{r} \Delta t\|w^{n+1}\| \|\nabla\times u^{n+1}\|\\ &\leq 4\nu_{r} \Delta t\|w^{n+1}\|^{2}+\nu_{r} \Delta t \|\nabla u^{n+1}\|^{2},\\ 2\Delta t(g(t_{n+1}),w^{n+1})&\leq 2\Delta t\|g(t_{n+1})\|\|\nabla w^{n+1}\|\\ &\leq \frac{\Delta t C}{\nu}\|g(t_{n+1})\|^2+\nu\Delta t\|\nabla w^{n+1}\|^2. \end{align*} $

将上述估计代入 (3.6), 我们有

$\begin{equation*} \begin{array}{l} ~~~|u^{n+\frac{1}{2}}\|^{2}-\|u^n\|^{2}+\|u^{n+\frac{1}{2}}-u^n\|^{2}+\nu\Delta t\|\nabla u^{n+\frac{1}{2}}\|^2\\ ~~+(j-4\nu_r\Delta t)\|w^{n+1}\|^2-(j+4\nu_r\Delta t)\|w^n\|^2+j\|w^{n+1}-w^{n}\|^{2}\\ ~~+c_{1}\Delta t\|\nabla w^{n+1}\|^{2}+2c_2\|\nabla\cdot w^{n+1}\|^2\\ \leq \frac{ C \Delta t}{\nu}\|f(t_{n+1})\|^2+\frac{ C \Delta t}{\nu}\|g(t_{n+1})\|^{2}. \end{array} \end{equation*}$

当 $ n=0,1,\cdots, N-1$ 时, 对上面的式子求和并用离散形式的 Gronwall 引理, 则有:

$\begin{equation*} \begin{array}{l} ~~~|u^{N+\frac{1}{2}}\|^2+(j-4\nu_r\Delta t)\|w^{N+1}\|^2+\sum\limits_{n=0}^{N-1}(\|u^{n+\frac{1}{2}}-u^n\|^{2}+j\|w^{n+1}-w^{n}\|^{2})\\ ~~+\Delta t \sum\limits_{n}^{N-1}(\nu\|\nabla u^{n+\frac{1}{2}}\|^2+c_1\|\nabla w^{n+1}\|^2+c_2\|\nabla\cdot w^{n+1}\|^2)\\ \leq \|u^0\|^2+(j+4\nu_{r}\Delta t)\| w^0\|^2+\frac{ C }{\nu}\|f\|^{2}_{l^\infty(L^2)}+\frac{ C }{\nu}\|g\|^{2}_{l^\infty(L^2)}. \end{array} \end{equation*}$

现在讨论算法 3.1 的速度和角速度的最佳误差估计. 为了进行收敛性分析, 方便起见, 引入以下符号:

$\begin{equation*} \begin{array}{l} e_u^{n+\frac{1}{2}}=u(t_{n+1})-u^{n+\frac{1}{2}}, \qquad e_u^{n+1}=u(t_{n+1})-u^{n+1},\qquad e_p^{n+1}=p(t_{n+1})-p^{n+1},\\ e_w^{n+1}=w(t_{n+1})-w^{n+1}, \qquad \xi_w^n=w(t_{n+1})-w^{n}. \end{array} \end{equation*}$

最后介绍了速度和角速度的误差估计.

定理 3.2 由正则性假设, 有 $0\leq N< M$, 则

$\begin{equation*} \begin{array}{l} \|e_u^{N+1}\|^2+j\|e_w^{N+1}\|^2+\sum\limits_{n=0}^{N}(\|e_u^{n+1}-e_u^n\|^2+j\|e_w^{n+1}-e_w^n\|^2)\\ +\Delta t\sum\limits_{n=0}^{N}(\nu_0\|e_u^{n+1}\|^2+c_1\|e_w^{n+1}\|^2+2c_2\|\nabla\cdot e_w^{n+1}\|^2)\leq C\Delta t^2. \end{array} \end{equation*}$

证 由 (3.1) 和 (3.2) 的第一个方程得:

$\begin{equation} \begin{array}{l} ~~~\frac{u^{n+1}-u^n}{\Delta t}-\nu_0 \Delta u^{n+1}+(u^n\cdot\nabla)u^{n+\frac{1}{2}}+\nabla p^{n+1}\\ =f(t_{n+1})+2\nu_{r}\nabla\times w^{n}. \end{array} \end{equation}$

令方程 (1.1) 中的 $t=t_{n+1}$, 有:

$\begin{equation} \left \{ \begin{array}{ll} \displaystyle~~~\frac{u(t_{n+1})-u(t_{n})}{\Delta t}-\nu_{0}\Delta u(t_{n+1})+(u(t_{n+1})\cdot\nabla)u(t_{n+1})+\nabla p(t_{n+1})\\ =2\nu_{r}\nabla\times w(t_{n+1})+f(t_{n+1})+R_{u}^{n},\\ \displaystyle~~~j\frac{w(t_{n+1})-w(t_{n})}{\Delta t}-c_{1}\Delta w(t_{n+1})+j(u(t_{n+1})\cdot\nabla)w(t_{n+1})-c_2\nabla\nabla\cdot w(t_{n+1})\\~~+4\nu_{r}w(t_{n+1})\\=2\nu_{r}\nabla\times u(t_{n+1})+g(t_{n+1})+R_{w}^{n},\\ \nabla\cdot u(t_{n+1})=0, \end{array} \right. \end{equation}$

这里

$\begin{equation*} R_u^{n}=-\frac{1}{\Delta t}\int_{t_n}^{t_{n+1}}(t-t_{n})u_{tt}(t)\text{d}t,\qquad R_w^{n}=-\frac{1}{\Delta t}\int_{t_n}^{t_{n+1}}(t-t_{n})w_{tt}(t)\text{d}t. \end{equation*}$

从 (3.8) 中的第一个方程减去 (3.7), 从 (3.8) 中的第二个方程减去 (3.3), 我们有

$\begin{equation} \left \{ \begin{array}{ll} \displaystyle~~~\frac{e_{u}^{n+1}-e_{u}^{n}}{\Delta t}-\nu_{0}\Delta e_{u}^{n+1}+(u(t_{n+1})\cdot\nabla)u(t_{n+1})-(u^n\cdot\nabla)u^{n+\frac{1}{2}}+\nabla e_{p}^{n+1}\\ =2\nu_{r}\nabla\times \xi_w^n+R_{u}^{n},\\ \displaystyle~~~j\frac{e_{w}^{n+1}- e_{w}^{n}}{\Delta t}-c_{1}\Delta e_{w}^{n+1}+j(u(t_{n+1})\cdot\nabla)w(t_{n+1})-(u^{n+1}\cdot\nabla)w^{n+1}\\ ~~-c_2\nabla\nabla\cdot e_{w}^{n+1}+4\nu_{r} e_{w}^{n+1}\\ =2\nu_{r}\nabla\times e_{u}^{n+1}+R_{w}^{n},\\ \nabla\cdot e_{u}^{n+1}=0. \end{array} \right. \end{equation}$

用 $2\Delta t e_{u}^{n+1}$ 与 (3.9)的第一个方程做内积, 用 $2\Delta t e_{w}^{k+1}$ 与第二个方程做内积, 我们有:

$\begin{matrix} &\quad\|e_{u}^{n+1}\|^2-\|e_{u}^{n}\|^{2}+\||e_{u}^{n+1}-e_{u}^{n}\|^{2}+2\nu_0\Delta t\|\nabla e_{u}^{n+1}\|^2\nonumber\\ &=2\Delta t((u^n\cdot\nabla)u^{n+\frac{1}{2}},e_{u}^{n+1})-2\Delta t((u(t_{n+1})\cdot\nabla)u(t_{n+1}), e_{u}^{n+1})\\ & +4\nu_{r}\Delta t(\nabla\times \xi_w^n, e_{u}^{n+1})+2\Delta t(R_{u}^{n},e_{u}^{n+1})\nonumber \end{matrix}$

和

$\begin{matrix} &~~~~j\|e_{w}^{n+1}\|^{2}-j\|e_{w}^{n}\|^{2}+j\||e_{w}^{n+1}-e_{w}^{n}\|^{2}-2c_{1}\Delta t(\Delta e_{w}^{n+1},e_{w}^{n+1})\nonumber\\ &~~~+2c_{2}\Delta t \|\nabla \cdot e_{w}^{n+1}\|^{2}+8\nu_r\Delta t\|e_w^{n+1}\|^2\\ &=2\Delta t(u^{n+1}\cdot\nabla)w^{n+1},e_{w}^{k+1})-2\Delta t((u(t_{n+1})\cdot\nabla)w(t_{n+1}),e_{w}^{n+1})\nonumber\\ &~~~+4\nu_{r}\Delta t(\nabla\times e_{u}^{n+1},e_{w}^{n+1})+2\Delta t(R_{w}^{n},e_{w}^{n+1}).\nonumber \end{matrix}$

将 (3.10) 右侧的非线性项的估计为:

$\begin{equation} \begin{array}{l} ~~~(u^n\cdot\nabla)u^{n+\frac{1}{2}}-(u(t_{n+1})\cdot\nabla)u(t_{n+1})\\ =((u(t_{n})-u(t_{n+1}))\cdot\nabla)u^{n+\frac{1}{2}}-(e_u^n\cdot\nabla)u^{n+\frac{1}{2}}- (u(t_{n+1})\cdot\nabla)e_u^{n+\frac{1}{2}},\\ \end{array} \end{equation}$

我们有

$\begin{matrix} & \|e_{u}^{n+1}\|^2-\|e_{u}^{n}\|^{2}+\||e_{u}^{n+1}-e_{u}^{n}\|^{2}+2\nu_0\Delta t\|\nabla e_u^{n+1}\|^2\nonumber\\ &=2\Delta tb((u(t_{n})-u(t_{n+1});u^{n+\frac{1}{2}},e_u^{n+1}) -2\Delta t b(e_u^n;u^{n+\frac{1}{2}},e_u^{n+1})\nonumber\\ & -2\Delta t b(u(t_{n+1});e_u^{n+\frac{1}{2}},e_u^{n+1})+4\nu_{r}\Delta t(\nabla\times \xi_w^n,e_u^{n+1})+2\Delta t(R_{u}^{n},e_u^{n+1})\nonumber\\ &=\sum\limits_{i=1}^{5} A_i\nonumber \end{matrix}$

和

$\begin{matrix} & j\|e_{w}^{n+1}\|^{2}-j\|e_{w}^{n}\|^{2}+j\||e_{w}^{n+1}-e_{w}^{n}\|^{2}+2c_{1}\Delta t\|\nabla e_{w}^{n+1}\|^{2}+2c_{2}\Delta t \|\nabla \cdot e_{w}^{n+1}\|^{2}+8\nu_r\Delta t\|e_w^{n+1}\|^2\nonumber\\ &=-2\Delta t \tilde{b}(e_{u}^{n+1};w(t_{n+1}),e_{w}^{n+1})+2\Delta t \tilde{b}(u^{n+1};e_{w}^{n+1},e_{w}^{n+1})\\ & +4\nu_{r}\Delta t(\nabla\times e_{u}^{n+1},e_{w}^{n+1})+2\Delta t(R_{w}^{n},e_{w}^{n+1})=\sum\limits_{i=1}^{4} B_i.\nonumber \end{matrix}$

利用引理 2.1, 柯西不等式和非线性项的结论, 我们讨论 (3.13) 的右侧项, 即

$\begin{equation*} \begin{array}{l} A_1=2\Delta tb((u(t_{n})-u(t_{n+1});u(t_{n+1}),e_u^{n+1}))\\ ~~~~=-2\Delta tb((u(t_{n})-u(t_{n+1});e_u^{n+1},u(t_{n+1})))\\ ~~~~\leq C\Delta t\|u(t_{n})-u(t_{n+1})\|\| \nabla e_u^{n+1}\|\|u(t_{n+1})\|_2\\ ~~~~\leq c\Delta t \|\int_{t_n}^{t_{n+1}}u_t \text{d}t\|\nabla e_u^{n+1}\|\\ ~~~~\leq \frac{c \Delta t}{\nu_0}^2\int_{t_n}^{t_{n+1}}\|u_t\|^2 \text{d}t+\frac{\nu_0 \Delta t}{4}\| \nabla e_u^{n+1}\|^2,\\ A_2=-2\Delta tb(e_u^{n};u^{n+\frac{1}{2}},e_u^{n+1})=-2\Delta tb(e_u^{n};e_u^{n+1},u^{n+\frac{1}{2}})\\ ~~~~\leq C\Delta t\|e_u^{n}\|\|\nabla e_u^{n+1}\|\|u^{n+\frac{1}{2}}\|_2\\ ~~~~\leq c\Delta t\|e_u^{n}\|\|\nabla e_u^{n+1}\|\\ ~~~~\leq \frac{c \Delta t}{\nu_0}\|e_u^{n}\|^2+\frac{\nu_0 \Delta t}{4}\| \nabla e_u^{n+1}\|^2,\\ A_3=-2\Delta t b(u(t_{n+1});e_u^{n+\frac{1}{2}},e_u^{n+1}) =2\Delta t b(u(t_{n+1});e_u^{n+1},e_u^{n+\frac{1}{2}})\\ ~~~~\leq C\Delta t\|u(t_{n+1})\|_2\|\nabla e_u^{n+1}\|\|e_u^{n+\frac{1}{2}}\| \leq c\Delta t\|\nabla e_u^{n+1}\|\|e_u^{n+\frac{1}{2}}\|\\ ~~~~\leq \frac{\nu_0 \Delta t}{4}\|\nabla e_u^{n+1}\|^2+\frac{c \Delta t}{\nu_0}\|e_u^{n+\frac{1}{2}}\|^2,\\ A_4=4\nu_{r}\Delta t(\nabla\times \xi_w^n, e_u^{n+1})=4\nu_{r}\Delta t(\xi_w^n,\nabla\times e_u^{n+1})\\ ~~~~\leq C\nu_r\Delta t\|\xi_w^n\|\|\nabla\times e_u^{n+1}\|\leq \frac{\nu_0\Delta t}{4}\|\xi_w^n\|^2+\frac{\nu_r C\Delta t}{\nu_0}\|e_u^{n+1}\|^2,\\ A_5=2\Delta t(R_{u}^{n},e_u^{n+1})\leq 2\Delta t\|R_u^n\|\|e_u^{n+1}\|_1 \leq \frac{\nu_0 \Delta t}{4}\|\nabla e_u^{n+1}\|^2+C\Delta t \|R_u^n\|^2_{-1}\\ ~~~~=\frac{\nu_0 \Delta t}{4}\|\nabla e_u^{n+1}\|^2+C\Delta t^{-1}\|\int_{t_n}^{t_{n+1}}(t-t_n)u_{tt}\|_{-1}^2\\ ~~~~\leq \frac{\nu_0 \Delta t}{4}\|\nabla e_u^{n+1}\|^2+C\Delta t^{-1}\|\int_{t_n}^{t_{n+1}}\|u_{tt}\|_{-1}^2\text{d}t\int_{t_n}^{t_{n+1}}(t-t_n)^2\text{d}t\\ ~~~~\leq \frac{\nu_0 \Delta t}{4}\|\nabla e_u^{n+1}\|^2+\frac{ C\Delta t^2}{\nu_0}\int_{t_n}^{t_{n+1}}\|u_{tt}\|_{-1}^2\text{d}t,\\ \end{array} \end{equation*}$
$\begin{equation*} \begin{array}{l} B_1=-2\Delta t \tilde{b}(e_{u}^{n+1};w(t_{n+1}),e_{w}^{n+1})) =2\Delta t \tilde{b}(e_{u}^{n+1};e_{w}^{n+1},w(t_{n+1}))\\ ~~~~\leq C\Delta t\|e_u^{n+1}\|\|\nabla e_{w}^{n+1}\|\|w(t_{n+1})\|_2\leq c\Delta t\|e_u^{n+1}\|\|\nabla e_{w}^{n+1}\|\\ ~~~~\leq \frac{c_1\Delta t}{4}\|\nabla e_{w}^{n+1}\|^2+\frac{c \Delta t}{c_1}\|e_u^{n+1}\|^2,\\ B_2=2\Delta t \tilde{b}(u^{n+1};e_{w}^{n+1},e_{w}^{n+1})=0,\\ B_3=4\nu_{r}\Delta t(\nabla\times e_{u}^{n+1},e_{w}^{n+1})\\ ~~~~\leq C\nu_r\Delta t\|e_u^{n+1}\|\|\nabla\times e_{w}^{n+1}\| \leq \frac{c_1\Delta t}{4}\|e_u^{n+1}\| ^2+C\Delta t\|\nabla e_{w}^{n+1}\|^2\\ \hspace{-4.6em} \text{和}\\ B_4=2\Delta t(R_{w}^{n},e_u^{n+1})\leq 2\Delta t\|R_w^n\|\|e_w^{n+1}\|_1 \leq \frac{c_1 \Delta t}{4}\|\nabla e_w^{n+1}\|^2+C\Delta t \|R_u^n\|^2_{-1}\\ ~~~~=\frac{c_1 \Delta t}{4}\|\nabla e_w^{n+1}\|^2+C\Delta t^{-1}\|\int_{t_n}^{t_{n+1}}(t-t_n)w_{tt}\|_1^2\\ ~~~~\leq \frac{c_1 \Delta t}{4}\|\nabla e_w^{n+1}\|^2+C\Delta t^{-1}\|\int_{t_n}^{t_{n+1}}\|w_{tt}\|_{-1}^2\text{d}t\int_{t_n}^{t_{n+1}}(t-t_n)^2\text{d}t\\ ~~~~\leq \frac{c_1 \Delta t}{4}\|\nabla e_w^{n+1}\|^2+C\Delta t^2\int_{t_n}^{t_{n+1}}\|w_{tt}\|_{-1}^2\text{d}t. \end{array} \end{equation*}$

将估计值 $A_i$ 和 $B_i$ 代入 (3.13), 我们有

$\begin{equation} \begin{array}{l} ~~~\|e_u^{n+1}\|^2-\|e_u^n\|^2+\|e_u^{n+1}-e_u^n\|^2+\nu_0\Delta t\|e_u^{n+1}\|^2+j(\|e_w^{n+1}\|^2\\ \quad-\|e_w^n\|^2+\|e_w^{n+1}-e_w^n\|^2)+\frac{c_1}{2}\Delta t\|e_w^{n+1}\|^2+2c_2\Delta t\|\nabla\cdot e_w^{n+1}\|^2\\ \leq \frac{c\Delta t^2}{\nu_0}\int_{t_n}^{t_{n+1}}\|u_{t}\|^2\text{d}t+\frac{c \Delta t}{\nu_0}\|e_u^{n}\|^2 +\frac{c\Delta t}{\nu_0}\| e_u^{n+\frac{1}{2}}\|^2+\frac{\nu_r\Delta t}{\nu_0}\| e_u^{n+1}\|^2\\ \quad+\frac{c \Delta t^2}{\nu_0}\int_{t_n}^{t_{n+1}}\|u_{tt}\|_{-1}^2\text{d}t+\frac{c \Delta t^2}{\nu_0}\int_{t_n}^{t_{n+1}}\|w_{tt}\|_{-1}^2\text{d}t. \end{array} \end{equation}$

当 $ n=0,1,\cdots, N-1$ 时, 对上面的式子求和并用离散形式的 Gronwall 引理, 则有:

$\begin{equation} \begin{array}{l} ~~~\|e_u^{N+1}\|^2+j\|e_w^{N+1}\|^2+\sum\limits_{n=0}^{N}(\|e_u^{n+1}-e_u^n\|^2+j\|e_w^{n+1}-e_w^n\|^2\\ \quad+\Delta t\sum\limits_{n=0}^{N}(\nu_0\|e_u^{n+1}\|^2+c_1\|e_w^{n+1}\|^2+2c_2\|\nabla\cdot e_w^{n+1}\|^2)\\ \leq C\Delta t^2+ C\Delta t\sum\limits_{n=0}^{N}(\|e_u^{n+1}\|^2+\|e_w^{n+1}\|^2). \end{array} \end{equation}$

用离散的 Gronwall 引理 2.1, 有:

$\begin{equation} \begin{array}{l} ~~~~\|e_u^{N+1}\|^2+j\|e_w^{N+1}\|^2+\sum\limits_{n=0}^{N}(\|e_u^{n+1}-e_u^n\|^2+j\|e_w^{n+1}-e_w^n\|^2)\\ ~~~+\Delta t\sum\limits_{n=0}^{N}(\nu_0\|e_u^{n+1}\|^2+c_1\|e_w^{n+1}\|^2+2c_2\|\nabla\cdot e_w^{n+1}\|^2)\\ \leq C\Delta t^2. \end{array} \end{equation}$

3.2 一阶全离散算子分裂法:

设 $ \tau_h=\{K\}$ 是 $\Omega$ 的一致正则三角网格, 且网格大小 $h=\max _{K\in\tau_h}\text{\{diam(K)\}}$ 中, 定义以下离散子空间:

$\begin{equation*} \begin{array}{l} X_h=\{v_h\in X\cap C^0(\overline \Omega)^d:v_h|_K\in P1b(K)^d,\quad \forall K\in \tau_h\},\\ M_h=\{q_h\in M \cap C^0(\overline \Omega)^d:q_h|_K\in P_1(K),\quad \forall K\in \tau_h\},\\ W_h=\{s_h\in W \cap C^0(\overline \Omega)^d:s_h|_K\in P1b(K)^d,\quad \forall K\in \tau_h\}. \end{array} \end{equation*}$

算法 3.2 一阶全离散算子分裂方法

步骤 1 : 已知 $\forall v_h\in X_h$, 求 $ u_h^{n+\frac{1}{2}}\in X_h$ 则有:

$\begin{equation} \begin{array}{l} (\frac{u_h^{n+\frac{1}{2}}-u_h^n}{\Delta t},v_h)+a( u_h^{n+\frac{1}{2}},v_h)+N(u_h^n; u_h^{n+\frac{1}{2}},v_h) =(f(t_{n+1}),v_h)+d(w^{n},v_h); \end{array} \end{equation}$

步骤 2 : 已知 $ u_h^{n+\frac{1}{2}}$, 求 $u_h^{n+1}\in X_h$ 和 $p_h^{n+1}\in M_h$ 则有:

$\begin{equation} \begin{array}{l} (\frac{u_h^{n+1}-u_h^{n+\frac{1}{2}}}{\Delta t},v_h)-\nu_{0} (\Delta(u_h^{n+1}-u_h^{n+\frac{1}{2}}),v_h)-b( p_h^{n+1},v_h)=0; \end{array} \end{equation}$

步骤 3 : 已知 $\forall s_h \in W_h$, 求 $w_h^{n+1}\in W_h$ 则有:

$\begin{equation} \begin{array}{l} ~~~(j\frac{w_h^{n+1}-w_h^{n}}{\Delta t},s_h)+\tilde a( w_h^{n+\frac{1}{2}},s_h)+\tilde N(u_h^{n+1};w_h^{n+1},s_h)+c(w_h^{n+1},s_h)+4\nu_{r} (w_h^{n+1},s_h)\\ =(g(t_{n+1}),s_h)+\tilde d( u_h^{n+1},s_h). \end{array} \end{equation}$

定理 3.3 (稳定性) 算法 3.2 有以下无条件稳定性:

$\begin{equation*} \begin{array}{l} ~~~\|u_h^{N+\frac{1}{2}}\|^2+(j+4\nu_r\Delta t)\|w_h^N\|^2++\sum\limits_{n=0}^{N-1}(\|u_h^{n+1}-u_h^n\|^{2}+j\|w_h^{n+1}-w_h^{n}\|^{2})\\ \quad+\Delta t \sum\limits_{n=0}^{N-1}(\nu\|\nabla u_h^{n+\frac{1}{2}}\|^2+c_1\|\nabla w_h^{n+1}\|^2+c_2\|\nabla\cdot w_h^{n+1}\|^2)\\ \leq \|u_h^0\|^2+(j+4\nu_{r}\Delta t)\| w_h^0\|^2+\frac{ C}{\nu}\|f\|^{2}_{l^\infty(L^2)}+\frac{ C}{\nu}\|g\|^{2}_{l^\infty(L^2)}. \end{array} \end{equation*}$

证 令$v_h=2\Delta t u_h^{n+\frac{1}{2}}$ 有:

$\begin{equation} \begin{array}{l} ~~~\|u_h^{n+\frac{1}{2}}\|^{2}-\|u_h^n\|^{2}+\|u_h^{n+1}-u_h^n\|^{2}+2\nu_0\Delta t\|\nabla u_h^{n+\frac{1}{2}}\|^2\\ =2\Delta t(f(t_{n+1}),u_h^{n+\frac{1}{2}})+4\nu_{r}\Delta t(\nabla\times w_h^{n},u_h^{n+\frac{1}{2}}). \end{array} \end{equation}$

类似地, 令 $s_h=2\Delta t w_h^{n+1}$, 我们有:

$\begin{equation} \begin{array}{l} ~~~j\|w_h^{n+1}\|^{2}-j\|w_h^{n}\|^{2}+j\|w_h^{n+1}-w_h^{n}\|^{2}+2c_{1}\Delta t\|\nabla w_h^{n+1}\|^{2}\\ \quad+2c_{2}\Delta t\|\nabla\cdot w_h^{n+1}\|^{2}+8\nu_{r}\Delta t\|w_h^{n+1}\|^{2}\\ =4\nu_{r} \Delta t(\nabla\times\tilde u_h^{n+1},w_h^{n+1})+2\Delta t(g(t_{n+1}),w_h^{n+1}). \end{array} \end{equation}$

将 (3.20) 和 (3.21) 相加, 得:

$\begin{matrix} &~~~\|u_h^{n+\frac{1}{2}}\|^{2}-\|u_h^n\|^{2}+\|u_h^{n+1}-u_h^n\|^{2}+2\nu_0\Delta t\|\nabla u_h^{n+\frac{1}{2}}\|^2 +j\|w_h^{n+1}\|^{2}-j\|w_h^{n}\|^{2}+j\|w_h^{n+1}\nonumber\\ &~~~-w_h^{n}\|^{2}+2c_{1}\Delta t\|\nabla w_h^{n+1}\|^{2}+2c_{2}\Delta t\|\nabla\cdot w_h^{n+1}\|^{2} +8\nu_{r}\Delta t\|w_h^{n+1}\|^{2}\\ &=4\nu_{r}\Delta t(\nabla\times w_h^{n},u_h^{n+\frac{1}{2}})+2\Delta t(f(t_{n+1}),u_h^{n+\frac{1}{2}}) +4\nu_{r} \Delta t(\nabla\times u_h^{n+1},w_h^{n+1})\nonumber\\&~~~+2\Delta t(g(t_{n+1}),w_h^{n+1}).\nonumber \end{matrix}$

我们有类似于定理 3.1 的右侧项估计:

$\begin{equation*} \begin{array}{l} ~~~\|u_h^{N+\frac{1}{2}}\|^2+(j+4\nu_r\Delta t)\|w_h^N\|^2++\sum\limits_{n=0}^{N-1}(\|u_h^{n+1}-u_h^n\|^{2}+j\|w_h^{n+1}-w_h^{n}\|^{2})\\ ~~~+\Delta t \sum\limits_{n=0}^{N-1}(\nu_0\|\nabla u_h^{n+\frac{1}{2}}\|^2+c_1\|\nabla w_h^{n+1}\|^2+c_2\|\nabla\cdot w_h^{n+1}\|^2)\\ \leq \|u_h^0\|^2+(j+4\nu_{r}\Delta t)\| w_h^0\|^2+\frac{ C}{\nu}\|f\|^{2}_{l^\infty(L^2)}+\frac{ C}{\nu}\|g\|^{2}_{l^\infty(L^2)}. \end{array} \end{equation*}$

该方法的稳定性估计似于算法 3.1 的稳定性估计.

4 数值实验

在本节中, 给出了几个数值实验. 考虑了二维和三维微极流体方程符合有限元逼近的收敛性测试. 给出了二维三维微极流体方程的能量稳定性. 为了表明所提方法的有效性和可靠性, 给出了圆柱绕流问题. 通过数值实验, 证明算子分裂方法的准确性和效率. 选择 $P1b-P1-P1b$ 有限元对来构造有限元空间.

4.1 收敛性验证

实例 1 : 二维微极流体方程的收敛性检验.

考虑了具有光滑解析解的问题. 方程 (1.1) 的右侧由精确解决定. 还给出了一阶运算方法的误差和收敛速度.

$\begin{equation*} \begin{array}{l} u_1(x,y,t)=-\sin(\pi y)\sin^2(\pi x)\cos(\pi y),\\ u_2(x,y,t)=\sin(\pi x)\sin^2(\pi y)\cos(\pi x),\\ p(x,y,t)=\cos(\pi y)\sin(\pi x),\\ w(x,y,t)=u_1(x,y,t)+u_2(x,y,t). \end{array} \end{equation*}$

表1   一阶格式的误差和收敛速度, $\Delta t=h^2$, $T$ =1.

新窗口打开| 下载CSV


实例 2 : 三维微极流体方程的收敛性检验.

在这个例子中, 给出了三维微极流体方程的收敛性测试, 并给出了以下精确解

$\begin{equation*} \begin{array}{l} u(x,y,z,t)=((x^4+z^2)\cos(t),(z^4+x^2)\cos(t),(x^4+y^2)\cos(t)),\\ p(x,y,z,t)=(2x-1)(2y-1)(2z-1)\cos(t),\\ w(x,y,z,t)=((\sin(y)+z)\cos(t),(\sin(z)+x)\cos(t),(\sin(x)+y)\cos(t)). \end{array} \end{equation*}$

表2   一阶格式的误差和收敛速度, $\Delta t=h^2$, $T$ =1

新窗口打开| 下载CSV


4.2 2D/3D 微极性流体方程的稳定性测试

当 $\Omega=[0,1]^d$, $d=2, 3$, 考虑了一阶时间离散算子分裂法的稳定性. 令 $f=g=0$ 并选择以下 2D/3D 初始条件:

$\begin{equation*} \begin{array}{l} u_1=2\pi\sin(\pi x)^2\sin(\pi y)\cos(\pi y),\quad u_2=-2\pi\sin(\pi x)\sin(\pi y)^2\cos(\pi x), \\ w_1=0.5\pi(\sin(2\pi y) - \sin(2\pi x) + \sin(2\pi(x - y))) \end{array} \end{equation*}$

和

$\begin{equation*} \begin{array}{l} u_1=2\cos(\pi y)\sin(\pi x)\cos(\pi z),\quad u_2=-\sin(\pi y)\cos(\pi x)\cos(\pi z),\\ u_3=-\cos(\pi y)\sin(\pi z)\cos(\pi x),\quad w_1=(1-\cos(2\pi x))\sin(2\pi y)\sin(2\pi z),\\ w_2=\sin(2\pi x)(1-\cos(2\pi y))\sin(2\pi z),\quad w_3=\sin(2\pi x)(1-\cos(2\pi z))\sin(2\pi y). \end{array} \end{equation*}$

这里, $T$=1, $h=1/60$. 给出了一阶半离散算子分裂方法的能量耗散过程图.对于能量系统, 有 $E_n:=\|u^n\|^2+\|v^n\|^2+\|w^n\|^2$, 给出了不同时间步长: $\Delta t=0.1, 0.05, 0.025, 0.0125$, $T$=1 时的能量耗散过程. 我们观察到, 所有个能量曲线都显示出所有时间步长的单调递减.

图1

新窗口打开| 下载原图ZIP| 生成PPT

图1   二维 (左) 和三维 (右) 微极流体方程不同时间步长的能量稳定性


4.3 圆柱绕流问题

我们的下一个数值试验是关于二维问题求解的圆柱体绕流. 是一个众所周知的基准问题. 这个流体流动问题不是湍流, 但确实有一些有趣的特点流动模式是由流体与壁的相互作用驱动的, 区域 $\Omega= [0,2.2]\times[0,0.41]$ 是流动区域[23,24]. 本实例中 $c_a=c_d=c_0=1$, $(\nu+\nu_r)=1/1000$, $\Delta t=0.002$. 与时间相关的流入和流出剖面是:

$\begin{equation} u_1(0,y,t)=u_1(2.2,y,t)=\frac{6}{0.41^2}\sin(\frac{\pi t}{8}y(0.41-y)), u_2(0,y,t)=u_2(2.2,y,t)=0. \end{equation}$

图2

新窗口打开| 下载原图ZIP| 生成PPT

图2   $T$=2, 4, 5, 6, 7 和 8 秒时的速度$u$.


如图所示, 随着流入量的增加, 两个涡流开始在圆柱体后面形成, (见 T=2 和 4s). 在 T=4 和 5s 之间, 涡流与圆柱体分开, 漩涡街形成. 漩涡在最后时间 T=8s.

5 总结

本文用算子分裂求解究微极性流体方程, 证明了该方法是无条件稳定性, 并给出了速度和角速度的误差分析. 数值实验表明算子分裂方法是有效的, 系统能量是耗散的, 并给出了圆柱绕流问题. 将在我们今后的工作中讨论中, 可以进一步考虑微极性流体方程的二阶全离散算子分裂方法及其理论分析.

参考文献

Dahler J S, Scriven L E.

Angular momentum of continua

Nature, 1961, 192: 36-37

DOI:10.1038/192036a0      [本文引用: 1]

Eringen A C. Microcontinuum Field Theories: I. Foundations and Solids. New York: Springer-Verlag, 1999, 116

[本文引用: 1]

Allen S J, Kline K A.

Lubrication theory for micropolar fluids

J Appl Mech, 1971, 38: 646-650

DOI:10.1115/1.3408868      URL     [本文引用: 1]

The equations governing the flow of a fluid with rigid, spherical substructure are summarized. A two-dimensional flow field is considered and applied to the geometry of a slider bearing. Order-of-magnitude arguments are used which reduce the governing equations to a system of coupled, linear, ordinary differential equations. The equations are solved subject to appropriate boundary conditions and the effects of substructure discussed with the help of a specific numerical example.

Tipei N.

Lubrication with micropolar liquids and its application to short bearings

J Lubr Technol, 1979, 101: 356-363

DOI:10.1115/1.3453375      URL     [本文引用: 1]

Microscopic effects, generated by micromotions of particles in suspension in a viscous fluid, drastically change the character of the flow between solid walls. To the modified momentum and continuity equations, an equation of angular (spin) particle momentum is added. A vectorial system of equations is presented, for variable material coefficients. General properties of this system are discussed and differential equations for pressure and velocity field are derived. For constant viscosity and micropolar coefficients across the lubricating film, important simplifications lead to easier workable expressions. Under this assumption, the short bearing performance has been analyzed. The fluid pressure increase (compared to Newtonian flow) is represented by a surface depending on two groups of micropolar parameters; the overall bearing characteristics also exhibit larger values with respect to simple Newtonian lubricants. However, for similar gap geometries, the friction coefficient has lower values. Some formulas, regarding the velocity field, friction stresses and side flow, are general and may be applied to any bearing length-diameter ratio.

Scriven J S.

Theory of structured continua I. General consideration of angular momentum and polarization

Proc R Soc, 1963, 275: 504-527

[本文引用: 1]

Łukaszewicz G. Micropolar Fluids:Theory and Applications. Birkhäuser Boston: MA, 1999

[本文引用: 1]

Szopa P.

On existence and regularity of solutions for 2-D micropolar fluid equations with periodic boundary conditions

Math Method Appl Sci, 2007, 30: 331-346

DOI:10.1002/mma.v30:3      URL     [本文引用: 1]

Dong B Q, Zhang Z F.

Global regularity of the 2D micropolar fluid flows with zero angular viscosity

J Differential Equations, 2010, 249: 200-213

DOI:10.1016/j.jde.2010.03.016      URL     [本文引用: 1]

Salgado A J.

Convergence analysis of fractional time-stepping techniques for incompressible fluids with microstructure

J Sci Comput, 2015, 64: 216-233

DOI:10.1007/s10915-014-9926-x      URL     [本文引用: 1]

Ortega-Torres E, Rojas-Medar M.

Optimal error estimate of the penalty finite element method for the micropolar fluid equations

Numer Funct Anal Optim, 2008, 29(5-6): 612-637

DOI:10.1080/01630560802099555      URL     [本文引用: 1]

Nochetto R H, Salgado A J.

The micropolar Navier-Stokes equations: A priori error analysis

Math Models Methods Appl Sci, 2014, 24: 1237-1264

DOI:10.1142/S0218202514500018      URL     [本文引用: 1]

Bo Y Y, Jiang Y L.

Analysis of two decoupled time-stepping finite-element methods for incompressible fluids with microstructure

Int J Comput Math, 2018, 95(4): 686-709

DOI:10.1080/00207160.2017.1294688      URL     [本文引用: 1]

Maimaiti H, Liu D.

Pressure-correction projection methods for the time dependent micropolar fluids

Int J Numer Meth Fluids, 2022, 94(4): 377-393

DOI:10.1002/fld.v94.4      URL     [本文引用: 1]

Jiang Y L, Bo Y Y.

Analysis of some protection methods for the incompressible fluids with microstructure

J Korean Math Soc, 2018, 55: 471-506

[本文引用: 2]

Shen J.

Velocity-correction projection methods for incompressible flows

SIAM J Numer Anal, 2003, 41: 112-134

DOI:10.1137/S0036142901395400      URL     [本文引用: 1]

Shen G J. On the error estimates for the rotational pressure-correction projection methods. Math Comput, 2004, 73: 1719-1737

[本文引用: 1]

Blasco J, Codina R.

A fractional-step method for the incompressible Navier-Stokes equations related to a predictor-multicorrector algorithm

Int J Numer Meth, 1998, 28: 1391-1419

DOI:10.1002/(ISSN)1097-0363      URL     [本文引用: 1]

Blasco J, Codina R.

Error estimates for an operator-splitting method for incompressible flows

Appl Numer Math, 2004, 51: 1-17

[本文引用: 1]

Dai X, Sun J, Cheng X.

Error estimates for an operator-splitting method for Navier-Stokes equations: Second-order schemes

J Comput Appl Math, 2009, 231: 696-704

DOI:10.1016/j.cam.2009.04.014      URL     [本文引用: 1]

Brenner S C, Scott L R. The Mathematical Theory of Finite Element Methods. New York: Springer-Verlag, 1994

[本文引用: 1]

Labovsky A, Layton W J, Manica C C, et al.

The stabilized extrapolated trapezoidal finite-element method for the Navier-Stokes equations

Computer Methods in Applied Mechanics and Engineering, 2009, 198(9-12): 958-974

DOI:10.1016/j.cma.2008.11.004      URL     [本文引用: 1]

Shen J.

On error estimates of some higher order projection and penalty-projection methods for Navier-Stokes equations

Numerische Mathematik, 1992, 62(1): 49-73

DOI:10.1007/BF01396220      URL     [本文引用: 1]

Lu X, Zhang L, Huang P.

A fully discrete finite element scheme for the Kelvin-Voigt model

Filomat, 2019, 33(18): 5813-5827

DOI:10.2298/FIL1918813L      [本文引用: 1]

In this paper, we study convergence of a fully discrete scheme for the two-dimensional nonstationary Kelvin-Voigt model. This scheme is based on a finite element approximation for space discretization and the Crank-Nicolson-type scheme for time discretization, which is a two step method. Moreover, we obtain error estimates of velocity and pressure. At last, the applicability and effectiveness of the present algorithm are illustrated by numerical experiments.

Turek S.

Benchmark computations of laminar flow around a cylinder//Hirschel E H. Flow Simulation with High-Performance Computers II

Wiesbaden: Vieweg, 1996: 547-566

[本文引用: 1]

/

〈 〉

版权所有 © 《数学物理学报》杂志社
地址：武汉市71010号信箱(430071)　电话: 027-87199206(中文版)　027-87199087(英文版)
E-mail: actams@wipm.ac.cn
本系统由北京玛格泰克科技发展有限公司设计开发