盐岩因其低渗透性和良好的自愈合能力，被国际公认为石油、天然气等常规能源储存的理想介质。中国近年来实行能源储备发展战略，开始在盐岩地层中建设能源地下储库群，把盐岩储气库工程作为能源战略储备计划的重要部分。盐穴储库通常采用单井油垫法建造，但中国层状盐岩地层中单层厚度薄（60～100 m） 、软弱夹层多，单井油垫法水溶造腔遇到了很多工程技术难题，如造腔管柱破损或弯曲导致施工受阻、储气库运营压力失稳、腔体形状很难控制^[1]。溶腔长期处于水环境中，流体的流动对腔体形状扩展起关键作用，从流体角度出发揭示溶腔流场机理，对优化工艺参数、提高溶腔生产效率具有指导性作用。

小井间距双井水溶造腔技术^[2]是近20年来提出的新的造腔方法，具有造腔工艺简单、造腔速度快、注采效率高等优点。小井间距双井水溶造腔的基本流程为在相聚20 m左右的区域向目标盐层钻取2口井，以一口井为目标井，从另一口井钻取水平井和目标井连通，根据造腔需求布置管柱，如图1所示。与单井水溶造腔不同的是，小井间距双井水溶造腔只需布置注水管和排卤管，造腔时通过进水管向目标盐层注入淡水，使淡水溶解盐岩形成卤水，通过注水压力，卤水从排卤管排出。一段时间后，改变注水方向，这样交替进行注水和排卤。

国内外开展了一系列盐岩单井溶腔流场相关研究。Saberian^[3]对近几年在盐岩溶解速率方面的成果进行了总结。Durie和Jessen^[4]研究了自然对流作用下的盐溶现象。Saberian等^[5]指出，随着注水流速的加快，溶解速率加快。Jessen^[6]研究了在湍流条件下盐岩溶解特性。Husband和Chang等^[7–8]针对盐岩和卤水的边界自然传质对流问题展开了研究。Kazemi等^[9]通过玻璃容器和半径方向的模型分析溶腔中的流动形式和浓度分布情况。Jiang等^[10]对损伤盐岩的溶解特性进行研究。霍琰^[11]通过试验结合数值模拟研究了单井建腔期流场运移规律。李银平等^[12]针对建槽期流场的不溶物运动进行研究。陈锋等^[13]对储库注气排卤期可排卤水进行分析。Liu等^[14]利用室内试验研究了岩盐钻井开采中淹没射流场及循环场的特点。吴乘胜等^[15]提出单井对流法水溶开采的数学模型。赵志成等^[16–17]根据物质平衡原理及Fick扩散定律，建立盐岩溶蚀物理模型；根据流体力学基本原理、对流扩散理论及物质平衡原理，建立溶腔溶质传输–流体流动数学模型。姜德义等^[18]利用量纲分析法研究了含夹层腔体的浓度场分布规律。任松等^[19]利用选定物理法探究了单井造腔期流场特征。陈结等^[20]采用PIV技术研究了盐穴建造期夹层与卤水的相互作用。Nolen等^[21]对单井水溶过程进行了计算机数值模拟，对卤水浓度变化、腔体体积增长速率及腔体形状进行了分析。总体而言，现有研究很少有针对小井间距溶腔流场特征方面的研究，本文通过物理模型试验，结合数值分析，探索小井间距溶腔流场及浓度场特性，从溶腔形态扩展特性角度探究小井间距水溶造腔方法的合理性及优化性，以期为小井间距水溶造腔方法在工程实际中的应用提供理论基础。

1 试验平台搭建

考虑到国内深部盐岩难以获取、离散性较大、加工困难，国外高纯度盐岩均质性更好，而国内外高纯度盐岩溶解速率差异性较小^[22]，试验选取巴基斯坦喜马拉雅山区盐岩。盐样埋深较大，呈浅红色，盐岩NaCl纯度高，密度为2 358 kg/m³，含有少量K₂SO₄和其他泥质等杂物，组份比例见表1。原盐加工成尺寸长宽高为200 mm $ \times $ 100 mm $ \times $ 50 mm的盐砖试件，磨平处理。由于盐岩遇水溶解快，为保证试验的可靠性，对试件5个面涂抹防水胶密封，确保单个面的溶蚀。图2所示标记面为溶蚀面，其余5个为密封面。

