MBR单根中空纤维膜丝的CFD模拟仿真与研究

CFD Simulation and Research of MBR Single Hollow Fiber Membrane

摘要: 在对MBR膜生物反应器的模拟仿真中，我们曾使用计算流体动力学(CFD)的相关软件建立MBR中空纤维膜的二维模型，并成功模拟出流体在中空纤维膜区域的流动现象。为了后期提高中空纤维膜的模拟效果，我们可以考虑先将单根中空纤维膜丝放进膜组件中，并建立三维模型进行模拟，计算和观察，为扩展到多根中空纤维膜丝的研究做好前期的工作。在模拟中，由于中空纤维膜管壁小孔的数量规模达十亿数量级，在建模阶段画出大量的小孔并不现实。所以，在模拟时，我们引入多孔介质的思想，将中空纤维膜的管壁定义为多孔介质区域，从而完成对中空纤维膜管壁小孔的模拟。在之后进行的CFD数值求解的计算中，我们采用比较精确的有限体积法求解离散方程，经迭代100次后，我们可以从残差曲线图中看到所有值均低于阈值，呈收敛状态。模拟结果表明，通过CFD相关软件对单根中空纤维膜丝的模拟[1]，基本符合流体的实际流动现象。

Abstract: By simulation analysis in MBR membrane bioreactors research, we have used computational fluid dynamics (CFD) related software to create a two-dimensional model and simulated the flow of fluid in the hollow fiber membrane. In order to improve the simulation effect and enhance precision in calculation, we can consider firstly putting a single fiber into the membrane module and establish a three-dimensional model for simulation, so as to do a good job for the study of the expansion of the 3000 hollow fiber membrane. It is not realistic to draw large numbers of holes in a simulation because the number of small holes in the hollow fiber membrane is the orders of magnitude of billion. Therefore, we use model of porous media to define the wall of hollow fiber membrane as a porous zone. In calculating the numerical solution of CFD, we use the more exact finite volume method (FVM) to solve the discrete equations. After 100 iterations, we can see from the residuals curve that all values are below the threshold and converge. The result shows that the simulation of single hollow fiber membrane by CFD software [1] is basically consistent with the actual flow phenomenon of the fluid.

1. 引言

MBR膜生物反应器是一种将膜分离单元和生物反应器技术相结合的新型污水处理技术。膜污染问题一直是阻碍其发展的重要因素，而使用中空纤维膜的MBR膜生物反应器，可以在更高的活性污泥浓度下保持高通量的稳定运行 [2] [3] 。近几年来，不少人使用过计算流体动力学(CFD)软件对中空纤维膜进行模拟仿真 [4] ，例如，使用Gambit建立二维模型，并使用有限差分法计算方程结果，模拟出中空纤维膜的速度场和浓度场 [5] ，并将多相流模型应用到中空纤维膜上，通过气、液、固三相流观察曝气量，以及液、固两相的流动情况。

随着CFD技术的发展，本文借助目前比较流行的专业CFD软件ICEM CFD建立三维模型，并在Fluent中使用多孔介质模型和有限体积法，求解离散方程。最后通过CFD-Post软件，生成压力云图和速度矢量图，供我们有效地观察和分析流动计算结果，掌握中空纤维膜周围流动的规律。为扩展到多根中空纤维膜丝的研究做好前期的工作 [1] 。本文在建立模型时，膜组件两侧水槽忽略不计，将中空纤维膜丝的两端直接接在出口处，这样的模型可能会在模拟结果中产生误差。同时，污水在进入膜组件之前，已经默认通过了预处理。经模拟结果表明，在忽略误差的前提下，通过CFD相关软件对单根中空纤维膜丝的模拟，基本符合流体的实际流动现象。

2. 单根中空纤维膜丝在CFD中的模拟

我们曾经使用过CFD对中空纤维膜进行模拟，但是在实验中也存在着较大的问题。在模拟的过程中，建立的模型是二维模型，并且在运用求解器进行CFD数值计算时，采用的是有限差分法，而这种方法仅当网格极其细密时，离散方程才满足积分守恒。

针对这些问题，我们首先将模型的维度从二维模型上升到三维模型。并采用CFD数值解法中的有限体积法 [6] 。由于有限体积法即使在粗网格情况下，也能够显示出准确的积分守恒，所以本实验在求解计算时用有限体积法来代替有限差分法，这样可以弥补有限差分法的不足，从而得到更精确的数值 [7] [8] 。

本实验模拟的中空纤维膜是天津膜天膜公司的UOF-4型膜组件中的8英寸膜，膜直径为203 mm，膜组件内浇注约3000~3500根膜丝。

2.1. 建立模型

我们借助CFD的前处理器ICEM CFD软件建立一个三维模型，该模型由一个膜组件和一根中空纤维膜组成。由于研究的主要对象是中空纤维膜丝，这里需要借助膜组件作为模拟流体流动时的壁面，所以为了观察的时候更方便，在不影响对中空纤维膜流体流动研究的基础上，只保留了膜组件的入口，出口和壁面，而膜组件两边出口处附近的两个水槽和出口的管道可忽略不计，让水直接流出。由于膜丝的过滤形式是外压式，所以水从入口处进入后，膜组件向中空纤维膜丝施加压力，水渗透过中空纤维膜并从两端的出口处流出。建立的三维模型如图1所示：

Figure 1. 3D model of membrane module and single membrane filament

图1. 膜组件和单根膜丝组成的三维模型

由于上图中的三维模型所在计算区域较为均匀，数据结构简单，同时考虑到网格的生成速度和质量，我们将模型划分成结构化网格，划分的网格如图2所示：

Figure 2. Structured grid

