非牛顿流体在不规则管道中的流动特性研究

科技与创新┃ScienceandTechnology&Innovation

·24·2020年第08期

文章编号：2095－6835（2020）08－0024－05

非牛顿流体在不规则管道中的流动特性研究＊

李亚飞，周懿，胡钺，高政

（船舶动力工程技术交通行业重点实验室，湖北武汉430063；武汉理工大学能源与动力工程学院，湖北武汉430063）

摘要：非牛顿流体被大量应用于工业生产中，如利用剪切稀化流体应用的酱类食品生产运输、高分子聚合物在

工业生产中的管道运输及钻井用洗井液完井液等。相比传统的牛顿流体，关于非牛顿流体在管道内的流动的研究还有很大发展空间，因此对其机理的研究具有重要的现实意义。借助多物理场仿真软件COMSOLMultiphysics，对非牛顿流体在不规则管道内的流动现象进行了数值模拟研究，分析了非牛顿流体在管道内的速度场分布、剪切应变速率分布以及黏度分布情况。结果表明，在同一不规则管道中，管内流动的非牛顿Carreau流体的剪切应变速率与黏度关系符合其本构方程特性，相比牛顿流体表现出了明显较差的流动性，且速度场的分布一直保持一致；实验组非牛顿流体的剪切应变速率量级比对照组牛顿流体小很多，解释了两种流体在速度场分布的差异性；验证了Carreau流体构型中剪切应变速率与黏度的负相关关系。关键词：非牛顿流体；管道；流动特性；幂律流体

中图分类号：O373文献标识码：ADOI：10.15913/j.cnki.kjycx.2020.08.009

非牛顿流体是剪切力与剪切应变率之间不是线性关系的流体。非牛顿流体的应用十分普遍，涉及采用剪切增稠液体的装甲、食品工厂设计果酱类食品的运输加工及高分子聚合物在化工生产中的运输等诸多方面。这类非牛顿流体在工业生产方面的应用研究，直接影响工业生产效率及国家经济，因此对非牛顿流体在管道中的流动机理具有较大的研究价值。对非牛顿流体的研究方法包括实验法、解析解法与数值解法。实验法直接，但成本高且实验普遍性不高。解析解法理论上为最理想的研究方法，通过建立合适的微分方程，使用纯数学方法得出方程的精确解；但对于非牛顿流体复杂的流动情况来说，求解的难度过大。数值解法则是应用计算机将物理场离散化，之后将流体微分方程组转化为代数方程并求出各个节点上的参数值，属于一种近似解法。由于数值解法容易获得且能保证足够的求解精度，其已经成为研究非牛顿流体问题最为常用的方法。本文选取一种注塑喷嘴模型作为研究对象，以COMSOLMultiphysics这一软件作为计算平台，采用有限元

法对非牛顿流体流经不规则管道时表现出的流动特性进行数值模拟研究。在不规则腔体内，对照研究非牛顿Carreau流体与幂律牛顿流体的流动特性，讨论非牛顿流体的剪切应变速率与黏度对于速度大小与速度场分布的影响。1管内流动数学模型

1.1控制方程组一切的流体流动过程，都以三个基本的物理学原理为基础：质量守恒定律，牛顿第二定律与能量守恒定律。将这些物理学用于构建流动模型，将会导出一组方程，即连续性方程、动量方程与能量方程。这些方程是上述物理学原理的数学描述，本文不讨论传热，因此不引入能量守恒方程：质量守恒定律（连续性方程）：

0zwyvxut（1）

式（1）表示瞬态三维可压流体流动的连续性方程，本文所分析的流体流动处于稳态且不可压缩，密度ρ不会随着时间的变化而改变，因此流体流动的数学描述为：

0zwyvxu（2）动量守恒方程：

xzxyxxxFzyxxρuudivtρu）（（3）

yzyyyxyFzyxyvudivtρv）（（4）

zzzyzzFzyxzwudivtρwx）（（5）

式（3）（4）（5）是对于任何流体都成立的动量守恒方程，是微元体内流体动量对于事件的变化率等于外界作用于该微元体上的各种力的和，简称动量方程，也称纳维斯托克斯（N-S）方程。1.2本构方程本构方程是反映物料宏观性质的数学模型，又被称为流——————————————————————————

＊［基金项目］国家自然科学基金青年科学基金项目（编号：51806157）.com.cn. All Rights Reserved.ScienceandTechnology&Innovation┃科技与创新2020年第08期

·25·变状态方程或是流变方程。在流变学中，本构方程是在某些

假定条件下，对流体或弹性体的材料力学行为的数学描述，

可用来区分流体类型。本构方程与连续性方程、运动方程一

起构成封闭的方程组，用于求解流体的流动特性。

在本文仿真中使用的幂律流体的本构方程式如下：nm（6）

式（6）中：m为黏稠系数，表示物料的黏稠程度，Pa·sn；n

为幂律流变指数（简称幂律指数），为无量纲量，表示非牛

顿流体的流动特性偏离牛顿流体的程度（n=1时为牛顿

流体）。

1.3非牛顿流体管道流动模型

剪切应变速率描述的是流体的剪切流动，定义为单位时间的剪切应变变化：

dtd（7）

式（7）中：γ为剪切应变。

值得注意的是，剪切应变速率常与速度梯度混淆。实际

上二者是不同的概念。速度梯度是流体的速度对空间坐标的

导数，用du/dy来表示。在数学上，二者的数值有时相等，这是因为一般速度梯度符合：

dtddtdydL

dydtdL

dydu（8）

但二者的物理意义并不相同，且数值上有时并不相等

（如流体在同轴圆筒之间的流动，此时有角速度的影响）。2物理模型及计算方法

在规则管道内，流体的速度场分布较为均匀，各流动特

性参数的变化梯度较小。本文在这一部分选取一不规则管道

腔体进行管道流动特性仿真，以更清楚地展示流体的流动特

性。另选取适当的非牛顿流体与牛顿流体进行变量控制，构

成实验组与对照组。通过对比两种流体的仿真结果，说明非

牛顿流体不同于牛顿流体的流动特性。

2.1物理模型

2.1.1几何模型

本组仿真的计算区域如图1所示，为一个呈圆柱对称的

不规则弯曲管道腔体。流体由腔体的圆形入口处流入，中途

流经一截面当量直径较小的喉口，最后由腔体的环形出口处

流出。腔体的总长度为42mm，圆形入口直径为18mm，环

形出口宽度为12mm，其余具体几何尺寸如图2所示。A—B为管道进口段；F—G为管道出口段；B—C—D—E—F与

G—H—I—J—A为管道壁面。使用此不规则腔体进行数值模

拟的原因：在腔体喉口处由于其的截面当量直径相对较小，

流体经过时会出现较大的压力梯度，进而显著地影响速度、

剪切率与动力黏度。所以在喉口当量直径最窄处选取一个截

面，计算分析两种流体各项流动参数在这个截面上的分布，以更精确地说明流体的流动特性。截面具体位置如图3所

示。喉口处截面的位置与图2中R12与R3圆心连线重合。

图1不规则管道腔体示意图

图2不规则管道腔体的图3喉口处截面的位置

具体尺寸（单位：mm）（单位：mm）2.1.2流体性质

2.1.2.1非牛顿流体性质（实验组）

本组仿真计算所使用的非牛顿流体为线性聚苯乙烯在1-氯苯中的溶液（Asolutionoflinearpolystyrenein

1-chloronaphthalene），其被广泛用于塑料制品制造中。该溶

液在图1中不规则腔体内的流动过程，可视为对塑料制品生

产中注塑工序的溶液流经注塑喷嘴工段过程的模拟。该聚苯

乙烯溶液属于非牛顿Carreau流体，其黏度μeff与流体剪切应

变速率γ相关，具体关系式如下：

212inf0infeff1n）（）（）（（9）

式（9）中：μ0为零剪切应变速率时的黏度；μinf为无限大剪

切应变速率时的黏度；λ为松弛时间，s；n为模型指数。Carreau流体的剪切应变速率的计算公式如下：



22z2rz2r422221ruvvuu）（）（）（（10）

Carreau流体在低剪切应变速率（/1）时，表现出

类似牛顿流体的流动特性（黏度与剪切应变速率成正比）；

在高剪切应变速率时（/1），则表现出非牛顿流体的流

动特性。

本溶液的相关非牛顿Carreau流体模型参数值与仿真需

要的其他物性如表1所示。对照组牛顿流体的物性参数值如

表2所示。

2.1.2.2牛顿流体性质（对照组）

根据控制单一变量的原则，本组仿真所使用的对照组牛.com.cn. All Rights Reserved.科技与创新┃ScienceandTechnology&Innovation

·26·2020年第08期

顿流体只进行黏度相关系数的调整，其他物性参数与实验组

非牛顿流体保持一致。对照组牛顿流体使用幂律流体模型，

幂律指数n=1。当幂律指数n=1时，流体的黏度μeff实际上

固定为1Pa·s。表1线性聚苯乙烯在1-氯苯中的溶液的物性参数值参数名参数值零剪切应变速率时黏度μ00无限大剪切应变速率时黏度μinf/（Pa·s）166松弛时间λ/s1.73×10-2Carreau模型指数n0.583密度ρ/（kg/m3）450表2对照组牛顿流体的物性参数值参数名参数值稠度指数m1幂律指数n1最低剪切应变速率min/s-10.01密度ρ/（kg/m3）4502.2边界条件

进口边界条件：进口为图2中A—B段。设定进口压力pin，取值从10kPa开始，以40kPa的梯度递增至210kPa，

模拟不同进口压力下实验组与对照组流体的流动特性。

出口边界条件：出口为图2中F—G段，设定出口压力

为0。

壁面边界条件：壁面为图2中B—C—D—E—F与G—H—I—J—A段。设定无滑移边界条件，即壁面处的流体速

度为0。

无黏性应力条件：n·[μ（▽u+（▽u）T）]。其他：流体不可压缩；所有管道区域内实验组非牛顿流

体处于层流状态，进口压力较大时，对照组牛顿流体会出现

湍流。

2.3求解器介绍

本次模拟使用的COMSOLMultiphysics是一款功能强

大的多物理场仿真软件，可以对多个领域的物理过程进行模

拟计算。本文研究的是非牛顿流体在管道内的流动特性，因此选

用的是COMSOLMultiphysics中的Laminarflow求解器模

块。这个模块内置了许多流动的基本物理参数，可以在这些

参数的基础上结合前文所述的数学模型定义更多仿真所需

要的变量。

2.4数值计算方法本文涉及的数值计算使用的是有限元分析法（Finiteelementmethod，FEM）。它的基本思路为：一个物体或系统

被分解为由多个相互联结、简单、独立的点组成的几何模型。

在这种方法中这些独立的点的数量是有限的，因此被称为有

限元。由实际的物理模型中推导出来的平衡方程式被使用到

每个点上，由此产生了一个方程组，这个方程组可以用线性代数的方法来求解。3仿真结果分析

3.1两种流体的速度场对比

当进口压力分别为10kPa、90kPa、170kPa时，从实

验数据分析可知，随着进口压力的增大，牛顿流体的速度场

分布情况发生了较为明显的变化。进口压力为10kPa时，

流体速度在腔体入口段达到最大值，喉口处的中心速度也达

到了局部峰值，喉口段之后的流体速度渐缓。进口压力增大

时，腔体入口段的速度场分布情况差异不大，但喉口处的局