土木建筑与环境工程  2014, Vol. 36 Issue (5): 44-49, 65   PDF    
多孔介质孔隙结构重构及水蒸气传递特性
郑佳宜, 戴晓丽, 陈振乾    
东南大学 能源与环境学院, 南京 210096
收稿日期：2014-01-06
基金项目：科技部国际科技合作技术交流专项（2011DFA60290）
作者简介：郑佳宜(1986-), 女, 博士生, 主要从事多孔介质传热传质研究, (E-mail)wwa_1986@163.com。
陈振乾(通信作者), 男, 教授, 博士生导师, (E-mail)zqchen@seu.edu.cn。
摘要：采用随机生长法分别构造了各向同性、各向异性的孔隙率相同，孔径分布和孔隙表面分形维数不同的多孔介质，探索了在多孔介质孔隙率相同的情况下，微观孔隙结构对孔隙内湿传递过程的影响。结果表明：决定多孔介质内部输运特性的孔隙结构参数，除孔隙率以外，孔径分布和孔隙表面分形维数等参数对多孔介质孔隙中水蒸气扩散过程的影响显著；孔隙率相同，孔径分布和孔隙表面分形维数不同的3种各向同性和3种各向异性多孔介质的水蒸气扩散浓度最大差别分别为30.6%和23.3%；各向同性多孔介质的小孔分布率越大，孔隙表面分形维数越大，说明孔隙的连通性和水蒸气在多孔介质孔隙中的扩散性能越好；各向异性多孔介质（水蒸气扩散的方向垂直于其生长率较大的方向）的大孔分布率越大，孔隙表面分形维数越小，水蒸气在多孔介质孔隙中的传递性能越好。
关键词：随机生长法    多孔介质    孔径分布    孔隙表面分形维数    湿传递    
Reconstruction of Pore Structures and Water Vapor Transfer Properties of Porous Medium
Zheng Jiayi, Dai Xiaoli, Chen Zhenqian    
School of Energy and Environment, Southeast University, Nanjing 210096, P. R. China
Received: 2014-01-06
Abstract: Isotropic and anisotropic porous medium with the same porosity, different pore size distributions and pore surface fractal dimensions were constructed by random growth method, and the impact of microscopic pore structure on moisture migration process in pores of porous medium was theoretically studied. The results show that the internal transport properties of porous medium are not only determinded by porosity, but also greatly influenced by pore size distribution and pore surface fractal dimension. The biggest water vapor diffusion concentration differences of three kinds of isotropic and three kinds of anisotropic porous medium with the same porosity, different pore size distributions and pore surface fractal dimensions are 30.6% and 23.3%, respectively. Better water vapor diffusion properties in pores of isotropic porous medium demands for bigger tiny pore size distribution rate and larger pore surface fractal dimension, which means the pore connectivity is better. While better water vapor diffusion properties in pores of anisotropic porous medium(water vapor diffusion direction is perpendicular to the direction of larger growth rate of anisotropic porous medium) demands for more big pore size distribution and smaller pore surface fractal dimension. The method lays a theoretical foundation for preparing indoor humidity control of building materials.
Key Words: random growth method    porous medium    pore size distribution    pore surface fractal dimension    water vapor miqration    

多孔介质传热传质现象普遍存在于自然界和人们的生产生活中，对国民经济发展、节约能源、改善人民生活水平具有重要的意义。土壤中的水、肥、污染物的迁移过程，地下岩层中地下水、天然气和石油等资源的开采，动植物的生命过程，纺织品、食品和建材的干燥，建筑围护结构使用寿命和对室内舒适性的影响等，都涉及到多孔介质中能量与物质的传递问题。

