本发明涉及油气井开采领域，尤其是一种根据目标井测井资料和目标井或其邻井的岩芯岩石力学测试资料预测横波时差剖面的方法。
背景技术：
油气井在进行钻井完井、分段分层、水力压裂等一系列施工设计中都需要利用到众多岩石力学参数作为设计依据，其中纵横声波时差对岩石力学参数的获取具有不可或缺的作用，其获取来源主要是油气井的测井资料和油气井所在区块的岩石力学测试资料。纵波时差剖面可由传统声波时差测井得到，成本较低，获取方便，但传统声波时差测井不能接收横波信号，因此横波时差剖面需要通过长源距声波、交叉偶极子声波等特殊的测井方式获取，测量周期长、成本高，油田上仅有数口重要的探井会进行横波测井，其余井的横波时差数据只能通过定深度取岩芯进行声波测试，且只能得到代表在实验围压条件下的岩芯对应油井深度处的横波数据，而不能得到整个井筒的横波时差剖面。
井场上，工程师往往将岩石力学测试的声波波速数据进行线性拟合的方法得到纵横波波速之间的一次线性关系，再将测井的纵波时差剖面代入该一次线性关系式中得到横波时差剖面。这种方法在一部分均质性较高的砂岩地层有可能会取得比较好的效果，而对于非均质性和各向异性较强的泥页岩储层或煤层，其纵横声波波速并不表现出很好的线性关系，因此通过此方法得到横波时差剖面缺乏较高的准确性。
对于纵横声波时差之间的关系，国内外学者做了大量研究。Xu-White考虑了岩石基质、孔隙大小及形状、孔隙内流体性质等因素的影响，提出一种计算横波时差的模型。Biot-Gassmann提出的横波时差计算模型组主要基于润湿颗粒估价的体积模量和剪切模量、孔隙空间模量、估价和孔隙流体的压缩模量等微观信息。张美玲提出一种以Gristensen经验方法为基础的横波时差计算模型，系数采用阻尼型高斯-牛顿法反演优化。然而诸多方法的实际应用效果由于其难以获得的参数或较为经验和局限适用对象等因素而受限，亟需一种适合普遍储层同时又满足精度要求的横波时差预测算法。
技术实现要素：
本发明的目的在于提供一种利用测井和岩芯资料获取横波时差剖面的方法，该算法原理可靠，操作简便，可克服现有技术的不足，能够解决求解过程中的病态问题，结果精度能够得到保证，为储层弹性模量、泊松比等岩石力学参数的获取提供了可靠的数据来源。
为达到以上技术目的，本发明提供以下技术方案：
(A)根据目标井岩性剖面对目标井按岩性划分层段，针对每个层段分别进行后续步骤的操作；
(B)收集目标井测井资料以及目标井或其邻井的对应目标层段的岩芯岩石力学测试资料，重点准备的数据包括纵波时差剖面、密度剖面，重点准备的岩芯岩石力学测试资料包括对应目标层段岩芯的纵波波速、密度、横波波速和实验围压条件；
(C)整理资料，保证在对应目标层段内取足三个以上岩芯的岩石力学测试数据，并按照下式分别构建观测输入矩阵、参数列向量和输出列向量：
Vs＝XA
Vs＝[vs(1) vs(2) … vs(m)]T
A＝[a1 a2 a3]T
x2(i)＝vp2(i)·Rc(i)
其中，Vs——m维输出列向量，T为转置符号；
m——对应目标层段岩石力学测试资料所涉及的岩芯数量，m>3；
vs(i)——编号为i的岩芯的横波波速，m/s，i取1,2,…,m；
X——m×3观测输入矩阵；
x1(i)、x2(i)、x3(i)——编号为i的岩芯的三种观测输入矩阵元素，i取1,2,…,m；
A——3维参数列向量；
a1、a2、a3——参数列向量的三个待定参数；
vp(i)——编号为i的岩芯的纵波波速，m/s，i取1,2,…,m；
Rc(i)——编号为i的岩芯的围压比，i取1,2,…,m；
ρ(i)——编号为i的岩芯的密度，g/cm3，i取1,2,…,m；
σcori(i)——编号为i的岩芯在原始地层条件下的围压，MPa，i取1,2,…,m；
σctes(i)——编号为i的岩芯在实验条件下的围压，MPa，i取1,2,…,m；
(D)构建三次Householder矩阵，并分别对观测输入矩阵和输出列向量做三次变换，该步骤的具体操作包括以下子步骤：
(D1)取根据m×3观测输入矩阵X第一列元素对m维列向量U1按照下式进行构建：
(D2)根据U1，按照下式构建m×m第一次Householder矩阵H1：
其中，I1——m×m单位矩阵；
||U1||——向量U1的2-范数；
(D3)对观测输入矩阵X和输出列向量Vs按下式做第一次Householder变换，得到第一次变换后的观测输入矩阵X1和第一次变换后的输出列向量Vs1：
Vs1＝H1Vs＝[vs'(1) vs'(2) … vs'(m)]T
其中，x1’(1)、x2’(i)、x3’(i)——第一次Householder变换后的观测输入矩阵X1中的元素，i取1,2,…,m；
vs’(i)——第一次Householder变换后的输出列向量Vs1中的元素，i取1,2,…,m；
(D4)取根据第一次变换后的m×3观测输入矩阵X1第二列后(m-1)个元素，对(m-1)维列向量U2按照下式进行构建：
(D5)根据U2，按照下式构建(m-1)×(m-1)第二次Householder矩阵H2：
其中，I2——(m-1)×(m-1)单位矩阵；
||U2||——向量U2的2-范数；
(D6)对第一次变换后的观测输入矩阵X1和输出列向量Vs1按下式做第二次Householder变换，得到第二次变换后的观测输入矩阵X2和第一次变换后的输出列向量Vs2：
其中，x1”(1)、x2”(1)、x2”(2)、x3”(i)——第二次Householder变换后的观测输入矩阵X2中的元素，i取1,2,…,m；
vs”(i)——第二次Householder变换后的输出列向量Vs2中的元素，i取1,2,…,m；
(D7)取根据第二次变换后的m×3观测输入矩阵X2第三列后(m-2)个元素，对(m-2)维列向量U3按照下式进行构建：
(D8)根据U3，按照下式构建(m-2)×(m-2)第三次Householder矩阵H3：
其中，I3——(m-2)×(m-2)单位矩阵；
||U3||——向量U3的2-范数；
(D9)对第二次变换后的观测输入矩阵X2和输出列向量Vs2按下式做第三次Householder变换，得到第三次变换后的观测输入矩阵X3和第一次变换后的输出列向量Vs3：
其中，x1”’(1)、x2”’(1)、x2”’(2)、x3”’(1)、x3”’(2)、x3”’(3)——第三次Householder变换后的观测输入矩阵X3中的元素；
vs”’(i)——第三次Householder变换后的输出列向量Vs3中的元素，i取1,2,…,m；
(E)分别取第三次变换后的观测输入矩阵X3和输出列向量Vs3的前三行元素构建线性方程组，并进行简单叠加运算求解参数列向量A：
(F)根据测井资料，按下式计算获得对应层段目标井横波时差剖面：
Δts＝106×vs-1
其中，Δts——横波时差，μs/m；
求得该对应目标层段的横波时差后，再按照以上流程求取下一个层段，直至所有要求层段求取完毕，最后拉通整个要求层段的横波时差剖面。
与现有技术比较，本发明紧密结合测井资料和岩芯资料，来源可靠，同时适用性广，对于均质的砂岩储层和具有非均质性和各向异性的泥页岩储层或煤层都能够计算得到较为准确的横波时差剖面，算法简便，易于实现程序化。
附图说明
图1是本发明一种利用测井和岩芯资料获取横波时差剖面的方法的技术流程图；
图2是某水平井储层段按本发明所述方法求取的横波时差剖面；
图3是某定向井储层段按本发明所述方法求取的横波时差剖面。
具体实施方式
下面结合附图对本发明做进一步说明，但本发明的保护范围不局限于以下所述。
位于鄂尔多斯盆地的某口水平井，其水平井段深度范围是1968m～3989m，水平井段所在的储层岩性为致密砂岩，现已获取该口井的测井资料以及其2口邻井同一目的层共计6块岩芯的岩石力学测试资料，但是测井资料缺乏所要求的水平井段的横波时差剖面，现欲对该口井按如图1所示的流程获取水平井段的横波时差剖面。具体操作如下：
(1)确定层段，要求层段为位于同一储层的整个水平井段，因此只需对整个水平井段进行横波时差剖面的预测计算；
(2)收集该水平井段的纵波时差、密度等测井资料，以及2口邻井的6块岩芯对应的深度、纵横波波速、密度、实验围压条件等岩石力学测试资料；
(3)整理岩芯资料，构建观测输入矩阵、参数列向量和输出列向量如下所示：X
(4)构建三次Householder矩阵，并对观测输入矩阵和输出列向量分别做三次变换，三次Householder矩阵分别为：
经过三次变换后的观测输入矩阵和输出列向量分别为：
(5)分别取三次Householder变换后的观测输入矩阵X3和输出列向量Vs3的前三行，构建线性方程组，解得参数列向量为：
(6)根据测井资料和参数向量A，计算横波时差剖面，结果如图2所示。
按照本发明所述方法计算得到的横波时差剖面曲线图(图2)，通过分析发现，在深度范围1968m～3150m内，由于密度曲线起伏不大，横波时差曲线与纵波时差曲线呈现比较良好的相关性；在深度范围3150m～3989m内，密度曲线起伏较大，横波曲线也受其影响偏离了与纵波时差的相关性。该口水平井目前已采用本发明的方法计算了横波时差剖面，并将所预测得到的横波时差剖面和其他测井资料用于预测弹性模量、泊松比等力学参数剖面，并投入到后续分段压裂施工应用中，从压裂效果上来看，多簇裂缝同时起裂延伸有效率≥80％，说明包括横波时差剖面在内的分段压裂设计依据准确性较高，进一步说明了本发明所述方法具有较高的可靠性。