小井间距水溶溶腔模型按照选定物理法则^[23]进行搭建（表2），通过可视化模型观察溶腔内流场运移。理想造腔腔体为近似竖直对称几何体，上顶板为油垫层，下壁为不溶物沉积层，两侧壁为褶皱状椎侧面，要实现全腔体的流场3维模型较难，腔体流场物理特征沿竖直中心平面具有高度对称性，因此试验模型采用500 mm $ \times $ 115 mm $ \times $ 650 mm的长方体，见图3。套管间距120 mm，模拟现实管距18 m。通过该近似腔体试验装置，以期获得实际小井间距水溶造腔中溶腔2维垂直断面流场特性。选用DK800–6玻璃转子流量计控制入口流量，流量为30 mL/min，套管为4 mm薄钢管。注水池高度高于装置顶面。正视面设置高分辨率CCD摄像机拍摄；背面设置取样小口，测取流场浓度。

2 试验结果及分析

根据拍摄的高清图像，结合已有流场研究和溶腔中流体各区域流动特征，把腔体内区域分为：1）边界溶蚀区：区域卤水既会发生横向扩散，也会由于重力作用下降。2）浮羽流区：注入淡水在浮力、重力等合力作用下加速度向上，初始速度向下，速度降为0后转而向上运动直至达到极限高度，然后向两侧运动；该区域是驱动流体运移的的关键区域，需要注意的是在不同流量下，尽管染色淡水上升高度有所不同，但运动到一定高度后，都不再向上运动，即都存在浮羽流极限高度。3）对流扩散区：根据是否受浮羽流影响分为强扩散区和弱扩散区。浮羽流效应形成的上升流与边界区的卤水下沉流对周围产生相反方向的卷吸作用，使得强扩散区低浓度与高浓度卤水相互扩散，卤水运移速度较快，以强迫对流为主；弱扩散区基本不受浮羽流影响，主要受浓度差驱动作用，卤水以自由扩散为主，区域卤水浓度在竖直方向上线性分布，横向浓度相差不大，盐岩溶蚀速度较慢。4）底部饱和区：此区域卤水近乎完全饱和，盐岩溶解极慢，可认为是静止区域，不溶物的主要沉积区，此区域盐岩几乎不溶解。几个区域相互联系，相互影响。如图4所示，通过染色图像处理，用箭头粗细代表流速大小，箭头越宽表示流速越大。

2.1 浮羽流区流场特征及微分模型

在浮羽流区不考虑黏性力时，染色剂流动受浮力羽流和对流扩散驱动。根据分层环境流体力学^[24]，淡水从注水管鞋流出，其密度远小于周围卤水密度，浮力远大于轻水自身重力，故淡水在管口以一定初速度流出，将形成小范围冲刷区域，如图5所示。浮力产生加速度逐渐减缓淡水速度，随后淡水开始上升，形成浮羽流。在初始动量和浮力的共同作用下，淡水流体微团弯曲向上，由于浓度差的存在，在上升过程中不断卷吸较重的周围流体，淡水微团本身逐渐变重，上升高度越高而浓度越低，密度越小，向上的浮力越来越小，乃至最后合力与速度方向相反，动量减小，直至达到浮羽流上升的极限高度。同时，流体微团在浮力羽流上升过程中受到卷吸作用，使得微团浓度升高，横向浓度差会产生横向推动力，使得分子向羽流中心线两边扩散。上升到一定高度，对流扩散高于浮力，浮羽流区域边界流体微团率先向两边扩散，浮羽流中轴线周边少部分流体继续上升。

图6为注水时不同流量与极限高度的关系曲线，根据流量与浮羽流极限高度的数学关系，可以划分出浮羽流区区域。

出水管冲刷区域远小于整个浮羽流区，可把羽流简化为有限空间点源浮羽流。忽略黏性切应力，采用柱坐标系，设羽流上升方向为x轴；径向方向为r轴；对应轴向速度为u；径向速度为v；半羽流特征厚度为b；取浮羽流区微团控制体，其密度为ρ；ρ_a表示周围环境密度（图5）。

根据质量守恒定理，控制体连续性方程为^[24]：

$\frac{{\partial u}}{{\partial {x}}} + \frac{1}{r}\frac{\partial }{{\partial r}}\left( {r\upsilon } \right) = 0$

(1)

引入卷吸假定，即认为羽流卷吸速度 ${\upsilon _{\rm e}}$ 与轴线流速 ${u_{\rm m}}$ 成比例，即 ${\upsilon _{\rm e}} \propto {u_{\rm m}}$ ；卷吸系数为 $\alpha $ ，则浮羽流连续性积分方程为：

$\frac{\rm d}{{{\rm d}x}}\int_0^{{b}} {\rho {{u}}{\rm d}r} = \alpha {\rho _{\rm{a}}}{{u}_{\rm{m}}}$

(2)