多孔介质传热传质理论的研究已有近百年的历史，现已形成比较完善的理论体系[1]。但是，多孔介质孔隙结构(包括孔隙率、孔径分布、孔隙连通性、孔容积、孔隙弯曲度等)的异常复杂和无规则性, 在过去成为学者们研究的瓶颈。毛细管束模型[2-4]、规则或不规则的网络模型[5-6]、分形模型[7-9]等简化模型被用于描述多孔介质的结构，但这与真实物质的孔隙形貌相差较大；随着数字图形信息技术的发展，人们通过CT扫描可以精确地获取三维多孔介质的孔隙形貌[10-13]，但此方法效率较低，花费较大，修正函数的关系式十分复杂；在多孔介质传热传质模型研究领域，大多数的数学模型仅以孔隙率作为多孔介质孔隙结构的表征[14-15]，而且对于非均质、各向异性多孔介质中通常有更多参数需要确定，从而带来求解困难并影响模型的广泛应用[16]。现有的模型鲜有从多孔介质物理模型本身出发来研究不同孔隙结构多孔介质中的热量或质量传递过程。

据此，本文采用随机生长法[17-18]构造出孔隙率相同的各向同性和各向异性多孔介质，研究了水蒸气在不同微观孔隙结构中的扩散规律。此方法简单高效，模拟结果直观，并且符合自然界物质本身的生长规律，从而为研究多孔介质中的热质传递过程奠定理论基础。

1 孔隙结构表征
1.1 随机生长法

随机生长法是一种描述多孔介质基质(或孔隙)在特定区域内生长的方法，通过在特定区域内随机分布基质生长核，调节预重构多孔介质孔隙结构的生长参数(生长概率、孔隙率等)，控制其孔隙形貌(各向同性或各向异性)，对其进行二维或三维重构。此方法能够准确地描述具有异常复杂和随机特性的多孔介质孔隙结构，解决了前人对多孔介质孔隙结构的评价参数都是建立在统计平均基础上的问题，符合自然界物质生长规律并能够准确反映多孔介质孔隙结构的几何特征。

1.2 多孔介质二维重构

图 1为多孔介质固体中心生长方向示意图，多孔介质的重构思路为：初始时，区域内全部为孔隙(ε=1)，以概率pro1在此区域内随机分布固相生长中心(i, j)，pro1必须小于预构造多孔介质中固相比例(1-ε)。然后，从固相生长中心以一定的概率向方向2(i, j+1)、4(i, j-1)、3(i-1, j)、1(i+1, j)、6(i-1, j+1)、5(i+1, j+1)、8(i+1, j-1)、7(i-1, j-1)生长。当向上、下、左、右的生长概率和向左上、右上、左下、右下的生长概率分别相等(p2=p4=p3=p1，p6=p5=p7=p8)时，多孔介质为各向同性；当向上、下和向左、右的生长概率分别相等又互相不同(p2=p4≠p3=p1)，且向左上、右上、左下、右下的生长概率相等(p6=p5=p7=p8)时，多孔介质为各向异性。按此规律生长直至孔隙率达到预设值，重构过程结束。采用以上方法重构了孔隙率均为0.45的3个各向同性多孔介质(图 2(a)~(c))和3个各向异性多孔介质(图 2(d)~(f))。

图 1 多孔介质固体中心生长方向示意图

图 2 各向同性((a)、(b)、(c))和各向异性((d)、(e)、(f))多孔介质

1.3 多孔介质的孔径分布

采用图像分析法统计出如图 2(a)~(f)所示的多孔介质二值图(黑色区域代表孔隙，白色区域代表固体)对应的孔径分布如图 3(a)~(f)所示，即通过如图 3所示的孔径区间与该区间内孔隙的累积面积占总孔隙面积的比例作柱状图。可以看出，虽然图 2中构造出的3个各向同性和3个各向异性多孔介质的孔隙率均为0.45，但其孔径分布各不相同。

图 3 各向同性和各向异性多孔介质分别对应的孔径分布

1.4 孔隙表面分形维数

采用计盒维数法[19-20]确定多孔介质的表面分形维数df，将构造的多孔介质二值图逐次栅格化，得到一系列不同边长(ε=1, 1/2, 1/22, 1/23, 1/24, …)的网格覆盖在分形多孔介质图形上，统计每一网格边长值对应的非空格子数(N(r))。将每次划分所得的非空格子数与对应的网格边长在双对数坐标下进行线性拟合所得拟合直线斜率的绝对值即为多孔介质的表面分形维数，图 4(a)~(f)为图 2(a)~(f)所示多孔介质对应的孔隙表面分形维数图。

