多孔地质材料中的反应扩散与孔隙结构演化：Ⅰ. 理论分析与数值模拟
姜元勇1，徐曾和1，李艳松2，李勤美3
(1. 东北大学 资源与土木工程学院，辽宁 沈阳，110819；
2. 河北工程大学 土木工程学院，河北 邯郸，056038；
3. 云南冶金集团股份有限公司，云南 昆明，650224)
摘要：考虑到多孔地质材料中的反应、扩散和孔隙结构演化同时发生在微观和细观2个不同的层次上，并且3个过程之间存在明显的相互作用，是耦合在一起的，研究反应固体和惰性固体共存情况下，多孔地质材料中的反应扩散与孔隙结构演化过程。同时，考虑到能够在实验室内进行实时测量的浮重力可以间接地反映多孔介质中的反应扩散和孔隙结构演化进程，研究多孔地质材料中的反应扩散和孔隙结构演化与浮重的关系。首先，在颗粒尺度上，给出同时考虑界面化学反应和边界层外扩散联合控制的颗粒动力学模型，描述固体反应颗粒的反应转化过程；接着，将颗粒动力学推广到包含足够多颗粒的多孔介质微元体上，得到描述多孔介质中反应扩散和孔隙结构演化的数学模型；最后，运用有限体积法对数学模型进行数值求解。研究结果表明：化学反应、溶质扩散和孔隙结构演化相互影响，密不可分；缓慢的反应和扩散过程，可以引起孔隙结构的显著改变，不容忽视；化学作用导致的孔隙结构改变远比力学作用的大；反应扩散与孔隙结构演化过程与浮重之间存在一一对应关系。本文所建立的数学模型可以反映非惰性多孔介质中的反应扩散和孔隙结构演化进程，对于进一步研究多孔介质中的物理化学流体动力学过程具有一定借鉴意义。
关键词：非惰性多孔介质；反应扩散；孔隙结构演化；微元体动力学；浮重
中图分类号：TD925 文献标志码：A 文章编号：1672-7207(2014)10-3561-10
Reaction-diffusion and pore evolution in porous geological material: Ⅰ. Theoretical analysis and numerical simulation
JIANG Yuanyong1, XU Zenghe1, LI Yansong2, LI Qinmei3
(1. College of Resources and Civil Engineering, Northeastern University, Shenyang 110819, China;
2. School of Civil Engineering, Hebei University of Engineering, Handan 056038, China;
3. Yunnan Metallurgical Group Company Limited, Kunming 650224, China)
Abstract: Considering that the chemical reaction, solute diffusion and pore structure evolution in porous geological material proceed simultaneously on micro-scale and meso-scale, and has interaction obviously, the three processes in the presence of inert solid and reactive solid coexistence were researched. In addition, considering that there are some one-one correspondences between microscopic reaction-diffusion and macroscopic buoyed weight, the relationship between reaction-diffusion and the buoyed weight was studied. A pellet kinetic model, which takes into consideration the resistances including the interface chemical reaction and the diffusion in boundary layer, was developed and then extended to representative element volume of porous media, solved by volume control method finally. By cases computation, some conclusions are as follows: The chemical reaction, solute diffusion and pore structure evolution in porous geological material are coupled and have interactions; Although the chemical reaction and solute diffusion in porous media take place slowly, the pore structure evolution due to them is obvious and can not be ignored; The change of porosity is more notable for chemical effect than single mechanical behavior; there are some one-one correspondences between microscopic reaction-diffusion and macroscopic buoyed weight. This model can be used to research the reaction-diffusion process in reactive porous media and has some reference value for further study in physicochemical hydrodynamics in porous media.
Key words: reactive porous media; reaction-diffusion; pore evolution; volume element kinetics; buoyed weight
多孔地质材料中的反应扩散与孔隙结构演化现象普遍存在，如流体污染物在地层中的聚集过程，煤层瓦斯的运移与抽放过程，化学法开采有用矿物的过程，喀斯特地貌的形成过程等。由于多孔地质材料形成过程的复杂性，以及其内部组分的多样性，使得研究其中的反应扩散和孔隙结构演化过程变得极其复杂。高红波等[1-4]针对岩盐矿床的水溶开采过程，对钙芒硝试样的溶解过程进行了大量的实验研究，考察了钙芒硝矿物溶解过程的物理化学特性，并采用CT技术获得了孔隙结构的图像。尹升华等[5-9]针对金属矿物的堆浸过程，采用柱浸实验研究了矿物浸出过程中的物理化学特性，采用CT技术观测了堆积矿块之间空隙的变化情况。Dadvar等[10-11]研究了多孔介质中存在化学反应情况下的孔扩散进程，但是他们忽略了化学反应所引起的孔隙结构演化。Romano等[12-13]研究了核工程中，核燃料燃烧和核废料处理所引起的多孔介质孔隙结构改变，但没有给出描述孔隙结构改变的理论表达形式。鉴于多孔地质材料中的反应扩散与孔隙结构演化过程是同时发生和相互影响的，并且此过程同时发生在微观和细观2个水平上，很难通过试验仪器直接测量。受到热重分析技术(TGA)的启发，考虑到随着多孔介质内部的反应扩散与孔隙结构演化的进行，置于流体中的多孔介质的浮重将发生变化，而浮重可以通过电子天平进行实时的准确测量，本文作者拟研究多孔地质材料中的反应扩散与孔隙结构演化过程，以及这一传质和反应耦合过程与多孔介质浮重的关系，以便为进一步深入研究多孔地质材料中的反应扩散和孔隙结构演化进程提供参考。
1 理论分析
考察如图1所示的圆片状多孔地质材料(图1中，R为圆片半径，m；H为圆片厚度，m)，假设圆片状多孔地质材料是由反应固体颗粒和惰性固体颗粒所构成，颗粒近似为球形，当反应流体扩散进入多孔地质材料时，其中的反应固体颗粒发生化学反应，反应固体与产物固体之间的体积差异将引起多孔介质孔隙结构的改变。
图1 圆片状多孔地质材料
Fig. 1 Disc porous geological material
图2所示为固体反应物颗粒示意图。图2中，rc为未反应颗粒半径，m；rs为颗粒初始半径，m；Cb为颗粒周围的流体产物质量浓度，kg/m3；Cc为未反应核表面的流体产物质量浓度，kg/m3。采用缩核模型描述单个反应颗粒的反应进程，考虑界面化学反应和产物流体外扩散联合控制，化学反应过程中无产物固体生成，即aA(l)+bB(s)=cC(l)型反应。
图2 固体反应物颗粒示意图
Fig. 2 Schematic diagram of reactive solid reactant pellet
未反应核表面的化学反应速率可以表示为
(1)
式中：kr为化学反应速率常数，m/s；Cs为流体产物的饱和质量浓度，kg/m3。
流体产物通过颗粒表面边界层向圆片孔隙内流体的传质速率为
(2)
式中：kg为外扩散传质系数，m/s。
由质量传递过程的拟稳态假设
(3)
联立式(1)~(3)，得到以产物流体C的生成表示的颗粒综合反应速率
(4)
由化学计量关系，固体反应物B的消耗速率和流体产物C的生成速率之间满足
(5)
式中：b和c为化学计量系数。
将颗粒动力学推广到含有足够多颗粒的微元体上，得到以反应固体B的消耗表示的微元体综合反应速率
(6)
式中：
为多孔介质的初始孔隙率；p为反应固体在多孔介质固相中的初始体积分数，可得
(7)
式中：
和
分别为反应固体B和惰性固体I的质量密度，kg/m3；XB0为初始时刻反应固体B占多孔介质固体骨架的质量分数。
根据几何关系，固体颗粒的反应速率还可以表示为
(8)
式中：t为时间，s。
联立式(4)，(5)和(8)，可以得到
(9)
考虑到t=0时，rc=rs，因此积分式(9)，可以得到
(10)
图3所示为流体中圆片状多孔介质的受力分析。图3中，Gs为固体骨架的重量，N；Gl为孔隙流体的重量，N；F1为浮力，N；F2为浮重，N。由图3可见：当多孔介质被反应流体包围时，对多孔圆片试件而言，假设圆片中的扩散过程仅发生在x方向，由质量守恒，多孔介质中的产物流体浓度满足
(11)
初始条件为
；
；
(12)
边界条件为
；
(13)
式中：
为多孔介质的孔隙率；Db为流体产物的孔扩散系数，m2/s；Cm为多孔介质周围流体产物的质量浓度，kg/m3。
图3 流体中圆片状多孔介质的受力分析
Fig. 3 Stress state of disc porous media in fluid
反应固体转化率的分布为
(14)
多孔介质中的孔隙率分布为
(15)
多孔介质中，反应固体的总转化率为
(16)
多孔介质的总孔隙率为
(17)
多孔介质的固体骨架质量为
(18)
式中：m0为初始时刻多孔介质的固体骨架质量，kg。
多孔介质孔隙中的流体质量为
(19)
式中：S为多孔介质的横截面积，m2，
；
为孔隙内流体的密度，kg/m3，针对不同的流体体系，可以得到流体质量密度与浓度之间的关系式
。
多孔介质所受的浮力为
(20)
式中：
为环境流体的密度，kg/m3；g为重力加速度，取g=10 m/s2。
由图3，根据力的平衡原理，可以得到
(21)
将式(18)~(20)代入式(21)，得到多孔介质的浮重表达式为
(22)
观察式(22)可以发现：多孔介质的浮重与固体骨架的质量、孔隙流体的密度和孔隙率等密切相关，而固体骨架的质量、孔隙流体的密度和孔隙率与多孔介质中的反应扩散和孔隙结构演化进程密切相关，也就是说多孔介质的浮重可以间接地反映发生在其中的反应扩散和孔隙结构演化进程。
2 数值方法
采用有限体积法[14]对数学模型式(6)~(13)进行数值求解，取控制容积的中心点作为控制节点，用控制容积积分法导出控制方程的离散形式，采用附加源项法处理边界条件。为了提高计算精度，扩散项采用具有二阶精度的中心差分格式进行离散，界面上的扩散系数采用两侧控制节点的调和平均表示，时间上采用全隐格式[15]。
图4所示为内部控制节点P离散示意图。对于任意内部控制节点P，可以得到离散后的控制方程为
(23)
式中：
；
；
；
；
。
图4 内部控制节点P离散示意图
Fig. 4 Discrete schematic diagram of inner control point P
图5所示为控制节点1离散示意图。对控制节点1，西界面便是x=0边界，根据附加源项法可以得到控制节点1的离散方程：
(24)
式中：
；
；
。
图5 控制节点1离散示意图
Fig. 5 Discrete schematic diagram of the first control point
图6所示为控制节点N离散示意图。对控制节点N，东界面便是x=0.5H边界，根据附加源项法可以得到控制节点N的离散方程：
(25)
式中：
；
；
。
图6 控制节点N离散示意图
Fig. 6 Discrete schematic diagram of the last control point N
综合以上离散结果，可以得到N个控制节点满足的矩阵方程：
(26)
由于数学模型的非线性以及数值离散过程中时间上采用全隐格式，导致式(26)的系数矩阵中含有待求参数，此处使用迭代方法进行求解，步骤如下：
1) 根据初始条件，
，
，假设第一时间水平多孔介质中固体颗粒的未反应核半径
，计算矩阵方程(26)，得到第一时间水平多孔介质中产物流体浓度的第1次计算值
；
2) 根据式(10)，计算第一时间水平固体反应物未反应核半径的第1次计算值
；
3) 比较
与
，C1与C0的差别，若同时满足精度要求，则进行第二时间水平的求解；若不能同时满足精度要求，则用代替，进行第一时间水平的第2次计算，根据矩阵方程(26)得到第一时间水平产物流体浓度的第二次计算值C2，然后根据式(10)，计算第一时间水平反应固体未反应核半径的第2次计算值
，比较与，C2与C1的差别，若同时满足精度要求，则进行第二时间水平的求解；若不能同时满足精度要求，则用代替，进行第一时间水平的第3次计算，…；
4) 在计算第n+1时间水平时，用第n时间水平得到的反应固体未反应核半径rc，作为初值进行计算，类比重复步骤3)，便可以得到所求时刻产物流体浓度和固体反应物未反应核半径的分布。
计算过程中对式(10)，(16)，(19)和(22)中的积分项均使用梯形法进行数值求解。
3 实例计算与分析
研究发生在一种多孔岩盐矿物中的反应扩散和孔隙结构演化过程，此岩盐矿物溶液浓度与密度之间满足
(27)
计算中所用的各类参数见表1。
在各种未开采的地下盐类矿床中，孔隙流体经过漫长地质年代的演化，常常是周围反应固体的饱和溶液，即流固反应停止。当矿床采用溶浸法开采时，随着反应流体进入孔隙，孔隙中的产物流体浓度降低(低于饱和浓度)，此时，才能发生进一步流固反应，固体反应物和流体反应物不断消耗，产物流体不断生成，如图7所示，孔隙中充满饱和产物流体的多孔圆片试件没入反应流体中，在0~60 min的时间段内，孔隙中的产物流体浓度、反应固体颗粒半径、孔隙率和固体转化率的分布和变化情况。
表1 计算参数表
Table 1 Parameters for computation
从图7可以看出：同一时刻，从多孔介质中心到表面，产物流体的浓度逐渐减小，反应固体颗粒半径逐渐减小，孔隙率逐渐增加，固体转化率逐渐增加。在t=0~60 min的过程中，随着时间的增加，多孔介质中产物流体的浓度逐渐降低，反应固体颗粒半径逐渐减小，孔隙率逐渐增加，固体转化率逐渐增加。容易想到，初始时刻多孔介质中充满产物流体的饱和溶液，而周围流体是纯反应流体，在浓度梯度的作用下，多孔介质中的产物流体通过孔扩散进入周围流体中，引起孔隙内产物流体浓度逐渐降低；与此同时，反应流体通过孔扩散进入多孔介质内部，与反应固体发生化学反应，导致反应固体颗粒半径减小，孔隙率增加，固体转化率增加。另外，由于流体反应物的扩散路径是从多孔介质外表面逐渐深入多孔介质内部的，这样靠近多孔介质表面的反应固体要先于多孔介质内部的反应固体发生化学反应，这使得靠近多孔介质表面的局部区域内反应固体颗粒半径最小，孔隙率最大，固体转化率最高。同理，由于产物流体的扩散路径是从多孔介质内部逐渐到达多孔介质的表面，产物流体的存在对流固化学反应起到了抑制作用，这也使得从多孔介质中心到表面，产物流体浓度逐渐减小，反应固体颗粒半径逐渐减小，孔隙率逐渐增加，固体转化率逐渐增加。
从图7还可以看出：与多孔介质的内部区域相比，多孔介质的表层区域内反应流体扩散迁移入多孔介质和产物流体扩散迁移出多孔介质的路径较短，因此扩散传质更容易，而产物流体扩散迁移出多孔介质和反应流体扩散迁移进多孔介质，都有利于流固化学反应的发生，即多孔介质表层区域内化学反应更容易发生，或者说化学反应的速率更快，因此，经历相同的时间，多孔介质表层附近区域内反应固体颗粒半径减少最大，孔隙率增加最大，固体转化率增加最大。
从图7还可以看出：由于流体的浓度扩散过程进行地非常缓慢，当时间t≤60 min时，多孔介质中的部分区域内，流体产物的浓度仍然为初始时的饱和浓度，导致流固化学反应在此区域内无法进行，因而这部分区域内反应固体颗粒半径、孔隙率、反应固体转化率没有发生改变，仍为初始值。然而，对比图7(b)，7(c)和7(d)可以发现：在0~60 min的时间段内，尽管多孔介质中反应固体仅有一小部分发生了转化，最高转化率低于9%(见图7(d))，反应固体颗粒半径变化小于4%(见图7(b))，但是，多孔介质的孔隙率最大却增加了50%(见图7(c))，这表明尽管多孔地质材料中的反应和扩散过程进行的非常缓慢，但是它们所引起的孔隙结构演化却是非常明显和不可忽视的。
图7 时间较短时，多孔介质中产物流体浓度、反应固体颗粒半径、孔隙率和固体转化率的分布和变化
Fig. 7 Liquid product concentration, solid reactant radius, porosity and conversion profiles in porous media for a short time
随着时间的增加，在浓度扩散的作用下，产物流体不断迁移出多孔介质，反应流体不断进入多孔介质，当时间足够长时，整个多孔介质中产物流体的浓度均低于饱和浓度，此时整个多孔介质中均有流固反应发生，如图8所示。从图8可以看出：从t=6 h开始，整个多孔介质中产物流体浓度均低于饱和浓度(见图8(a))，整个多孔介质中都有流固反应发生(见图8(d))。随着时间的进一步增加(t=12，18，24 h)，化学反应更明显地在多孔介质全长上发生，导致整个多孔介质中产物流体浓度进一步降低，反应固体颗粒半径进一步减小，孔隙率进一步增加，反应固体的转化率进一步增加。
从图8还可以看出：当t=24 h时，多孔介质表层的反应固体已经接近完全转化，多孔介质中心处的反应固体也有约50%的转化(见图8(d))，此时多孔介质表层的孔隙率接近初始孔隙率的7倍，多孔介质中心处的孔隙率约为初始孔隙率的4倍。对于多孔地质材料而言，如果在单纯力学作用下，孔隙率发生如此大地改变，几乎是不可能的。可见，与单纯的力学作用相比，化学反应所引起的多孔介质孔隙率改变更加明显。
综合图7和图8可以看出：多孔介质中产物流体浓度梯度最大的区域(即产物流体浓度扩散最快的区域)，与反应固体转化率最大的区域(即平均化学反应速率最大的区域)，以及孔隙率最大的区域，是重合的。
图8 时间较长时，多孔介质中产物流体浓度、反应固体颗粒半径、孔隙率和固体转化率的分布和变化
Fig. 8 Liquid product concentration, solid reactant radius, porosity and conversion profiles in porous media for a long time
也就是说扩散与反应和孔隙结构演化存在明显的对应关系，即浓度扩散引起产物流体浓度降低，加速化学反应；而化学反应导致反应固体颗粒半径减小，引起多孔介质孔隙结构的改变，使得多孔介质的孔隙增大，增大的孔隙进一步促进了产物流体的扩散传质。因此，多孔介质中的反应扩散与孔隙结构演化是相互影响，密不可分的。
通过对图7和图8的分析可以知道：浸没于反应流体中的非惰性多孔介质，反应流体通过孔扩散进入多孔介质内部与多孔介质骨架中的反应固体发生化学反应，产物流体通过孔扩散迁移出多孔介质，进入周围流体中，孔扩散和流固化学反应不断进行，反应固体逐渐变为产物流体，引起固体骨架空间结构的改变，直接导致多孔介质总孔隙率的变化，以及固体骨架质量和多孔介质浮重的改变，如图9所示。
从图9可以看到：随着时间的增加，反应固体总转化率逐渐增加，多孔介质总孔隙率逐渐增加，固体骨架的质量逐渐减小，浮重逐渐减小。随着反应流体扩散进入多孔介质，产物流体扩散出多孔介质，流固反应持续进行，反应区域不断地向多孔介质内部扩展，越来越多的反应固体被消耗，即反应固体总转化率增加；孔隙空间逐渐增加，即多孔介质总孔隙率逐渐增加；固体骨架逐渐消耗，即固体骨架质量逐渐减小。根据式(22)可以得到：由于反应流体密度小于产物流体密度，同时反应流体密度和产物流体的密度均小于反应固体的密度，因此消耗掉的反应固体的质量大于同体积的反应流体(或产物流体)的质量，而多孔介质所受到的浮力在整个过程中保持不变，这最终导致多孔介质的浮重逐渐减小。可以看到：当流固反应进行至24 h时，固体总转化率已经超过70%，总孔隙率增加至初始孔隙率的5倍多，固体骨架的质量减少至初始质量的50%以下，多孔介质浮重的相对变化也超过50%，这说明随着多孔介质中的反应扩散与孔隙结构演化过程的进行，多孔介质的物理性质发生了非常明显的改变。
从图9还可以看到：尽管多孔地质材料中的反应扩散与孔隙结构演化引起多孔介质总孔隙率、固体骨架质量和浮重的明显改变，但是多孔介质总孔隙率和固体骨架质量很难在实验室内进行实时测量。相反，浮重可以通过电子天平进行实时测量，并且目前高精度电子天平的精度已经可以到达0.01 mg，对于变化范围在克数量级上的浮重而言，具有足够的测量精度，并且电子天平与计算机联合使用，可以实现测量数据的实时采集和存储。因此，通过测量多孔介质的浮重，可以间接地获得多孔介质中反应扩散和孔隙结构演化的相关信息。
图9 多孔介质固体总转化率、总孔隙率、骨架质量和浮重随时间的变化
Fig. 9 Variation of whole conversion, whole porosity, solid skeleton mass and buoyed weight for porous media
由于多孔介质中的反应扩散与孔隙结构演化之间存在相互作用，这使得固体总转化速率、总孔隙率变化速率、固体骨架质量变化速率和浮重变化速率随时间呈现出明显的非线性特征，如图10所示。
从图10可以看出：随着时间的增加，固体总转化速率和总孔隙率的变化速率都先增大后减小，而固体骨架质量的变化速率和浮重的变化速率都先减小后增加。由微元体综合反应速率表达式(6)可以看出：流固反应的快慢程度与孔隙内产物流体的浓度成反比，只有当孔隙内产物流体的浓度低于饱和浓度时，流固反应才可能发生。由于初始时刻多孔介质中充满产物流体的饱和溶液，同时产物流体通过孔扩散进入周围流体中的速度是极其缓慢地，这使得一定时间内化学反应仅仅发生在多孔介质表层的局部区域内，反应固体仅有很少一部分参与化学反应，导致时间较短时反应固体的总转化速率较低，总孔隙率的变化速率较低，骨架质量变化率较低，浮重变化率较低。随着时间的增加，多孔介质中有更多的产物流体通过浓度扩散进入周围流体中，致使多孔介质中有更多区域内产物流体浓度低于饱和浓度，这样发生流固反应的区域不断扩大，直至在整个多孔介质内部全部发生。此时，反应固体的转化速率达到最大，固体质量减小速率也达到最大，总孔隙率增加速率最大，浮重减小速率最大。当时间进一步增加时，不再有新的未反应固体颗粒参与反应，相反，原来参与反应的固体颗粒半径却在逐渐减小，导致有效反应比表面积减小，化学反应的速率降低，反应固体转化速率减小，总孔隙率的增加速率减小，骨架质量的减小速率减小，浮重的减小速率减小。
由图10还可以看出：尽管固体总转化速率、总孔隙率变化速率、骨架质量变化速率和浮重变化速率均呈现出明显的非线性变化的特点，但是它们增加或减小的步伐却是一致的，非线性变化的时间段也是一一对应的。因此，尽管多孔介质中的反应和扩散过程发生在微观和细观水平上，很难进行准确测量，但它们与宏观水平上可以测量的浮重之间却存在一一对应关系，若能够通过宏观可测量的浮重，则获得多孔介质中微观和细观水平上的反应扩散信息，将是一个更加有意义的研究课题。
图10 多孔介质固体总转化速率、总孔隙率变化速率、骨架质量变化速率和浮重变化速率
Fig. 10 Variation rate of whole conversion, whole porosity, solid skeleton mass and buoyed weight for porous media
4 结论
1) 多孔地质材料中的反应扩散与孔隙结构演化是相互影响，密不可分的。
2) 尽管多孔地质材料中的反应和扩散过程进行的非常缓慢，但是它们所引起的孔隙结构演化却是非常明显和不可忽视的。
3) 与单纯的力学作用相比，化学反应所引起的多孔介质孔隙率改变更加明显。
4) 多孔地质材料中的反应扩散与孔隙结构演化过程与浮重之间存在一一对应关系。
参考文献：
[1] 高红波, 徐素国, 韩传磊, 等. 钙芒硝盐岩溶解机理细观研究[J]. 地下空间与工程学报, 2012, 8(4): 722-725.
GAO Hongbo, XU Suguo, HAN Chuanlei, et al. Study of meso-dissolution mechanism of glauberite[J]. Chinese Journal of Underground Space and Engineering, 2012, 8(4): 722-725.
[2] 刘中华, 胡耀青, 徐素国, 等. 钙芒硝溶解重结晶过程中孔隙演化规律试验研究[J]. 岩石力学与工程学报, 2011, 30(增1): 2743-2748.
LIU Zhonghua, HU Yaoqing, XU Suguo, et al. Experimental study of pore evolution law during dissolution-recrystallization for glauberite[J]. Chinese Journal of Rock Mechanics and Engineering, 2011, 30(Suppl.1): 2743-2748.
[3] 徐素国, 梁卫国, 赵阳升. 钙芒硝岩盐溶解的温度效应[J]. 矿业研究与开发, 2010, 30(3): 42-44.
XU Suguo, LIANG Weiguo, ZHAO Yangsheng. Temperature effect of glauberite salt rock dissolution[J]. Mining Research and Development, 2010, 30(3): 42-44.
[4] 刘中华, 徐素国, 胡耀青, 等. 钙芒硝盐岩溶蚀试验研究[J]. 岩石力学与工程学报, 2010, 29(增2): 3616-3621.
LIU Zhonghua, XU Suguo, HU Yaoqing, et al. Experimental study of corrosion of glauberite salt rock[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(Suppl.2): 3616-3621.
[5] 尹升华, 吴爱祥, 王少勇, 等. 溶液对斑岩铜矿废石的浸蚀机理及其力学特性[J]. 中南大学学报(自然科学版), 2013, 44(1): 309-314.
YIN Shenghua, WU Aixiang, WANG Shaoyong, et al. Erosion mechanism of waste rocks and its mechanical properties during dump leaching[J]. Journal of Central South University (Science and Technology), 2013, 44(1): 309-314.
[6] 杨保华, 吴爱祥, 刘金枝. 基于数字图像处理的堆浸散体孔隙连通性分析[J]. 矿冶工程, 2013, 33(1): 18-21.
YANG Baohua, WU Aixiang, LIU Jinzhi. Analysis of pore connectivity of granular ore media after heap leaching based on digital image processing[J]. Mining and Metallurgical Engineering, 2013, 33(1): 18-21.
[7] YANG Baohua, WU Aixiang, WANG Chunlai, et al. Three-dimensional simulation of pore scale fluid flow in granular ore media with realistic geometry[J]. Transactions of Nonferrous Metals Society of China, 2012, 22(12): 3081-3086.
[8] YANG Baohua, WU Aixiang, YIN Shenghua. Simulation of pore scale fluid flow of granular ore media in heap leaching based on realistic model[J]. Journal of Central South University of Technology, 2011, 18(3): 848-853.
[9] YIN Shenghua, WU Aixiang, LI Xiwen, et al. Mathematical model for coupled reactive flow and solute transport during heap bioleaching of copper sulfide[J]. Journal of Central South University of Technology, 2011, 18(5): 1434-1440.
[10] Dadvar M, Sahimi M. The effective diffusivities in porous media with and without nonlinear reactions[J]. Chemical Engineering Science, 2007, 62(5): 1466-1476.
[11] Valdes-Parada F J, Alvarez-Ramirez J. On the effective diffusivity under chemical reaction in porous media[J]. Chemical Engineering Science, 2010, 65(13): 4100-4104.
[12] Romano A, Horvath M I, Restani R. Evolution of porosity in the high-burnup fuel structure[J]. Journal of Nuclear Materials, 2007, 361(1): 62-68.
[13] Zheng L, Samper J, Montenegro L, et al. Multiphase flow and multicomponent reactive transport model of the ventilation experiment in opalinus clay[J]. Physics and Chemistry of the Earth, 2008, 33(1): 186-195.
[14] 陶文铨. 数值传热学[M]. 2版. 西安: 西安交通大学出版社, 2004: 28-43.
TAO Wenquan. Numerical heat transfer[M]. 2nd ed. Xi’an: Xi’an Jiaotong University Press, 2004: 28-43.
[15] 陶文铨. 计算传热学的近代进展[M]. 北京: 科学出版社, 2000: 233-248.
TAO Wenquan. Advances in computational heat transfer[M]. Beijing: Science Press, 2000: 233-248.
(编辑 杨幼平)
收稿日期：2013-10-13；修回日期：2013-12-08
基金项目(Foundation item)：国家自然科学基金资助项目(50974030，51304035)；教育部博士点基金资助项目(20130042120034)；博士后科学基金资助项目(20110491512)；中央高校基本科研业务费资助项目(N120401008) (Projects(50974030, 51304035) supported by National Natural Science Foundation of China; Project(20130042120034) supported by Specialized Research Fund for the Doctoral Program of Higher Education; Project(20110491512) supported by China Postdoctoral Science Foundation; Project(N120401008) supported by Fundamental Research Funds for Central Universities of China)
通信作者：姜元勇(1978-), 男，山东即墨人，博士，讲师，从事多孔相变材料中的流动与反应研究；电话：15002483816；E-mail: jiangyuanyong@mail.neu.edu.cn
