引言

随着中国天然气事业及储运配套设施建设的步伐推进，盐岩储气库逐渐引起业界关注，其建设进度不断加快^[1-3]。但由于盐岩本身所具有的层状分布、多泥质夹层、不溶杂质含量高等物性特征，腔体建成后在底部会堆积有一定厚度的不溶物^[4-7]。相较国内，国外盐层纯度高，产生堆积物数量小，其研究也多集中在水溶建腔方面，对于堆积物的提及较少^[8-9]。而国内盐层该问题较为突出，部分地区不溶物含量甚至可达50%，由此产生的影响也较为明显^[10-12]。堆积物的存在，不仅浪费了相当一部分库容体积，同时也为后续注气排卤及井下作业带来了安全隐患。因此有必要对堆积物分布特征及空隙体积进行研究，以为科学评价腔体库容体积、合理布置施工处理措施及提升储库整体建设质量提供理论基础与技术借鉴。

1 堆积物样品制备与测定 1.1 样品制备

由于腔底堆积物蓬松且浸满卤水，加之各方面条件限制，实现对其取样较为困难。根据相似性原理及现场工艺流程，自主设计模拟实验装置，利用腔体不同层位岩芯分别进行相同条件下的室内水溶试验以获取堆积物样品。

实验装置见图 1c，对岩芯进行中部钻孔(模拟井眼)，分别下入同心圆管(模拟溶腔内、外管)，在外管与孔壁间加注防护液(模拟溶腔油垫)，控制溶蚀高度。在实验过程中，通过控制流量、调整油垫与管柱相对位置、改变循环方式，控制溶腔形态发展。

结合现场盐层与不溶物含量数据以及溶腔阶段发展动态，选用典型层位岩芯进行试验。以金坛矿区某井为例，试验选用岩芯截面约10 cm，对其进行处理打磨，使重量相等，为1.17 kg(以最小者为基准)，层位及描述如表 1。

水溶实验后，对形成的不溶堆积物烘干、分析，并进行分布测定。

1.2 粒度分布测定