图 4 各向同性和各向异性多孔介质分别对应的表面分形维数

2 数学模型

主要研究孔隙率不变时，不同的微观孔隙结构(孔径分布和孔隙分形维数)对多孔介质孔隙中水蒸气传递性能的影响，因此对模型进行如下简化：1)湿传递过程为等温过程；2)湿传递过程中无相变传质，仅为单纯的水蒸气在多孔介质孔隙中扩散；3)无化学反应，固体骨架和气相均不可压缩。基于以上假设，建立如式(1)所示的控制方程。图 5为多孔介质孔隙中水蒸气扩散过程示意图，左边界为定浓度边界，右边界为自然对流边界，上、下边界为绝湿边界。

图 5 多孔介质中水蒸气扩散示意图

多孔介质孔隙中水蒸气扩散方程为

$ \frac{{\partial c}}{{\partial t}} = \frac{\partial }{{\partial x}}\left( {{D_{\rm{M}}}\frac{{\partial c}}{{\partial x}}} \right) + \frac{\partial }{{\partial y}}\left( {{D_{\rm{M}}}\frac{{\partial c}}{{\partial y}}} \right) $ (1)

初始条件：

$ t = 0, c = {c_0} $ (2)

边界条件：

$ t > 0, c = {c_{{\rm{bo}}}}\;\;0 \le y < H, x = 0 $ (3)
$ - {D_{\rm{M}}}\frac{{\partial c}}{{\partial y}} = 0\;\;\;\;0 \le x < L, y = 0和y = H $ (4)
$ - {D_{\rm{M}}}\frac{{\partial c}}{{\partial x}} = {k_{\rm{c}}}\left( {{c_{\rm{b}}} - c} \right)\;\;\;\;0 \le y < H, \mathit{x} = \mathit{L} $ (5)
2.1 扩散系数确定

为了确定多孔介质孔隙中水蒸气扩散系数，首先需要判断多孔介质孔隙中水蒸气分子的扩散类型，为此需要计算Knudsen数，Knudsen数被定义为水蒸气分子平均自由程λ与多孔介质的孔径d之比。多孔介质孔隙中水蒸气分子的扩散类型包括分子扩散(Kn≤0.1)、Knudson扩散(Kn≥10)、过渡扩散(0.1≤Kn≤10)。

分子运动的平均自由程表示每两次连续碰撞之间，一个分子自由运动所走过的平均路程。根据分子动力学理论，平均自由程可表示为

$ \lambda = \frac{{3.2\mu }}{P}{\left( {\frac{{{R_{\rm{s}}}T}}{{2{\rm{ \mathsf{ π} }}M}}} \right)^{0.5}} $ (6)

根据气体动理学，平均自由程可表示为

$ \lambda = \frac{{kT}}{{\sqrt 2 {\rm{ \mathsf{ π} }}{\mathit{d}^{\rm{2}}}P}} $ (7)

根据扩散传质理论，当环境温度T=298 K，P=101 325 Pa时，计算得出Kn≤0.01，这说明水蒸气在多孔介质孔隙中的运动为符合fick定律的分子扩散。假设多孔介质的孔隙为圆柱形，多孔介质孔隙中的水蒸气分子扩散系数为

$ {D_{\rm{M}}} = \frac{\varepsilon }{\tau }{D_{{\rm{va}}}} $ (8)

其中，水蒸气在空气中的扩散系数[21]为

$ {D_{{\rm{va}}}} = \frac{{0.016\;7{T^{1.5}}}}{{P{{\left( {V_\rm{V}^{0.5} + V_{^\rm{a}}^{0.5}} \right)}^2}}} \cdot {\left( {\frac{1}{{{M_{\rm{V}}}}} + \frac{1}{{{M_{\rm{a}}}}}} \right)^{0.5}} \cdot \alpha $ (9)
$ \alpha = 1 - 2.25 \times {10^{ - 2}}\left( {\mu _{\rm{V}}^{0.5} + \mu _{\rm{a}}^2} \right) $
2.2 数值计算

考虑到多孔介质孔隙结构的异常复杂和无规则性，应用自由三角形网格对多孔介质的孔隙进行网格划分，对局部区域进行网格加密，同时采用不同的网格尺寸检测网格的独立性，从而保证所得到的数值解都是网格独立的解。

采用有限元法对多孔介质孔隙中水蒸气传递的控制方程进行数值求解，在求解区域内迭代求解水蒸气传递方程，在每个时间节点，当相邻两个迭代步长之间，各节点的水蒸气浓度残差值均小于1×10-6，可以认为该时间节点上迭代收敛，从而得到求解区域内的水蒸气浓度分布。

3 结果与讨论

模型中整个计算区域的初始水蒸气浓度为0.001 5 mol/m3，温度恒定，左边界为浓度边界0.15 mol/m3，右边界为通量边界，也叫自然对流边界，上、下边界为绝湿边界。图 6(a)、(b)分别为各向同性和各向异性多孔介质右上角点水蒸气浓度随时间变化图。从中可以看出，各向同性多孔介质右上角点水蒸气浓度达到稳定的时间约为50 s，而各向异性多孔介质右上角点水蒸气浓度达到稳定的时间约为200 s；浓度曲线趋于稳定时，各向同性多孔介质右上角点水蒸气浓度值普遍高于各向异性多孔介质右上角点水蒸气浓度值，这是因为各向同性多孔介质基质的生长率在正4向(p2=p4=p3=p1)和对角4向(p6=p5=p7=p8)分别相等，所以水蒸气扩散速率和水蒸气扩散浓度达到稳定分布的时间相对稳定，而各向异性多孔介质中水蒸气浓度梯度方向垂直于其生长率较大的方向，导致水蒸气穿过多孔介质孔隙时的阻力增加和浓度减小。这种结构上的差异同时也导致各向异性多孔介质右上角点水蒸气浓度达到稳定的时间远大于各向同性多孔介质右上角点水蒸气浓度达到稳定所需的时间。各向同性多孔介质的右上角点水蒸气浓度随时间变化曲线是不同的，小孔分布率越大(如图 3(c)的20~30 μm的孔径分布最密集)，孔隙表面分形维数越大(如图 4(a)~(c))所示)，右上角点水蒸气浓度值越大(如图 6(a)中的3#所示)，说明孔隙的连通性和水蒸气在多孔介质孔隙中的扩散性能越好；各向异性多孔介质的右上角点水蒸气浓度随大孔分布率增大而增大(图 3(f)的80 μm的孔径分布率达70%)，随孔隙表面分形维数的增加而减小(如图 4(d)~(f))所示)，这是因为水蒸气扩散方向垂直于多孔介质生长率较大的方向，孔隙表面分形维数增加将导致阻碍水蒸气传递的固体生长率增加，从而导致穿过孔隙的水蒸气浓度降低，所以大孔分布率越大，各向异性多孔介质的水蒸气传递性能越好(如图 6(b)中的3#所示)。

图 6 多孔介质右上角点的水蒸气浓度变化

图 7(a)~(c)为运用COMSOL Multiphysics 4.1模拟出的3个孔隙率为0.45的各向同性多孔介质孔隙中水蒸气浓度分布图，图 7(d)~(f)为运用COMSOL Multiphysics 4.1模拟出的3个孔隙率为0.45的各向异性多孔介质孔隙中水蒸气浓度分布图。从图 7(a)~(c)中可以直观地看出，图 7(c)的孔隙结构中水蒸气浓度分布最均匀，只有右边界处出现少量的浓度阶跃，最低浓度分布值为0.087 mol/m3；图 7(b)的孔隙结构中水蒸气浓度在多孔介质x轴断面约1/3处出现明显的阶跃现象，说明此处孔径较细，甚至孔隙不连通，此多孔介质的孔隙连通性最差，最低水蒸气浓度值为0.041 mol/m3；由此可以计算出各向同性多孔介质右边界的水蒸气扩散浓度相对于左边界恒定浓度的最大差别为(0.087 mol/m3-0.041 mol/m3)/ 0.15 mol/m3=30.6%，而图 7(a)的孔隙结构中水蒸气浓度分布介于图 7(b)、(c)之间，在约1/4的右上角部分水蒸气浓度出现了阶跃现象，最低浓度分布值为0.052 mol/m3，说明在此处孔隙的连通性变差，出现了孔径较小或者不连通的孔隙。图 7(d)~(f)所示的各向异性多孔介质孔隙中水蒸气浓度分布普遍没有各向同性多孔介孔隙质中水蒸气浓度分布均匀，这是由于各向异性多孔介质孔隙中水蒸气扩散方向垂直于其生长率较大的方向，而各向异性多孔介质中生长率较大方向的孔隙连通性好于生长率较小方向的孔隙连通性。沿水蒸气扩散方向，各向异性多孔介质中出现水蒸气浓度阶跃现象的位置可以初步判断水蒸气在多孔介质中扩散情况，图 7(e)出现水蒸气浓度阶跃的位置在约4×10-7 m处，图 7(d)出现水蒸气浓度阶跃的位置在约5.5×10-7 m处，图 7(f)出现水蒸气浓度阶跃的位置在约5.5×10-7 m处，所以图 7(e)最早出现水蒸气浓度阶跃，水蒸气浓度值为最低0.006 mol/m3；从模拟的数值结果可以进一步判断图 7(d)、(f)的水蒸气浓度分布情况，图 7(d)、(f)的最低水蒸气浓度值分别为为0.023、0.041 mol/m3，说明图 7(f)的孔隙连通性最好，图 7(d)次之，图 7(e)最差。由此可以计算出各向异性多孔介质右边界的水蒸气扩散浓度相对于左边界恒定浓度的最大差别为(0.041 mol/m3-0.006 mol/m3)/ 0.15 mol/m3=23.3%。

图 7 各向同性与各向异性多孔介质的水蒸气浓度分布图

4 结论

1) 随机生长法能够有效重构孔隙形貌呈各向同性或各向异性的多孔介质，从而直观地观察多孔介质中流体流动或传热传质过程，实现实验难以达到的效果。

2) 各向同性多孔介质的小孔分布率和孔隙表面分形维数越大，说明孔隙连通性和水蒸气在多孔介质孔隙中的扩散性能越好；各向异性多孔介质(水蒸气扩散方向垂直于多孔介质生长率较大方向)的大孔分布率越大，孔隙表面分形维数越小时，说明水蒸气在多孔介质孔隙中的传递性越好。

3) 孔隙率为0.45，孔径分布和孔隙表面分形维数不同的3种各向同性和3种各向异性多孔介质的水蒸气扩散浓度最大差别分别为30.6%和23.3%。

符号说明

c：水蒸气浓度，mol·m-3；

M：摩尔质量，1×10-3 kg·mol-1；

c0：水蒸气的初始浓度，mol·m-3;

p：压力，Pa;

cbo：边界水蒸气浓度，mol·m-3;

Rs：气体常数，N·m·mol-1·K-1;

DM：多孔介质孔隙中水蒸气的分子扩散系数，m2·s-1;

T：热力学温度，K;

d：分子有效直径，m;

x、y：分别为横坐标和纵坐标，m;

H：多孔介质试样的高度，m;

α：考虑分子间作用力对气体扩散系数影响的修正系数;

kc：水蒸气的质扩散系数，m2·s-1;

μ：动力黏度，Pa·s;

k：玻尔兹曼常数，1.38×10-23J·K-1;

μv、μa：分别为水蒸气和空气的偶极距，C·m;

L：多孔介质试样的长度，m。

下角标

v：水蒸气;

a：空气。