沿轴向羽流流量等于单位长度被卷吸流量，因此单位卷吸流量为：

${Q_{\rm e}} = \frac{\rm d}{{{\rm d}x}}\int_0^\infty u 2\text{π} r{\rm d}r = 2\text{π} b\alpha {u_{\rm m}}$

(3)

引入鲍辛奈斯克（Boussinesq）近似，密度差与浓度差联系起来：

${\rho _{\rm a}} - \rho = \rho \beta (c_{\rm m} - c)$

(4)

其中： $\;\beta $ 为水的体积膨胀系数， $\;\beta = 2.08 \times {10^{ - 4}}$ 。质量力只有重力，和射流一样，采用边界层方程，则羽流微元控制体运动方程为：

$\begin{aligned}[b]u\displaystyle\frac{{\partial u}}{{\partial x}} & + \upsilon \displaystyle\frac{{\partial u}}{{\partial r}} = \displaystyle\frac{{{\rho _{\rm a}} - \rho }}{\rho }g - \displaystyle\frac{1}{r}\displaystyle\frac{\partial }{{\partial r}}\left( {r\overline {u'\upsilon '} } \right) = \\ &\beta (c_{\rm m} - c)g - \displaystyle\frac{1}{r}\displaystyle\frac{\partial }{{\partial r}}\left( {r\overline {u'\upsilon '} } \right)\end{aligned}$

(5)

相似性假定羽流各断面流速分布、浓度分布存在自相似性，并服从高斯分布，即

$\frac{\mu }{{{\mu _{\rm m}}}} = \exp \left[ { - {{\left( {\frac{r}{b}} \right)}^2}} \right]$

(6)

$\frac{{\Delta \rho }}{{\Delta {\rho _{\rm m}}}} = \exp \left[ { - {{\left( {\frac{r}{{\lambda b}}} \right)}^2}} \right]$

(7)

其中， $\lambda $ 为浓度分布于流速分布的厚度比。对于式（5），当 $r = 0$ ， $r = \infty $ 时， $\overline {\mu '\upsilon '} = 0$ ， $\upsilon = 0$ 。将式（6）和（7）代入式（5），积分整理得运动微分方程为：

$\frac{\rm d}{{\rm d}x}\left( {\frac{\text{π} }{2}u_{\rm m}^2{b^2}} \right) = \text{π} \frac{{\Delta {\rho _{\rm m}}}}{\rho}g{\lambda ^2}{b^2}$

(8)

用C表示浮羽流中含有物的浓度，含有物质量守恒关系为圆柱坐标的定常紊流中的扩散关系式为：

$u\frac{{\partial C}}{{\partial x}} + \upsilon \frac{{\partial C}}{{\partial r}} = - \frac{1}{r}\frac{\partial }{{\partial r}}\left( {r\overline {u'C'} } \right)$

(9)

羽流中卤水浓度C与周围环境卤水浓度 ${C_{\rm a}}$ 之差为 $\Delta C$ ，则浓度扩散方程为：

$u\frac{{\partial \Delta C}}{{\partial x}} + \upsilon \frac{{\partial \Delta C}}{{\partial r}} = - \frac{1}{r}\frac{\partial }{{\partial r}}\left( {r\overline {u'\Delta C'} } \right)$

(10)

2.2 对流扩散区流场特征

从图4可以看出，根据流体特征，以浮羽流极限高度水平面为界限，对流扩散区又可分为强扩散区、弱扩散区。当注水后流体流出浮羽流区，强扩散区的强迫对流占主导，浮羽流中轴线周边流体会以较大速度上升，浮羽流中轴线靠外的流体微团由于上升过程中受卷吸作用影响较大，质量增加较快，流体垂直向上速度变小，因此示踪染色剂在强扩散区宏观表现为缓慢上升。弱扩散区流体微团在经过浮羽流影响区域后，浓度差使横向对流扩散成主力，流体微团沿浮羽流中轴线两边近似对称扩散，遇上岩层，表现出类似烟气流的顶棚射流效应，向端部迅速扩散。微团流向边界溶解交换岩层，边界层流体浓度高，加上浮羽流向上卷吸，使微团向下曲折环形流动，呈现出“回”型流动，直至到达出水管鞋水平液面。试验显示，在流量为30、50、70 mL/min时，强扩散区染色流体上升速度分别为1.1 $ \times $ 10^–4、1.9 $ \times $ 10^–4、4.8 $ \times $ 10^–4m/s，弱扩散区染色示踪粒子下降速度分别为9.5 $ \times $ 10^–4、13.9 $ \times $ 10^–4、27.8 $ \times $ 10^–4m/s，弱扩散区流体流动速度大于强扩散区速度。