粒度分布测定方法较多，其中利用筛分法既可测定粒度，还可以绘制累积粒度分布曲线^[13]。对各组岩芯所获取堆积物岩样分别采用筛分法进行测定，进而累加，得到堆积物样品总体粒径数据及各分组数据(以3#岩芯样品为例)如图 2、图 3。

各分组所获取的堆积物样品粒度测定数据，其粒径范围及分级、累积质量百分数有所区别但差距有限，限于篇幅及数据重复性考虑，不将各组数据全部列出，仅以3#岩芯获取的不溶物岩样为例(下文处理同)，测定结果如图 4、图 5。由图 4、图 5可知，样品粒径多集中在30~40 mm及10~20 mm，而各分组粒径集中范围有所差异，如3#样品粒径相对较大，有30.8%部分大于40 mm。

2 分布特征研究 2.1 分布模型选定

为进行后续特征描述、评价，研究与之相符的分布函数。

由粒径分级质量图(图 4、图 5)可看出，样品总体符合双正态分布，但分组粒径却各有差别。鉴于分形理论近年来的发展创新及独特优势^[14-16]，采用分形函数对堆积物粒度特征进行描述分析。

根据Mandelbrot分形理论，多分散颗粒系统若符合分形分布，应满足^[17]

$ n_{(>x)}=C\cdot x^{-D} $

(1)

式中：

$n_{(>x)}$—系统中粒度大于$x$的颗粒数，无因次；

$x$—粒度，mm；

$C$—常数，无因次；

$D$—分形维数，用以表征颗粒组成集中性和均匀性的特征参数。

定义粒度$x$到$x+\textrm{d}x$间颗粒数$m(x)$

$ m(x)=\lim {\dfrac{{ - \Delta n_{(>x)}}}{{\Delta x}}} = - \dfrac{{\textrm{d}n_{(>x)}}}{{\textrm{d}x}} $

(2)

整理，得

$ m(x)=CD{x^{ - (D + 1)}} $

(3)

根据分布矩概念，将累积分布与分形维数联系起来，累积质量为

$ w_{(>x)}=\int_x^{{x_{\max }}} {\rho {C_{\textrm{V}}}{x^3}} m(x)\textrm{d}x \\~~~~~~~~ = \rho {C_{\textrm{V}}}C\dfrac{D}{{3 -D}}\left (x_{\max }^{3 -D} -{x^{3 -D}}\right) $

(4)

式中：$w_{(>x)}$—累积质量，g；

$C_{\textrm{V}}$—体积形状系数，无因次；

$\rho$—颗粒的密度，g/cm$^3$。

注意到总累积质量

$ {w_{\textrm{T}}}=\int_0^{{x_{\max }}} {\rho {C_{\textrm{V}}}{x^3}} m(x)\textrm{d}x = \rho {C_{\textrm{V}}}C\dfrac{D}{{3 - D}}x_{\max }^{3 - D} $

式中：$w_{\textrm{T}}$—总累积质量，g；

则大于特定粒径值的累积质量，即粒度上累积分布，可用函数表示为

$ N(x)=\dfrac{{w_{(>x)}}}{{{w_{\textrm{T}}}}} \times 100 = 100\left[{1-{{\left (\dfrac{x}{{{x_{\max }}}}\right )}^{3-D}}} \right] $

(6)

式中：$N(x)$—粒度上累积分布函数。

令$T=3-D$，转换，可得下累积分布函数$M(x)$

$ M(x){\textrm{ = 100}}{\left( {\dfrac{x}{{{x_{\max }}}}} \right)^T} $

(7)

式中：$M(x)$—粒度上累积分布函数。

2.2 分布拟合及特征参量计算

将粒度值和累积颗粒数分别取对数，进行分布拟合。3#岩芯不溶物岩样及样品总体拟合结果图如图 6、图 7。

对比实际数据线(实线)与拟合趋势线(虚线)可知，拟合偏差$R^2$在0.97~0.99，两者在双对数坐标系下呈线性相关，即其分布具有分形特点。特征参量即分形维数$D$计算结果如表 2。

根据计算结果，对比样品粒径数据可知，分形维数越大，颗粒集中性和均匀程度越差；相反，分形维数越小时，颗粒均匀性越好，颗粒尺寸也越集中。

2.3 分形维数与溶解度关系

由于腔底堆积物是各岩层段不溶物沉积累加的结果，为避免采用某单一特征参数表征其分布特性，并考虑颗粒级配产生的协同效应^[18]，分析分形维数与溶解度数据之间关系。岩芯各层分形维数和溶解度数据如表 3。

对分形维数与溶解度进行拟合，结果如图 8。由图 8可知，不溶物分形维数与溶解度呈线性相关。

为对两者关系的独立性做进一步分析研究，并验证与粒径组别范围的无关性，将5组样品不溶物中任意两组进行质量组合累加，得到10个新的组合，计算各组分形维数及溶解度并进行拟合(图 9)。

由图 9可看出，不溶物分形维数与溶解度成线性关系，且关系曲线可描述为

$ y=0.0091x+1.8186 $

式中：$x$—溶解度，无因次；

$y$—堆积物分形维数。

相关度在0.98以上。由此，堆积物分布分形维数可由其溶解度数据计算得到，且不溶物溶解度越大，形成的堆积物分形维数越大，亦即颗粒尺寸越不均匀，即堆积物分布特征可用分形函数进行表征。

3 空隙体积数值研究 3.1 空隙率数值模拟

受实验条件及因素限制，采用离散元仿真软件模拟不溶物堆积过程，并计算堆积空隙率。根据前述实验测定的堆积物分布结果，利用Origin拟合得出其分布特征参数，采用三维离散元方法模拟其堆积过程，并对模拟结果进行数值分析，从而得到各组不溶物岩样及总体堆积空隙率，结果如表 4。

3.2 堆积空隙率与分形维数关系

根据相关研究^[19]，随着颗粒粒度分形维数增加，其堆积空隙率线性减小，即其堆积结构越好，填充能力增强。据此对不溶物岩样数据进行空隙率与分形维数关系验证，验证结果如图 10。

由图 10可知，堆积物空隙率与其粒度分形维数成线性关系，且随分形维数增加，堆积空隙率下降。分形维数越大即颗粒越分散，粒级越多，这时颗粒间会互相填隙，空隙率从而变小。

为进一步探索堆积空隙率与分形维数间线性关系的独立性，研究相同分形维数不同粒度范围的颗粒堆积空隙率变化。所用各组样品颗粒分布如表 5。

采用同样处理方法，模拟得到各组堆积空隙率值，如表 6。

由表 6可知，分形维数相同粒度范围不同时，空隙率值也基本一致，即堆积空隙率与粒度范围无关，取决于其分形维数即颗粒分布特征。

4 空隙体积模型 4.1 可压缩堆积模型

腔底堆积物在形成过程中，颗粒间是相互影响的，小颗粒的存在会使大颗粒无法达到最紧密堆积，即松动效应；而大颗粒的存在会使小颗粒在靠近大颗粒壁面附近堆积，不能最密实，即附壁效应。同时，颗粒与颗粒间还存在互相压实的作用。

de Larrard等在考虑上述影响的基础上，为克服一般模型单一粒径假设的局限性，提出了可压缩堆积模型^[20-22]。根据该研究，在多粒级颗粒堆积系统中，假设第$i$级颗粒占主导地位，则系统各堆积参数为

$ \left\{ \begin{array}{l} {\gamma _i} = {\beta _i}/\left[ {1 - \sum\limits_{j = 1}^{i - 1} {\left[ {1 - {\beta _i} + {b_{ij}}{\beta _i}\left( {1 - 1/{\beta _j}} \right)} \right]{y_j} - \sum\limits_{j = i + 1}^n {\left( {1 - {a_{ij}}{\beta _i}/{\beta _j}} \right){y_j}} } } \right]\\ \gamma {\rm{ = }}\min \left( {{\gamma _i}} \right)\\ k = \frac{{{u_i}/{\beta _i}}}{{1/\mathit{\Phi - 1/}\gamma }} \end{array} \right. $

(8)

式中：$\gamma _i$—第$i$级颗粒为主时的虚拟堆积密实度，无因次；

$\beta _i$—第$i$级颗粒的剩余堆积密实度，无因次；

$y_j$—第$j$级颗粒所占体积分数，无因次；

$\gamma$—系统实际堆积密实度，无因次；

$K$—颗粒压实指数，无因次；

$\varPhi$—为系统实际堆积密实度，无因次；

$a_{ij}$，$b_{ij}$—第$j$级颗粒对第$i$级颗粒产生的松动效应和附壁效应系数，无因次，且$a_{ij} = \sqrt {1 - {\left( 1 - {d_j}\right)}^{1.02}}$，$j=i+1$，$\cdots$，$n$，$ b_{ij} = 1 - {\left( 1 - {d_i}\right)}^{1.50}$，$j=1$，$\cdots$，$i-1$；

$d_i$—第$i$级颗粒的粒径，mm。

4.2 模型计算分析

应用时，剩余堆积密实度可采用相关文献^[23-24]所得数据。根据Larrard的试验校准研究，压实指数在简单倾倒形式堆积时可取4.1，盐岩不溶物在堆积过程中主要受重力沉降作用并在此基础上压实形成^[25-27]，故按该情况处理。

由此，在确定堆积物各粒级尺寸及质量分数后，通过为各粒级剩余堆积密实度取值，求出系统虚拟堆积密实度，进而根据压实指数公式求得系统实际堆积密实度，最终确定堆积空隙率。

4.3 计算模型

上述堆积空隙率计算方法所针对的是盐岩溶腔各层位，要计算整个腔体的堆积物空隙率，需对各分层数据进行累加求解。假设共$n$个溶腔层位，第$i$层位盐岩含量即溶解度为$x_i$，层高$h_i$，溶腔半径$l_i$，则腔体不溶物总量$D_m$为

$ {D_m}{\textrm{ = }}\sum\nolimits_{i = 1}^n {(1 - {x_i}) \cdot {h_i} \cdot \pi } {l_i}^2 $

(9)

对于溶腔各层位，溶解度可通过岩性测井数据转换求解，进而求得各层位堆积物空隙率，再进行累加求解，即可得到腔体总体空隙率，表达式如下

$ \phi {\textrm{ = }}\dfrac{{\sum\nolimits_{i = 1}^n {{\phi _i} \cdot \lambda \cdot (1 - {x_i}){h_i} \cdot \pi } {l_i}^2}}{{{D_m} \cdot \lambda }}\\~~~~~= \dfrac{{\sum\nolimits_{i = 1}^n {{\phi _i} \cdot (1 - {x_i}){h_i} \cdot \pi } {l_i}^2}}{{{D_m}}} $

(10)

式中：$\lambda$—不溶物堆积膨胀系数，无因次。

由此，利用相关数据可对溶腔现场腔底堆积物空隙体积进行计算。

4.4 应用实例

实例井某溶腔段腔体测井解释数据如表 7。

采用可压缩堆积模型计算堆积空隙率，需首先确定堆积物各粒级组成，受条件限制，无法获取盐岩储库腔底不溶物具体颗粒组成。但根据前述分形研究结论，在分形维数相同时，空隙率值与粒度范围关系很小。故可根据上述基础数据，确定各层堆积后的分形维数，并选定50 mm(此时筛下累积质量为100%)作为最大粒径，以此确定各盐岩层位溶解后粒径分布，具体数据如表 8。

采用前述计算方法，对以上数据进行实例井计算，该溶腔阶段腔体堆积物空隙率计算结果如表 9。

根据计算结果，实例井某溶腔阶段堆积物总体空隙率为48.23%，现场通过后续声纳测腔分析得出的空隙率44.26%，本计算方法结果值稍大，但偏差不超过8.97%，故计算方法可行。

5 结论

(1) 通过室内溶腔实验获取了盐岩储气库腔底堆积物样品，采用筛分法对其粒度分布进行了实验测定，从而得到了样品粒径分级质量数据及总体累积质量数据。

(2) 腔底堆积物粒度分布可用分形函数进行表征，且分形维数作为特征参量能很好地描述颗粒分布的集中性与均匀程度。