参考文献
[1] Incropera F P, Lavine A S, de Witt D P.Fundamentals of heat and mass transfer[M].John Wiley & Sons Incorporated, 2011.
[2] Haring R E, Greenkorn R A. A statistical model of a porous medium with nonuniform pores[J]. AICHE Journal, 1970, 16(3): 477–483. DOI:10.1002/(ISSN)1547-5905
[3] Wei Y, Durian D J. Effect of hydrogel particle additives on water-accessible pore structure of sandy soils:A custom pressure plate apparatus and capillary bundle model[J]. Physical Review E, 2013, 87(5): 053013. DOI:10.1103/PhysRevE.87.053013
[4] Jackson M D. Multiphase electro-kinetic coupling:Insights into the impact of fluid and charge distribution at the pore scale from a bundle of capillary tubes model[J]. Journal of Geophysical Research, 2010, 115(B7): B07206.
[5] Sahimi M.Flow and transport in porous media and fractured rock[M].Wiley.com, 2012.
[6] Blunt M J, Jackson M D, Piri M, et al. Detailed physics, predictive capabilities and macroscopic consequences for pore-network models of multiphase flow[J]. Advances in Water Resources, 2002, 25(8): 1069–1089.
[7] Xu P, Yu B. Developing a new form of permeability and Kozeny-Carman constant for homogeneous porous media by means of fractal geometry[J]. Advances in Water Resources, 2008, 31(1): 74–81. DOI:10.1016/j.advwatres.2007.06.003
[8] Yun M, Yu B, Cai J. A fractal model for the starting pressure gradient for Bingham fluids in porous media[J]. International Journal of Heat and Mass Transfer, 2008, 51(5): 1402–1408.
[9] Zhang L Z. A fractal model for gas permeation through porous membranes[J]. International Journal of Heat and Mass Transfer, 2008, 51(21): 5288–5295.
[10] Narsilio G A, Buzzi O, Fityus S, et al. Upscaling of navier-stokes equations in porous media:Theoretical, numerical and experimental approach[J]. Computers and Geotechnics, 2009, 36(7): 1200–1206. DOI:10.1016/j.compgeo.2009.05.006
[11] Du D X, Beni A N, Farajzadeh R, et al. Effect of water solubility on carbon dioxide foam flow in porous media:an X-ray computed tomography study[J]. Industrial & Engineering Chemistry Research, 2008, 47(16): 6298–6306.
[12] Zhang T, Lu D, Li D.Porous media reconstruction using a cross-section image and multiple-point geostatistics[C]//Advanced Computer Control, 2009.ICACC'09. International Conference on.IEEE, 2009:24-29.
[13] Zhou N, Matsumoto T, Hosokawa T, et al. Pore-scale visualization of gas trapping in porous media by X-ray CT scanning[J]. Flow Measurement and Instrumentation, 2010, 21(3): 262–267. DOI:10.1016/j.flowmeasinst.2010.05.002
[14] Bear J.Dynamics of fluids in porous media[M].Dover Publications.com, 2013.
[15] Woods J, Pellegrino J, Burch J. Generalized guidance for considering pore-size distribution in membrane distillation[J]. Journal of Membrane Science, 2011, 368(1): 124–133.
[16] Bear J, Bachmat Y.Introduction to modeling of transport phenomena in porous media[M].Springer, 1990.
[17] Wang M, Pan N. Numerical analyses of effective dielectric constant of multiphase microporous media[J]. Journal of Applied Physics, 2007, 101(11): 114102. DOI:10.1063/1.2743738
[18] Wang M, Pan N. Predictions of effective physical properties of complex multiphase materials[J]. Materials Science and Engineering:R:Reports, 2008, 63(1): 1–30. DOI:10.1016/j.mser.2008.07.001
[19] Ge M, Lin Q. Realizing the box-counting method for calculating fractal dimension of urban form based on remote sensing image[J]. Geo-spatial Information Science, 2009, 12(4): 265–270. DOI:10.1007/s11806-009-0096-1
[20] Li J, Du Q, Sun C. An improved box-counting method for image fractal dimension estimation[J]. Pattern Recognition, 2009, 42(11): 2460–2469. DOI:10.1016/j.patcog.2009.03.001
[21] 陈昭宜, 刘蕊梅. 一个新的气体扩散系数计算公式[J]. 科学通报, 1989, 11: 29.
Chen Z Y, Liu R M. A new formula about gas diffusion coefficient[J]. Chinese Science Bulletin, 1989, 11: 29. (in Chinese)
X
    多孔介质孔隙结构重构及水蒸气传递特性
    郑佳宜, 戴晓丽, 陈振乾