对流扩散区在质量传输过程中速度场和浓度场相互作用、相互影响。注水管淡水流入腔体，强迫对流使溶质分子分散，对流扩散区呈非稳态浓度场。同时，重力作用下，溶液密度随着腔体深度增大，从上层至下层密度连续变化，因此溶液中溶质自然沉降，呈沉降扩散平衡分布，如图7（a）所示。若无注水、出水，溶腔内流体浓度平衡分布类似自然界中海水含盐量随深度的分布，具有垂直线性^[25]。淡水注入与排卤产生紊流，多尺度结构的扰动使得注水管鞋水平面以上浓度表现出宏观非线性分层。从试验分析可得，流体流动并不影响底部区溶液浓度的线性分层。

图7（b）为溶腔不同高度对应的浓度分布曲线。对流扩散区在溶腔平衡时，从顶层到管鞋平面浓度依次增大。对流扩散区浓度分为3个区域：浓度激增区、浓度平稳区、顶层浓度区。浓度激增区为管鞋平面以上小部分区域，约占腔体总体积9.6%，溶腔动态稳定时，此区内流体流动宏观运动现象明显，旋度大，流体的大尺度扰动使得区域内垂直方向上浓度变化率大，完成了对流扩散区低浓度到底部区高浓度的过渡。浓度平稳区是对流扩散区的主要区域，约占腔体总体积63.7%。溶腔溶解稳定时，对流扩散区流动结构尺度较小，只有注水管、出水管周边小部分区域结构尺度大；小尺度、低旋度的流体环境使得溶腔内物质交换速率低，强迫对流影响微小；自然对流扩散和沉降对流扩散相互作用，使得该区域浓度变化均匀，相对浓度激增区变化率小。顶层浓度区为溶腔顶部少部区域，约占腔体总体积的9.3%。结合图6，淡水流出注水管，羽流区中心区域（速度 ${u_{\rm m}}$ 区域）会冒出液面，冒出的液体浓度相对周边流体较低，在重力和液面表面张力的约束下，形成液体表面波，即涟波。涟波向四周散开，强制扰动该区域流体，同时冒出流体浓度在散出过程中对顶层区域流体稀释，加上本身的自由沉降扩散，使得顶层区域浓度降低。因此，图7（b）中顶层浓度最低，且变化率比浓度稳定区小。

2.3 边界层区流场特征

在边界层区域，盐层跟溶腔内流场流体发生质量交换，侧壁溶解使得边界层卤水浓度增大，密度升高，流体沿侧壁向下流动。在浮羽流极限高度平面附近，横向扩散速度较快，流体的扩散在小部分边界层形成紊动壁面浮射流区，若不考虑盐壁粗糙度，则边界层流体2种运移方式如图8所示。

普通边界层，即竖直盐壁面溶解过程区域，服从经典边界层理论。该边界溶蚀区内浓度和速度分布可用式（11）、（12）表示。

紊动壁面浮射流边界层溶解较为复杂，但根据格劳特壁面射流理论，近似服从经典边界层理论，同时，该区域具有壁面效应。

$C = {C_{\rm{d}}}{\left( {1 - \frac{y}{\delta }} \right)^2}$

(11)

$U = {U_1}\frac{y}{\delta }{\left( {1 - \frac{y}{\delta }} \right)^2}$

(12)

壁面浮射流区内外边界层存在不同的流动特征，对于无紊动射流影响的垂直壁面普通边界层，微观局部区域与单井溶解相同，区别主要在宏观运动路径上。可以引借吴乘胜等^[15]研究取得的单井水溶采矿中盐岩倾角侧壁盐岩溶解方程：

(13)

式中，w为溶解速度，g为重力加速度，D为扩散常数， $ \textit{z}$ 为盐壁高度， $\theta $ 为溶腔侧壁与竖直方向夹角， $\beta $ 为溶积膨胀系数， $\upsilon $ 为动力黏度系数， $C_{\rm d}$ 为淡水与盐表面直接接触的含盐量和淡水浓度之间的差值， $B = 1.44{w_0}$ 。

2.4 底部区流场特征

底部流体在浓度分层稳定后，浓度分布平稳，在轴向上近似呈线性分布。根据沉降扩散平衡理论，垂直方向上的浓度差导致垂直向上的扩散运动，重力作用下溶液垂直向下沉降运动。底部区盐溶液近似服从线性关系，即

$\frac{{{\rm d}{\rho _{\rm a}}}}{{{\rm d}y}} = cont$

(14)

