CN110532727A - 可用于常见非牛顿流体的数值模拟方法 - Google Patents
可用于常见非牛顿流体的数值模拟方法 Download PDFInfo
- Publication number
- CN110532727A CN110532727A CN201910847006.5A CN201910847006A CN110532727A CN 110532727 A CN110532727 A CN 110532727A CN 201910847006 A CN201910847006 A CN 201910847006A CN 110532727 A CN110532727 A CN 110532727A
- Authority
- CN
- China
- Prior art keywords
- fluid
- newtonian
- formula
- external force
- numerical simulation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明涉及流体动力学计算机模拟运算领域内的一种可用于常见非牛顿流体的数值模拟方法,首先基于格子玻尔兹曼方法建立含外力项的多松弛参数格子玻尔兹曼方法;之后根据三种常见非牛顿流体的流变方程,将它们的非牛顿特征转换成离散外力项,代入计算后,完成数值模拟过程。该方法可以有效改善应用格子玻尔兹曼方法模拟非牛顿流体时的稳定性和准确性,结合介观格子玻尔兹曼方法原理简单、计算方便、易于实现的特点,可指导非牛顿流体的流动仿真。
Description
技术领域
本发明涉及流体动力学计算机模拟运算领域,特别涉及一种非牛顿流体的数值模拟方法。
背景技术
随着流体动力学的不断发展,针对复杂的流体流动,涌现出一些不同于传统有限元计算的方法,如光滑粒子流方法、分子动力学方法以及格子玻尔兹曼方法等,从微观或介观角度更细致地诠释流体的流动过程,解决了利用传统有限元软件仿真,求解复杂微分方程组容易陷入局部解而无法求得最终解的难题,其中格子玻尔兹曼方法具有原理简单、计算方便、易于实现的特点,因此越来越广泛地用于解决复杂的流体流动。
中国专利数据库中,公开了一种移动粒子半隐式算法中自由表面流动模型的构建方法,其申请号:CN201210349290.1,申请日:20120919,公开号:CN102867094A,公开日:20130109,该方法引入的湍流模型包括静态Smagorinsky模型和拉格朗日形式的动态Smagorinsky模型;采用变黏度牛顿流体模型处理本构方程形如μ=f(|γ|)的非牛顿流体;引入三次样条核函数,采用光滑粒子流体动力学方法的散度离散方案离散非牛顿流体的剪切应力,处理非牛顿流体自由表面流动,其可用于非牛顿流体的计算。
现有技术的不足之处在于:对于非牛顿流体仿真,在利用普通格子玻尔兹曼方法时,会出现稳定性较弱,误差偏大的情况。
发明内容
本发明的目的是提供一种可用于常见非牛顿流体的数值模拟方法,可以有效改善格子玻尔兹曼方法模拟非牛顿流体流动时稳定性弱和误差偏大的不足。
本发明的目的是这样实现的:一种可用于常见非牛顿流体的数值模拟方法,包括步骤:
(1)基于格子玻尔兹曼方法,建立含有外力项的多松弛参数格子玻尔兹曼方法;
(2)根据三种常见非牛顿流体的流变方程,将它们的非牛顿特征转换成离散外力项,代入(1)中计算,完成数值模拟过程。
所述步骤(1)包括如下分步骤:
S2.1定义各物理量初值,及划分网格区域;
S2.2计算平衡态分布函数
式中feq(r,t)-流体在时刻t,位置r时平衡态分布函数;
ωi-权重系数,具体定义为
ρ-流体的密度;
ei-离散速度,具体定义为
u-流体的流动速度;
cs-格子声速,具体定义为
c-格子速度,具体定义为c=δx/δt;
δx-格子步长;
δt-时间步长;
S2.3计算碰撞步,具体公式为
式中f+(r,t)-流体执行碰撞步后在时刻t,位置r时的分布函数;
f(r,t)-流体在时刻t,位置r时的分布函数;
-中间变量,具体定义为
M-转换矩阵,具体定义为
-松弛过程相关的主对角线矩阵,具体定义为其中s8=1/τ;
τ-松弛时间,
ρ-动力粘度;
F’-计算外力项,具体定义为
I-单位矩阵;
-外力项分项,对于牛顿流体
S2.4计算应变率张量和剪切率,涉及公式如下:
式中Sαβ-应变率张量;
eiα,eiβ-分别为x方向和y方向离散速度ei的分量;
DII-应变率张量第二不变量;
l-仿真维度,这里l=2;
-剪切率;
S2.5计算迁移步,具体公式为f(r+eiδt,t+δt)=f+(r,t)
式中f(r+eiδt,t+δt)-流体在时刻t+δt,位置r+eiδt时的分布函数。
S2.6采用非平衡反弹形式进行壁面及边界处理;
S2.7计算相对误差若其小于10-4,则跳转至S2.9,否则执行S2.8;
S2.8跳至步骤S2.2,执行下一次循环计算;
S2.9计算所需的宏观量,所述宏观量至少包括速度和压力分布,并利用图形呈现具体的流线图、压力场和速度场分布图,完成模拟过程。
进一步地,所述内容(2)根据三种常见非牛顿流体的流变方程,将它们的非牛顿特征转换成特殊的外力项,代入(1)中计算,完成数值模拟过程,具体内容如下:
A幂律流体
对于幂律流体,其离散外力项可以表达为
式中μp0-幂律流体流变方程的粘度系数。
B宾汉流体
对于宾汉流体,其离散外力项可以表达为
式中μb0-宾汉流体流变方程的粘度系数;
τb0-宾汉流体流变方程的初始屈服应力;
m-一个与流体粘度、应力相关的参数;
-剪切率。
C Herschel-Bulkley流体
对于Herschel-Bulkley流体,其离散外力项可以表达为
式中μh0-Herschel-Bulkley流体流变方程的粘度系数;
τh0-Herschel-Bulkley流体流变方程的初始屈服应力;
m-一个与流体粘度、应力相关的参数;
-剪切率。
进一步地,本方法可通过编程实现,在利用MATLAB软件编程时,对于S2.2-S2.7步可利用Cell元胞数组计算代替for循环,可以有效提高运算效率。
与现有技术相比,本发明的有益效果在于:该方法可以有效改善应用格子玻尔兹曼方法模拟非牛顿流体时的稳定性和准确性,结合介观格子玻尔兹曼方法原理简单、计算方便、易于实现的特点,可指导常见非牛顿流体的流动仿真。
附图说明
图1为本发明的流程图。
具体实施方式
下面结合实施例及附图对本发明作进一步详细的描述,但本发明的的实施方式不限于此。
如图1所示,一种可用于常见非牛顿流体的数值模拟方法,其特征在于包括步骤:
(1)基于格子玻尔兹曼方法,建立含有外力项的多松弛参数格子玻尔兹曼方法;
(2)根据三种常见非牛顿流体的流变方程,将它们的非牛顿特征转换成离散外力项,代入(1)中计算,完成数值模拟过程。
所述步骤(1)包括如下分步骤:
S2.1定义各物理量初值,及划分网格区域;
S2.2计算平衡态分布函数
式中feq(r,t)-流体在时刻t,位置r时平衡态分布函数;
ωi-权重系数,具体定义为
ρ-流体的密度;
ei-离散速度,具体定义为
u-流体的流动速度;
cs-格子声速,具体定义为
c-格子速度,具体定义为c=δx/δt;
δx-格子步长;
δt-时间步长;
S2.3计算碰撞步,具体公式为
式中f+(r,t)-流体执行碰撞步后在时刻t,位置r时的分布函数;
f(r,t)-流体在时刻t,位置r时的分布函数;
-中间变量,具体定义为
M-转换矩阵,具体定义为
-松弛过程相关的主对角线矩阵,具体定义为其中s8=1/τ;
τ-松弛时间,
ρ-动力粘度;
F’-计算外力项,具体定义为
I-单位矩阵;
-外力项分项,对于牛顿流体
S2.4计算应变率张量和剪切率,涉及公式如下:
式中Sαβ-应变率张量;
eiα,eiβ-分别为x方向和y方向离散速度ei的分量;
DII-应变率张量第二不变量;
l-仿真维度,这里l=2;
-剪切率;
S2.5计算迁移步,具体公式为f(r+eiδt,t+δt)=f+(r,t)
式中f(r+eiδt,t+δt)-流体在时刻t+δt,位置r+eiδt时的分布函数。
S2.6采用非平衡反弹形式进行壁面及边界处理;
S2.7计算相对误差若其小于10-4,则跳转至S2.9,否则执行S2.8;
S2.8跳至步骤S2.2,执行下一次循环计算;
S2.9计算所需的宏观量,所述宏观量至少包括速度和压力分布,并利用图形呈现具体的流线图、压力场和速度场分布图,完成模拟过程。
步骤(2)中针对幂律流体,其离散外力项表达为
式中μp0-幂律流体流变方程的粘度系数。
步骤(2)中针对宾汉流体,其离散外力项表达为
式中μb0-宾汉流体流变方程的粘度系数;
τb0-宾汉流体流变方程的初始屈服应力;
m-一个与流体粘度、应力相关的参数;
-剪切率。
步骤(2)中针对Herschel-Bulkley流体,其离散外力项表达为
式中μh0-Herschel-Bulkley流体流变方程的粘度系数;
τh0-Herschel-Bulkley流体流变方程的初始屈服应力;
m-一个与流体粘度、应力相关的参数;
-剪切率。
利用MATLAB软件编程,对于S2.2-S2.7步利用Cell元胞数组计算代替for循环,以提高运算效率。
本发明并不局限于上述实施例,在本发明公开的技术方案的基础上,本领域的技术人员根据所公开的技术内容,不需要创造性的劳动就可以对其中的一些技术特征作出一些替换和变形,这些替换和变形均在本发明的保护范围内。
Claims (6)
1.一种可用于常见非牛顿流体的数值模拟方法,其特征在于包括步骤:
(1)基于格子玻尔兹曼方法,建立含有外力项的多松弛参数格子玻尔兹曼方法;
(2)根据三种常见非牛顿流体的流变方程,将它们的非牛顿特征转换成离散外力项,代入(1)中计算,完成数值模拟过程。
2.根据权利要求1所述的可用于常见非牛顿流体的数值模拟方法,其特征在于,所述步骤(1)包括如下分步骤:
S2.1定义各物理量初值,及划分网格区域;
S2.2计算平衡态分布函数
式中feq(r,t)-流体在时刻t,位置r时平衡态分布函数;
ωi-权重系数,具体定义为
ρ-流体的密度;
ei-离散速度,具体定义为
u-流体的流动速度;
cs-格子声速,具体定义为
c-格子速度,具体定义为c=δx/δt;
δx-格子步长;
δt-时间步长;
S2.3计算碰撞步,具体公式为
式中f+(r,t)-流体执行碰撞步后在时刻t,位置r时的分布函数;
f(r,t)-流体在时刻t,位置r时的分布函数;
-中间变量,具体定义为
M-转换矩阵,具体定义为
-松弛过程相关的主对角线矩阵,具体定义为其中s8=1/τ;
τ-松弛时间,
ρ-动力粘度;
F’-计算外力项,具体定义为
I-单位矩阵;
-外力项分项,对于牛顿流体
S2.4计算应变率张量和剪切率,涉及公式如下:
式中Sαβ-应变率张量;
eiα,eiβ-分别为x方向和y方向离散速度ei的分量;
DII-应变率张量第二不变量;
l-仿真维度,这里l=2;
-剪切率;
S2.5计算迁移步,具体公式为f(r+eiδt,t+δt)=f+(r,t)
式中f(r+eiδt,t+δt)-流体在时刻t+δt,位置r+eiδt时的分布函数。
S2.6采用非平衡反弹形式进行壁面及边界处理;
S2.7计算相对误差若其小于10-4,则跳转至S2.9,否则执行S2.8;
S2.8跳至步骤S2.2,执行下一次循环计算;
S2.9计算所需的宏观量,所述宏观量至少包括速度和压力分布,并利用图形呈现具体的流线图、压力场和速度场分布图,完成模拟过程。
3.根据权利要求2所述的可用于常见非牛顿流体的数值模拟方法,其特征在于,所述步骤(2)中针对幂律流体,其离散外力项表达为
式中μp0-幂律流体流变方程的粘度系数。
4.根据权利要求2所述的可用于常见非牛顿流体的数值模拟方法,其特征在于,所述步骤(2)中针对宾汉流体,其离散外力项表达为
式中μb0-宾汉流体流变方程的粘度系数;
τb0-宾汉流体流变方程的初始屈服应力。
5.根据权利要求2所述的可用于常见非牛顿流体的数值模拟方法,其特征在于,所述步骤(2)中针对Herschel-Bulkley流体,其离散外力项表达为
式中μh0-Herschel-Bulkley流体流变方程的粘度系数;
τh0-Herschel-Bulkley流体流变方程的初始屈服应力;
m-一个与流体粘度、应力相关的参数;
-剪切率。
6.根据权利要求2-5任一项所述的可用于常见非牛顿流体的数值模拟方法,其特征在于,利用MATLAB软件编程,对于S2.2-S2.7步利用Cell元胞数组计算代替for循环,以提高运算效率。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910847006.5A CN110532727B (zh) | 2019-09-09 | 2019-09-09 | 可用于常见非牛顿流体的数值模拟方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910847006.5A CN110532727B (zh) | 2019-09-09 | 2019-09-09 | 可用于常见非牛顿流体的数值模拟方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110532727A true CN110532727A (zh) | 2019-12-03 |
CN110532727B CN110532727B (zh) | 2023-05-23 |
Family
ID=68667660
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910847006.5A Active CN110532727B (zh) | 2019-09-09 | 2019-09-09 | 可用于常见非牛顿流体的数值模拟方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110532727B (zh) |
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111174648A (zh) * | 2020-01-13 | 2020-05-19 | 内蒙古工业大学 | 一种基于压板结构的降低药柱气孔率的浇注方法 |
CN112287622A (zh) * | 2020-12-28 | 2021-01-29 | 中国人民解放军国防科技大学 | 基于链接方向人工压缩的快速湍流数值模拟方法和装置 |
CN112765841A (zh) * | 2020-12-31 | 2021-05-07 | 西北工业大学 | 一种用以解决流体变物性计算的预处理格子玻尔兹曼方法 |
CN113051842A (zh) * | 2021-03-05 | 2021-06-29 | 清华大学 | 一种非牛顿流体模拟方法及装置 |
CN113468829A (zh) * | 2021-06-24 | 2021-10-01 | 西南石油大学 | 基于孔隙网络模型的非稳态非牛顿两相流体驱替模拟方法 |
-
2019
- 2019-09-09 CN CN201910847006.5A patent/CN110532727B/zh active Active
Non-Patent Citations (2)
Title |
---|
吴伟伟等: "Herschel-Bulkley流体的修正多松弛参数格子玻尔兹曼方法", 《工程热物理学报》 * |
王世博等: "微多孔介质非牛顿流体格子Boltzmann模拟", 《工程热物理学报》 * |
Cited By (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN111174648A (zh) * | 2020-01-13 | 2020-05-19 | 内蒙古工业大学 | 一种基于压板结构的降低药柱气孔率的浇注方法 |
CN111174648B (zh) * | 2020-01-13 | 2022-05-06 | 内蒙古工业大学 | 一种基于压板结构的降低药柱气孔率的浇注方法 |
CN112287622A (zh) * | 2020-12-28 | 2021-01-29 | 中国人民解放军国防科技大学 | 基于链接方向人工压缩的快速湍流数值模拟方法和装置 |
CN112287622B (zh) * | 2020-12-28 | 2021-03-09 | 中国人民解放军国防科技大学 | 基于链接方向人工压缩的快速湍流数值模拟方法和装置 |
CN112765841A (zh) * | 2020-12-31 | 2021-05-07 | 西北工业大学 | 一种用以解决流体变物性计算的预处理格子玻尔兹曼方法 |
CN112765841B (zh) * | 2020-12-31 | 2024-05-07 | 西北工业大学深圳研究院 | 一种用以解决流体变物性计算的预处理格子玻尔兹曼方法 |
CN113051842A (zh) * | 2021-03-05 | 2021-06-29 | 清华大学 | 一种非牛顿流体模拟方法及装置 |
CN113468829A (zh) * | 2021-06-24 | 2021-10-01 | 西南石油大学 | 基于孔隙网络模型的非稳态非牛顿两相流体驱替模拟方法 |
CN113468829B (zh) * | 2021-06-24 | 2024-04-26 | 西南石油大学 | 基于孔隙网络模型的非稳态非牛顿两相流体驱替模拟方法 |
Also Published As
Publication number | Publication date |
---|---|
CN110532727B (zh) | 2023-05-23 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110532727A (zh) | 可用于常见非牛顿流体的数值模拟方法 | |
Rosti et al. | Numerical simulations of emulsions in shear flows | |
CN104182560B (zh) | 飞行器颤振预测分析方法和装置 | |
Zhang et al. | Diffuse interface simulation of ternary fluids in contact with solid | |
Xu et al. | A level-set continuum method for two-phase flows with insoluble surfactant | |
Lee et al. | Two-dimensional Kelvin–Helmholtz instabilities of multi-component fluids | |
Aland et al. | Two-phase flow in complex geometries: A diffuse domain approach | |
Xie et al. | A balanced-force control volume finite element method for interfacial flows with surface tension using adaptive anisotropic unstructured meshes | |
Jiang et al. | A sharp-interface immersed smoothed finite element method for interactions between incompressible flows and large deformation solids | |
Filali et al. | Some experiences with the numerical simulation of Newtonian and Bingham fluids in dip coating | |
CN106372320B (zh) | 一种采用亚滤波尺度模型对公路隧道湍流进行大涡模拟的方法 | |
Guillaument et al. | An original algorithm for VOF based method to handle wetting effect in multiphase flow simulation | |
Reuther et al. | Incompressible two-phase flows with an inextensible Newtonian fluid interface | |
Aland et al. | Modeling and numerical approximations for bubbles in liquid metal | |
He et al. | Numerical simulation of interaction between multiphase flows and thin flexible structures | |
Rahmati et al. | Application of a modified pseudopotential lattice Boltzmann model for simulation of splashing phenomenon | |
Horton et al. | Benchmarking of computational fluid methodologies in resolving shear-driven flow fields | |
Li et al. | A mass-conserving lattice Boltzmann method for bubble behavior estimation | |
CN115879231B (zh) | 一种微型飞行器流动转捩信息的获取方法及装置 | |
Wang et al. | Large eddy simulation of stably stratified turbulent open channel flows with low-to high-Prandtl number | |
Zhang et al. | A comparative study of interface-conforming ALE-FE scheme and diffuse interface AMR-LB scheme for interfacial dynamics | |
Lee et al. | Full 3d simulations of two-phase core–annular flow in horizontal pipe using level set method | |
Garg et al. | A partitioned solver for compressible/incompressible fluid flow and light structure | |
Zuhal et al. | Core Spreading Vortex Method for Simulating 3D Flow around Bluff Bodies. | |
Wan et al. | Overset-rans computations of two surface ships moving in viscous fluids |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |