番茄酱、血液等非牛顿流体会展现出与应变速率相关的特性，这使得模拟通过多孔介质的流动变得复杂。多孔材料的复杂结构,如通道、停滞区和孤立孔隙会导致流体特性各不相同。在没有通用理论的情况下，常常需要通过测量或孔隙尺度仿真的方法开发出特定的流体-材料组合模型。在这篇博客中，我们演示了如何利用仿真或测量结果开发出一种均质化方法，用于模拟非牛顿流体在多孔结构中的流动。
非牛顿流方法
日常生活中某些流体会表现出反直觉特性。一个众所周知的例子就是番茄酱，刚开始它会顽固地附着在瓶壁，然后突然在整个盘子上铺满红色的酱汁。当黏度随着剪切速率的增加而降低时，就会出现这种被称为“剪切稀化”的现象。类似地，大多数聚合物也会出现这种现象。较少见的是剪切增稠流体，即黏度随剪切速率的增加而增加。一个著名的例子是 Oobleck（一种玉米淀粉悬浮液），其黏度会随剪切速率增加，甚至可使人在其表面快速行走时如履平地。
剪切速率表征相邻层间流体的相对运动速度，取决于流体速度场 \mathbf{u} 和流道曲率引起的梯度变化，其数学表达式为：
其中，\mathbf{S} = \frac{1}{2}\left(\nabla \mathbf{u} + (\nabla \mathbf{u})^T\right) 为应变率张量 。在多孔介质中，错综复杂的孔隙分叉和汇合导致速度和剪切速率快速变化，从而引起孔隙空间内黏度的显著变化，如下图所示的 Carreau 流体流经不规则孔隙结构的示例。
流线图展示了 Carreau 流体的黏度分布：剪切速率较高的狭窄通道展现出较低的黏度，而剪切速率较低的较宽区域则维持较高的黏度，流线突出显示了孔隙结构对流体行为的影响。
这种动态黏度变化给非牛顿流体流动仿真带来挑战。为实现更大尺度的模拟，我们需要有效的升尺度方法，将孔隙尺度行为转化为宏观流动特征。
孔隙尺度模型的升尺度转换
对于大尺度应用，在孔隙尺度上模拟非牛顿流体的流动是行不通的，因此必须进行升尺度处理。其中一种方法涉及定义表观黏度 \mu_\textrm{app}，本质上是会产生与非牛顿流体相同压降的（恒定）黏度（参考文献 1）。
根据达西定律，表观黏度可以表示为 :
(1)
式中， \kappa 是多孔介质的渗透率， v 是速度， \nabla p 是压力梯度。这些数值可以通过测量获得，也可以像我们的示例一样，通过代表性体积单元（RVE）的孔隙尺度模型获得。（本征）渗透率是多孔基质的一种性质，仅取决于孔隙大小、形状和连通性，其计算方法可参考这篇博客及相应模型。
在已知结构渗透率的前提下，令 Carreau 流体通过该多孔结构，其表观黏度或有效黏度的定义如下
(2)
其中，参数 \mu_0 和 \mu_\textrm{inf} 分别为零剪切速率和无限剪切速率情况下的黏度，\lambda 为弛豫时间，n 为幂律指数。为了进一步证明升尺度方法的必要性，我们将相同的 Carreau 流体通过具有相同孔隙率和渗透率的不同多孔结构，并施加相同的压降。不出所料，剪切率受流道曲率的影响。由于不同结构的流道不同，剪切速率也不同，因此，我们测量或使用方程(1)中的方法进行计算的结果也不同。
在相同的边界条件下，具有相同孔隙率和渗透率的不同多孔结构展现出不同的表观黏度，这突出表明了流体和结构对剪切速率的影响。
在构建多孔结构的升尺度模型时，需要一个剪切速率的表达式，因为不能仅根据达西速度直接计算剪切速率。我们将这个表达式称为表观剪切速率\dot\gamma_\textrm{app} ，然后将其用于方程 (2)中。
对于不同应用场景，存在多种计算公式。一种常用的形式是:
(3)
式中，孔隙结构的影响用校正因子\alpha 表示。虽然可采用现有经验关系式确定 \alpha，但将 \alpha 与测量值或模拟数据进行拟合可以获得更精确的结果。
计算和验证校正因子
校正因子可通过 COMSOL® 系统内置的方法确定。首先，利用方程 (2) 计算表观剪切速率 \dot\gamma_\textrm{app}, 同时计算 \frac{|\mathbf{u}|}{\sqrt{\kappa\epsilon_\textrm{p}}} 项，我们称之为归一化速度。再根据方程 (3) 计算出 \alpha 。虽然通常可以假定 \alpha 在一定的压降范围内保持不变，但必须注意的是，它可能会因流动模式的变化（如湍流的发生）而变化。
本案例模拟了 50—5000 Pa/m 压降范围内的流动，COMSOL® 将计算结果输出到数据表，随后通过 最小二乘 拟合函数 对方程 (3) 进行参数拟合以确定 \alpha 值。
该方法不仅支持孔隙尺度模拟数据，也可直接采用实验测量数据进行拟合验证。
为了验证这种方法，我们可以根据 Brinkman 方程建立一个均质化模型，并对结果进行比较。如果结果一致，则证明了近似方法的准确性和可靠性。
Brinkman 方程 接口和 达西定律 接口现在支持使用表观剪切速率法（包括热效应）模拟非牛顿流体流动，也就是说您可以在 6 种常用的非牛顿本构关系中进行选择。此外，您还可以输入自定义校正因子或选择毛细管束方法，该方法提供了另一种适用于毛细管束状多孔介质的 \alpha 描述。下图是 COMSOL Multiphysics® 用户界面的截图，显示了流体属性设置窗口。
下图显示的结果表明，表观剪切速率法提供了准确可靠的预测。
孔隙尺度模型与均质化方法的表观黏度-归一化速度关系对比。
下一步
由于流体特性与孔隙结构之间存在复杂的相互作用，模拟多孔介质中的非牛顿流体流动具有独特的挑战性。将孔隙尺度模拟或测量结果与升尺度技术相结合，可以得到精确的宏观模型。多孔介质流动接口目前支持表观剪切速率法，提供了一种考虑动态黏度变化和结构影响的可靠方法。
想自己建立非牛顿流体在多孔介质中的流动均质化模型吗？COMSOL 案例库中提供了 MPH 文件和分步说明，欢迎下载：
在后续博客中，我们将以狭窄肺动脉中的血液流动为例，展示表观剪切速率法在真实生物力学场景中的应用。
参考文献
- N. Zamani, I. Bondino, R. Kaufmann, A. Skauge, “Computation of polymer in-situ rheology using direct numerical simulation”, Journal of Petroleum Science and Engineering, Volume 159, 2017; https://doi.org/10.1016/j.petrol.2017.09.011.
评论 (0)