图2. 结构化网格

结构化网格划分成功后，为了保证网格导入到fluent中求解方程时能够收敛，检查网格的相关参数Quality，如图3所示：

Figure 3. Quality

图3. 网格质量图

从图3中可以看出，网格的参数Quality趋向于1，而Quality的质量已达到0.4，网格质量已经基本符合要求。

2.2. Fluent求解

将画好的结构化网格导入到Fluent中，进行计算求解，求解步骤如下：

1) 查看模型的网格信息，并检查网格的最小体积和最小面积是否为负数。

2) 由于ICEM CFD软件在建模时是无量纲的，所以需要设置模型的单位。

在CFD中，判断流动是层流还是湍流，由雷诺数Re的值决定，Re是一种用来表征流体流动情况的无量纲数，具体公式如下：

$\mathrm{Re}=du\rho /\mu$ (1)

其中，d为8寸膜的直径，u为流体的流速， $\rho$ 为待处理废水的密度， $\mu$ 为待处理废水的黏度。流速u的计算公式如下：

$u=A/S$ (2)

其中，A为流体的流量，S为膜面积。

当雷诺数Re的值小于临界雷诺数2000时，即Re < 2000，流动是层流，否则为湍流。因此，代入参数A = 2 m³，S = 40 m²，d = 203 mm，水的密度 $\rho$ 为 = 1000 kg/m³，废水的黏度 $\mu$ 为 = 0.8949 × 10⁻³ kg/m-s，求得Re = 3.15，由于Re的值远小于2000，所以判断流动为层流。

由于实验中不会涉及到温度的变化，所以这里不考虑能量守恒方程，只并打开动量守恒的设置。

3) 选择材料为水，并输入密度 $\rho$ = 1000 kg/m³，黏度 $\mu$ = 0.8949 × 10⁻³ kg/m-s。

4) 设置多孔介质模型和参数，由于在步骤(3)中已判断流动为层流，所以可以将惯性阻力系数忽略不计，只需要计算黏性阻力系数。计算黏性阻力系数的相关公式为：

$\Delta p=\mu \Delta n/\alpha$ (3)

其中，p为压降，n为中空纤维膜的厚度， $\mu$ 为黏性阻力系数。代入参数p = 0.15 MP，n = 0.5 mm， $\mu$ = 0.8949 × 10⁻³ kg/m-s，求得 $\alpha$ = 167,616,493,462.9568。

5) 设置条件边界 [9] 。入口为速度入口，速度为步骤(3)中得到的流速，设置出口为压力出口。

6) 根据实际工作经验得到，如果速度是100 m/s，那么临界值一般为Fluent中的默认值10⁻³，本课题中的流速为 $1.39\times {10}^{-5}\text{m}/\text{s}$ ，所以这里将临界值设为10⁻⁹。

7) Fluent默认为有限体积法进行计算，所以使用默认值，参数初始化使用入口的速度值进行初始化，迭代次数为1000次，最后生成残差曲线图如图4所示：

Figure 4. Residual curve

图4. 残差曲线图

从图4中可以看到，计算约至780步左右的时候，所有值均小于阈值，残差达到收敛标准，计算结束。

3. 模拟与结果分析

为了能够更直观有效的观察和分析流动计算结果，我们将Fluent的计算结果导入到CFD-POST中，并分别创建截面z-01 (xy平面z = 1 mm处)，z-05，z-06，z-18 (xy平面z = −8 mm处)，分别显示这4个截面的压力云图，结果如图5所示：

Figure 5. Pressure cloud chart

图5. 残差曲线图

如图5所示，左边的颜色条表示压力值从下到上依次增加，最小值为0.201082 Pa，最大值为0.201087 Pa。在这里，压力的变化主要由膜污染程度决定，水流速度在一定范围内，影响可以忽略，最大值到最小值相差5也在正常预估范围内。

另外，创建截面y-0，z-6，并生成速度矢量图，结果如图6所示：

Figure 6. Velocity vector chart

图6. 速度矢量图

从图6中，我们可以看出，水在出口的流速明显小于入口处的流速，这是因为在边界条件中设置了压力出口的边界条件，并且在通过多孔介质区域时形成了阻力，影响了流体继续往出口方向流动。再结合图5的压力云图可以看出，出口处的压力要大于其他区域的压力，这是因为水流在出口处受到压力出口边界条件的影响，将水不断的向多孔介质区域挤压，所以在出口处，水的流量逐渐增大。

经过此次的模拟结果，我们可以分析出，该模拟基本实现了水流在膜组件中通过单根中空纤维膜丝的流动现象。

4. 结论与展望

为了能够提高流体在MBR中空纤维膜区域的模拟效果，本文先以单根中空纤维膜丝为基础 [10] ，并借助目前比较流行的专业CFD软件进行模拟和计算。在模拟中，我们引入了多孔介质模型，将中空纤维膜的管壁定义为多孔介质区域。而且，我们引入CFD数值解法中的有限体积法，使离散方程能够在粗网格情况下也能显示出准确的积分守恒。结合以上两点并建立三维模型，使模拟结果更接近流体在MBR中空纤维膜中的流动现象。模拟结果表明，通过CFD相关软件对单根中空纤维膜丝的模拟图像是符合实际流体的流动现象的。

在接下来的研究中，我们将逐次在膜组件中递增中空纤维膜丝的数量，分别观察每次添加膜丝后流体的流动现象，直至实现完整的MBR在实际中的运行情况。

基金项目

国家自然科学基金项目(51378350)；国家自然科学基金青年科学基金资助项目(50808130)；国家自然科学基金青年科学基金项目(21506159)。

参考文献