底部区浓度高，区域内盐壁溶解速度慢，流速慢；造腔期完成后，最终接近饱和状态。在实际小井间距水平井中，不溶物沉淀在底部，溶腔内流场稳定后，底部区腔体形状固定，扩展量基本保持不变。值得一提的是，流体在溶腔内的湍流运动，随着时间的增加，其动能不断衰减，同时质量的传递会使整个溶腔内流体浓度升高，最终浓度场在深度方向出现非线性分层稳定状态。各分区之间没有特定分界面，主要取决于进入腔体内流体的初始特性。

2.5 腔体形态扩展分析

图9为小井间距水溶造腔方法下的腔体扩展实物。溶腔在腔底形成船体形的空间区域，该底部构造能够沉积大量的残渣而不占用腔体过多的高度。这种优势在后期腔体利用中能够发挥巨大的作用。由于卤水浓度在竖直方向上存在一定的浓度差，从而在腔底形成了55°～65°的溶蚀倾角。

图10为单井油垫法水溶造腔下的腔体扩展实物。可以看出该方法下形成的溶腔具有很好的对称性，这是由于造腔管柱位于溶腔的中心位置，使得同一水平面的卤水浓度具有很好的均一性，所以溶腔在同一水平面的盐壁会以相同的速度向外扩展。溶腔最后的形状为上大下小的近似梨形。

图11比较了2种造腔方法下的溶腔轮廓，图12比较了2种造腔方法下单位时间溶蚀的腔体体积。可以看出2条曲线的变化规律具有很好的相似性。在造腔前期，单位时间溶蚀的腔体体积随溶蚀时间的增加而增大；在造腔后期，单位时间溶蚀的腔体体积趋于平稳，基本上为一个稳定的数值。从2种方法比较来看，小井间距双井水溶造腔法在溶蚀效率上具有一定的优势，且优势在造腔前期更加明显。

图13比较了2种造腔方法下的溶腔立体形状。可以看出两者区别很大，单井油垫法水溶技术下的溶腔底部为陀螺形，在底部形成尖角，而小井间距双井水溶技术下的溶腔底部为平台形。相比较而言，平台形的底部能够存放更多的不溶物残渣。

结合流场分析结果可看出，在小井间距双井水溶造腔过程中，注入淡水形成的浮羽流对腔体竖直方向扩展起主导作用，紊动壁面浮射流区直接影响实际盐腔最大直径，约束腔体的形状。又因为卤水浓度在竖直方向上存在一定的浓度差，从而在腔底形成了55°～65°的溶蚀倾角。腔体形态扩展主要受流场运移和腔体卤水浓度影响。

与单井油垫法水溶造腔相比，小井间距双井水溶造腔技术的主要优势体现在造腔的前期和中期，在造腔后期，两者在这个参数上非常接近。从形状控制上来看，单井油垫法对溶腔形状的控制更加有利。小井间距双井水溶技术由于要来回倒井的原因，对油垫的控制更加困难，所以很容易形成腔体两端不对称的情况，在腔体形态上没有单井油垫法形成的规则。

3 结　论

1）根据流场区域现象，小井间距溶腔可分为4个主体区域：边界层、浮羽流区、对流扩散区、底部区，各个区域相互影响。浮羽流区是驱动整个腔体流场运移的关键区域。流体流入腔体，经过浮羽流区到达极限高度，部分向周边扩散，部分向顶部冒出。利用卷吸假定、Boussunesq近似及相似性假定，得出羽流微分方程为 $\displaystyle\frac{\rm d}{{{\rm d}x}}\left( {\displaystyle\frac{\text{π} }{2}u_{\rm m}^2{b^2}} \right) = \text{π} \displaystyle\frac{{\Delta {\rho _{{\rm m}}}}}{\rho }g{\lambda ^2}{b^2}$ 。

2）边界层区域可分为普通边界区和紊动壁面浮射流区。普通边界层符合基本边界层浓度分布和速度分布，紊动壁面浮射流区影响腔体最大直径，溶解较为复杂。底部区示踪特征不明显，浓度呈线性分布。

3）单井造腔中同一水平面的卤水浓度具有很好的均一性，溶腔在同一水平面的盐壁会以相同的速度向外扩展，所以形状比较规则。小井间距溶腔底部可以形成船体形的空间区域，这种底部构造能够沉积大量的残渣而不占用腔体过多的高度，可以更高效地储存石油和天然气。

4）在造腔前期和中期，小井间距双井水溶法的造腔效率比单井油垫法高，在造腔后期二者基本接近。需要注意的是，由于来回倒井和油垫难以控制的问题，小井间距双井水溶造腔容易形成两端不对称的腔体，这也是下一步准备解决的问题。