一种二维水平定向钻进轨迹计算方法、装置及介质

1.本发明涉及岩土工程领域，具体来说，涉及一种二维水平定向钻进轨迹计算方法、装置及介质。


背景技术：

2.水平定向钻进指的是施工过程中在较浅深度范围的地表土层，按照预先设计的钻孔轨迹先行钻进一个较小的导向孔，然后卸下导向钻头换上扩孔钻头反向对已经形成的导向孔进行扩孔，与此同时扩孔钻头后牵扯待铺管线进入孔内完成管线铺设。作为一种近些年发展迅速非开挖施工技术方法，水平定向钻进施工灵活，成本低，对周边环境扰动小，正逐渐取代传统的开挖施工方法成为城市内地下管线安装、更换的主流施工手段。水平定向钻进穿越距离的不断提高，意味着穿越区段内地质条件更加多变，面临的工程问题也随之更加复杂，重复无效的钻进对于工程成本的提高难以估计，因对长距离穿越的水平定向钻进轨迹的设计优化问题进行研究，找到兼具安全稳定和经济效益的最优钻孔轨迹具有重要的工程意义。
3.在实际工程中，绝大多数仍是在钻进施工时依靠随钻测量技术根据现场情况实时更新钻进轨迹。考虑具体的轨迹模型和参数对钻孔轨迹进行设计优化，主流方法仍是以曲线-直线组合的空间曲线模型为主，国内众多学者针对不同的优化目标提出了相应的优化模型和方法展开了深入研究(具体参考文献1、2)；国外亦有学者提出了以悬链线轨迹形式为主的钻孔轨迹设计优化方法，获得了回拖力较小，管道应力分布合理，管线寿命较长的轨迹形式。
4.参考文献1：杨先亢等，悬链法在水平定向钻轨迹设计中的应用，非开挖技术，2010年4月第1期。
5.参考文献2：陈军红，非开挖导向钻进钻轨迹优化设计，山西建筑，2007年1月第33卷第2期。
6.本实施例提供的背景描述用于总体上呈现本公开的上下文的目的。除非本实施例另外指示，在该章节中描述的资料不是该申请的权利要求的现有技术并且不要通过包括在该章节内来承认其成为现有技术。


技术实现要素：

7.针对相关技术中的上述技术问题，本发明提出一种二维水平定向钻进轨迹计算方法，其包括如下步骤：
8.s1，根据获取的水平定向钻进轨迹起点与终点的水平穿越长度a，纵向高差δz，利用悬链线钻孔轨迹模型建立所述水平定向钻进轨迹的第一轨迹方程，利用“五段法”钻孔轨迹模型建立所述水平定向钻进轨迹的第二轨迹方
9.^程并建立所述第二轨迹方程的约束条件；根据所述第一轨迹方程的纵坐标zi和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数；根据最大泥浆压力map和最小泥浆
压力mdp建立所述第二轨迹方程的泥浆压力约束条件；
10.s2，将所述第二轨迹方程的轨迹参数转化为改进的径向移动算法irmo
11.^的基础变量矩阵，根据所述第一轨迹方程的纵坐标zi和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数并结合irmo算法实现“五段法”钻孔轨迹与悬链线轨迹的最佳拟合，确定两者的最佳拟合区间；
12.s3，将最佳拟合区间确定为轨迹长度优化的最佳优化区间，结合所述最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件并采用irmo算法在最佳优化区间内搜索求解轨迹长度多维非线性目标函数的最优值，得到最佳优化区间内的钻进轨迹长度的最优钻进轨迹长度。
13.具体的，所述“五段法”钻孔轨迹模型的六个轨迹参数如下：入口段直线段长度l1；入口段直线段倾角θ1；入口段过渡圆弧半径r2；中心轨迹直线段长度l3；出口段过渡圆弧半径r4；出口段直线段倾角θ2。
14.具体的，所述步骤s2还包括：
15.s2.1根据所确定的六个轨迹参数，确定每个轨迹参数的初始取值范围x
jmax
和x
jmin，
其中，j＝nod＝6；
16.s2.2在步骤s2.1的取值范围内随机生成nop个初始位置点，由nop个初始位置点建立初始种群，所述初始位置点的数值信息通过计算式得出；
17.x(i，j)＝x
jmin
+rand(0，1)(x
jmax-x
jmin
)
18.其中，x(i，j)为生成的第i个初始位置点的第j个轨迹参数；rand(0，1)为在0到1之间的随机数；
19.s2.3确定生成新的预位置点的更新条件，具体的如下，
[0020][0021][0022]
其中，y
i，jk
是指第k代新生成的第i个位置点的第j个参数变量；centre
jk
是指第k代中心位置的第j个参数变量；wk为惯性权值，随迭代次数递减，用于决定算法的收敛速度；k是指当前迭代次数；g是指最大迭代次数；
[0023]
中心位置随着当代最优位置的rbestxk和全局最优位置gbestxk的移动而移动，具体的，
[0024]
centre
k+1
＝centrek+0.4(gbestx
k-centrek)+0.5(rbestx
k-centrek),
[0025]
其中，centrek为第k代的中心位置；centre
k+1
为第k-1代的中心位置；
[0026]
s2.4采用更新条件生成新的预位置点y
ik
，并通过拟合程度评估函数f1(x)计算出各预位置点y
ik
对应的适应度值；
[0027]
s2.5更新位置信息，具体的，将第k代预位置点的适应度值fitness(y
ik
)与第k-1代位置点的适应度值fitness(x
ik-1
)进行比较，若fitness(y
ik
)《fitness(x
ik-1
)，则需要更新位置点信息，令fitness(x
ik
)＝fitness(y
ik
)，x
ik
＝y
ik
；否则，令fitness(x
ik
)＝fitness(x
ik-1
)，x
ik
＝x
ik-1
；其中，x
ik
为第k代第i个位置点；x
ik-1
为第k-1代第i个位置点；y
ik
为第k代
生成的第i个预位置点；
[0028]
s2.6确定当代最优位置rbestxk和全局最优位置gbestxk，具体的，将更新后的第k代中的fitness(x
ik
)中的最小值的位置点作为当代最优位置rbestxk；将更新后的前k代中的fitness(x
ik
)中的最小值的位置点作为全局最优位置gbestxk；
[0029]
s2.7更新全局最优位置，具体的，比较当代最优位置的rbestxk的适应度值fitness(rbestxk)与全局最优位置gbestxk的适应度值fitness(gbestxk)；若fitness(rbestxk)《fitness(gbestxk)，则需要更新全局最优位置，令fitness(gbestxk)＝fitness(rbestxk)，gbestxk＝rbestxk；否则，令fitness(gbestxk)＝fitness(gbestx
k-1
)，gbestxk＝gbestx
k-1
。
[0030]
具体的，钻孔轨迹长度适应度函数f2(x)＝l1+r2·
θ1+l3+r4·
θ2+l5；
[0031][0032]
其中，l1为入口段直线段长度、l3为中心轨迹直线段长度、l5为出口段直线段长度，r2为入土圆弧段的半径、r4为出口段过渡圆弧半径，θ1为入口段直线段倾角、θ2为出口段直线段倾角。
[0033]
具体的，所述步骤s3具体包括：
[0034]
获得六个轨迹参数的最佳优化区间x
bjmax
和x
bjmin
，选择钻孔轨迹长度适应度函数f2(x)代替步骤s2.4中的拟合程度评估函数f1(x)，利用irmo算法对最佳优化区间内轨迹的最短长度进行优化搜索，其算法搜索步骤同s2.1-s2.7。
[0035]
第二方面，本发明的另一实施例公开了一种二维水平定向钻进轨迹计算装置，其包括如下单元：
[0036]
模型建立单元，用于根据获取的水平定向钻进轨迹起点与终点的水平穿越长度a，纵向高差δz，利用悬链线钻孔轨迹模型建立所述水平定向钻进轨迹的第一轨迹方程，利用“五段法”钻孔轨迹模型建立所述水平定向钻进轨迹的第二轨迹方程并建立所述第二轨迹方程的约束条件；根据所述第一轨迹
[0037]
^
[0038]
方程的纵坐标zi和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数；根据最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件；
[0039]
最佳拟合区间计算单元，用于将所述第二轨迹方程的轨迹参数转化为改进的径向移动算法irmo的基础变量矩阵，根据所述第一轨迹方程的纵坐标和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数并结合irmo算法实现“五段法”钻孔轨迹与悬链线轨迹的最佳拟合，确定两者的最佳拟合区间；
[0040]
轨迹长度计算单元，用于将最佳拟合区间确定为轨迹长度优化的最佳优化区间，结合所述最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件并采用irmo算法在最佳优化区间内搜索求解轨迹长度多维非线性目标函数的最优值，得到最佳优化区间内的钻进轨迹长度的最优钻进轨迹长度。
[0041]
具体的，所述“五段法”钻孔轨迹模型的六个轨迹参数如下：入口段直线段长度l1；入口段直线段倾角θ1；入口段过渡圆弧半径r2；中心轨迹直线段长度l3；出口段过渡圆弧半
径r4；出口段直线段倾角θ2。
[0042]
具体的，所述最佳拟合区间计算单元还包括：
[0043]
第一计算单元，用于根据所确定的六个轨迹参数，确定每个轨迹参数的初始取值范围x
jmax
和x
jmin
，其中，j＝nod＝6；
[0044]
第二计算单元，用于在第一计算单元获取的取值范围内随机生成nop个初始位置点，由nop个初始位置点建立初始种群，所述初始位置点的数值信息通过计算式得出；
[0045]
x(i，j)＝x
jmin
+rand(0，1)(x
jmax-x
jmin
)
[0046]
其中，x(i，j)为生成的第i个初始位置点的第j个轨迹参数；rand(0，1)为在0到1之间的随机数；
[0047]
第三计算单元，用于确定生成新的预位置点的更新条件，具体的如下，
[0048][0049][0050]
其中，y
i，jk
是指第k代新生成的第i个位置点的第j个参数变量；centre
jk
是指第k代中心位置的第j个参数变量；wk为惯性权值，随迭代次数递减，用于决定算法的收敛速度；k是指当前迭代次数；g是指最大迭代次数；
[0051]
中心位置随着当代最优位置的rbestxk和全局最优位置gbestxk的移动而移动，具体的，
[0052]
centre
k+1
＝centrek+0.4(gbestx
k-centrek)+0.5(rbestx
k-centrek),
[0053]
其中，centrek为第k代的中心位置；centre
k+1
为第k-1代的中心位置；
[0054]
第四计算单元，用于采用更新条件生成新的预位置点y
ik
，并通过拟合程度评估函数f1(x)计算出各预位置点y
ik
对应的适应度值；
[0055]
第五计算单元，用于更新位置信息，具体的，将第k代预位置点的适应度值fitness(y
ik
)与第k-1代位置点的适应度值fitness(x
ik-1
)进行比较，若fitness(y
ik
)《fitness(x
ik-1
)，则需要更新位置点信息，令fitness(x
ik
)＝fitness(y
ik
)，x
ik
＝y
ik
；否则，令fitness(x
ik
)＝fitness(x
ik-1
)，x
ik
＝x
ik-1
；其中，x
ik
为第k代第i个位置点；x
ik-1
为第k-1代第i个位置点；y
ik
为第k代生成的第i个预位置点；
[0056]
第六计算单元，用于确定当代最优位置rbestxk和全局最优位置gbestxk，具体的，将更新后的第k代中的fitness(x
ik
)中的最小值的位置点作为当代最优位置rbestxk；将更新后的前k代中的fitness(x
ik
)中的最小值的位置点作为全局最优位置gbestxk；
[0057]
第七计算单元，用于更新全局最优位置，具体的，比较当代最优位置的rbestxk的适应度值fitness(rbestxk)与全局最优位置gbestxk的适应度值fitness(gbestxk)；若fitness(rbestxk)《fitness(gbestxk)，则需要更新全局最优位置，令fitness(gbestxk)＝fitness(rbestxk)，gbestxk＝rbestxk；否则，令fitness(gbestxk)＝fitness(gbestx
k-1
)，gbestxk＝gbestx
k-1
。
[0058]
具体的，钻孔轨迹长度适应度函数f2(x)＝l1+r2·
θ1+l3+r4·
θ2+l5；
[0059][0060]
其中，l1为入口段直线段长度、l3为中心轨迹直线段长度、l5为出口段直线段长度，r2为入土圆弧段的半径、r4为出口段过渡圆弧半径，θ1为入口段直线段倾角、θ2为出口段直线段倾角。
[0061]
具体的，所述轨迹长度计算单元具体包括：
[0062]
获得六个轨迹参数的最佳优化区间x
bjmax
和x
bjmin
，选择钻孔轨迹长度适应度函数f2(x)代替第四计算单元中的拟合程度评估函数f1(x)，利用irmo算法对最佳优化区间内轨迹的最短长度进行优化搜索，其算法搜索步骤同第一计算单元至第七计算单元。
[0063]
第三方面，本发明的另一个实施例公开了一种非易失性存储器，所述存储上存储有指令，所述指令被处理器执行时，用于实现上述任一项所述二维水平定向钻进轨迹计算方法。
[0064]
本发明根据水平定向钻进轨迹的悬链线轨迹模型和“直线-圆弧-直线-圆弧-直线”五段式轨迹模型，在掌握工程要求和区域地质条件的情况下，分别建立拟合度和轨迹长度的求解方程，然后利用irmo算法将拟合度求解方程作为算法的适应度函数获得钻孔轨迹的最佳优化区间，进一步将钻孔轨迹长度作为算法的适应度函数在最佳优化区间内搜索求解轨迹长度最优值，同时加入了以泥浆压力理论为基础的轨迹可行性检验条件保证轨迹优化结果的合理性和可行性。该方法能够有效优化水平定向钻进轨迹的长度，精确计算出最优结果下的轨迹参数。
附图说明
[0065]
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
[0066]
图1是本发明实施例提供的一种二维水平定向钻进轨迹计算方法流程图；
[0067]
图2是本发明实施例提供的悬链线模型示意图；
[0068]
图3是本发明实施例提供的“五段法”钻孔轨迹参数及其几何约束示意图；
[0069]
图4是本发明实施例提供的轨迹可行性检验流程示意图；
[0070]
图5是本发明实施例提供的irmo算法实现流程图；
[0071]
图6是本发明实施例提供的series4条件下连续搜索20次的拟合结果sos波动示意图；
[0072]
图7是本发明实施例提供的series4条件下连续搜索20次拟合轨迹图示意图；
[0073]
图8是本发明实施例提供的irmo算法连续搜索20次的钻进距离最优轨迹示意图；
[0074]
图9是本发明实施例提供的连续搜索20次轨迹长度优化前后对比示意图；
[0075]
图10是本发明实施例提供的一种二维水平定向钻进轨迹计算装置示意图；
[0076]
图11是本发明实施例提供一种二维水平定向钻进轨迹计算设备的示意图。
具体实施方式
[0077]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员所获得的所有其他实施例，都属于本发明保护的范围。
[0078]
实施例一
[0079]
参考图1，本实施例公开了一种二维水平定向钻进轨迹计算方法，其包括如下步骤：
[0080]
s1，根据获取的水平定向钻进轨迹起点与终点的水平穿越长度a，纵向高差δz，利用悬链线钻孔轨迹模型建立所述水平定向钻进轨迹的第一轨迹方程，利用“五段法”钻孔轨迹模型建立所述水平定向钻进轨迹的第二轨迹方程并建立所述第二轨迹方程的约束条件；根据所述第一轨迹方程的纵坐标和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数；根据最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件；
[0081]
参考图2，水平钻进时，钻杆钻进水平距离很长，其直径相对来说很小，可以近似将水平定向钻进钻杆看成一条悬链线，基于悬链线的标准形式方程，获得水平定向钻进轨迹设计中的悬链线轨迹设计模型。本实施例以悬链线左端点为坐标原点，左右端点的横纵坐标差分别为a(a为地表起终点的水平距离)和δz(地表起终点的竖向距离)，轨迹顶点坐标为(x0,z0)，得到的悬链线轨迹方程和图像如下：
[0082][0083]
其中，q为管线重量，单位为n/m，np为回拖力，单位为kfg。
[0084]
基于“直线-曲线”组合的基本钻孔轨迹形式，依据其二维平面内的几何关系，确定水平定向钻进轨迹的优化参数，从而获得“五段法”钻孔轨迹的基本模型，表达出钻孔轨迹上任意一点的坐标，在此基础上与悬链线轨迹进行拟合并计算出钻孔轨迹的长度。
[0085]
参考图3，在同一竖直平面内，已知圆弧的张角γ和圆弧半径r即可确定圆弧内任意一点的坐标参数，公式(2)。以入土段设计为例(图3以入土段为例进行说明，则此时入土段的直线段长度为l1)，若已知任一相邻直线段的长度l、水平倾角θ(水平倾角是指钻进方向与水平面之间的夹角；图3以入土段为例进行说明，则此时入土段的水平倾角为θ1)及某一端点的坐标，便可根据坐标关系获得另一直线上任意一点的坐标参数，公式(3)、(4)。
[0086][0087]
α＝θ+γ
ꢀꢀ
(3)
[0088]
[0089]
其中，i为是标号，为第i个水平微分段轨迹，xi为微分段的端点水平坐标，α为另一端直线的水平夹角，在以入土段为例中可指中心水平段的水平倾角，d
x
为钻孔轨迹沿水平方向x微分每一小段长度(水平微分长度)，l1为入口段直线段的长度；γ为圆弧张角，δγ为圆弧张角变化量。
[0090]
基于公式(4)-(6)，即可通过以下四个轨迹参数：直线段倾角θ(顺时针为正，逆时针为负)、圆弧的张角γ、圆弧半径r以及直线段长度l确定由任意“曲线-直线”组成的钻孔轨迹。
[0091]
本实施例仅考虑“直线-圆弧-水平直线-圆弧-直线”的五段法轨迹形式，因此中心轨迹的设计通常保持水平(倾角α＝0
°
)，直线段倾角θ与圆弧张角γ可以在几何条件上得到统一：θ1＝γ1、θ2＝γ2，由此得到本实施例所需的“五段法”轨迹的所有优化参数如下：
[0092]
·
入口段直线段长度l1；
[0093]
·
入口段直线段倾角θ1；
[0094]
·
入口段过渡圆弧半径r2；
[0095]
·
中心轨迹直线段长度l3；
[0096]
·
出口段过渡圆弧半径r4；
[0097]
·
出口段直线段倾角θ2；
[0098]
本实施例统一悬链线轨迹和“五段法”轨迹的平面坐标系，已知轨迹起点与终点的水平穿越长度为a，纵向高差为δz，在满足水平定向钻孔轨迹的基本要求上(式7-8)，分别得到了拟合的计算公式和轨迹长度的计算公式。
[0099]
中心轨迹上任一点的埋深可由式6确定，需满足：h＝zi≥h；
[0100]
其中h为规范要求的最小埋深。
[0101]
轨迹终点(出土点)约束范围：
[0102]
其中，b为水平长度误差控制范围，c为纵向高差误差控制范围。
[0103]
为保证“五段法”轨迹与悬链线轨迹的最大相似，采用数理统计中的最小二乘法，将悬链线轨迹纵坐标(式1)与“五段法”轨迹纵坐标zi(式4)的差值平方和sos作为评价拟合程度的计算公式，sos的计算结果越小，表明两者的拟合程度越高，式5：
[0104][0105]
依据“五段法”六个轨迹参数，可直接获得轨迹长度length的计算公式，式6-7：
[0106]
length＝l1+r2·
θ1+l3+r4·
θ2+l5ꢀꢀ
(6)
[0107][0108]
其中，l是五段法中每一个直线段的长度，即本实施例中的直线-圆弧-水平直线-圆弧-直线，其中l1是指五段法中的第一段直线即入口段直线段长度、l3指五段法中的中心轨迹直线段长度、l5指五段法中出口段直线段长度，r是五段法中的圆弧的半径、r2是第二段(入土圆弧段)的半径；r4是出口段过渡圆弧半径。
[0109]
本实施例在轨迹优化的过程中，考虑了工程实际中井喷或塌孔等不利情况，对优化后的轨迹是否合理、可行进行判定。基于“五段法”钻孔轨迹的坐标公式，结合钻孔的泥浆压力理论，逐段对钻孔轨迹上各点的可行性进行检验，保证所设计优化出的钻孔轨迹符合工程实际。
[0110]
基于公式6，取dx＝1.0m，计算并比较轨迹上各点的最大泥浆压力map和最小泥浆压力mdp，式8：
[0111][0112]
最小泥浆压力值(mdp)通常由静态钻孔泥浆压力与泥浆压力损失组成，可由式9计算：
[0113]
mdp＝ps+p
l
ꢀꢀꢀ
(9)
[0114]
式中：ps为静态钻孔泥浆压力(pa)；p
l
为泥浆压力损失(pa)；
[0115]
以土体剪切破坏模式为基础的最大泥浆压力(map)计算方法主要基于经典孔扩张理论模型，即delft方程，式：
[0116]
map＝p
lim
ꢀꢀꢀ
(10)
[0117]
式中：p
lim
为delft方程计算的最大泥浆压力(pa)；
[0118]
一旦轨迹中任意一点的g1(zi)≤g2(zi),则视为该轨迹存在井喷风险，将其视为无效轨迹并舍弃。将此过程作为irmo算法优化过程中的判定条件之一，仅保留满足钻孔稳定性检验的有效优化轨迹进入下一次优化，最终获得最优的有效轨迹。
[0119]
s2，将所述第二轨迹方程的轨迹参数转化为改进的径向移动算法irmo的基础变量矩阵，根据所述第一轨迹方程的纵坐标和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数并结合irmo算法实现“五段法”钻孔轨迹与悬链线轨迹的最佳拟合，确定两者的最佳拟合区间；
[0120]
通过l1、θ1、r2、l3、r4、θ2六个轨迹参数即可确定一个完整的“五段法”钻孔轨迹。因此，将任意一条钻孔轨迹中的六个轨迹参数整合为irmo算法中的一个解向量(m＝6)，n条钻孔轨迹共同组成irmo算法中的一个基础变量矩阵，即可实现irmo算法与水平定向钻进轨迹的衔接。由六个轨迹参数组成的基础变量矩阵如式11所示。
[0121][0122]
依据以上基础变量矩阵，对各个位置的变量施加约束条件，在matlab平台上即可实现基于irmo算法的水平定向钻进轨迹设计优化。在实例计算过程中，拟合和钻孔轨迹长度的优化均使用了irmo算法进行，因此需要依据使用场景的不同更换算法的适应度函数。
[0123]
·
拟合适应度函数：f1(x)＝sos
[0124]
·
钻孔轨迹长度适应度函数：f2(x)＝length
[0125]
依据所确定的基础变量矩阵，选择适应度函数f1(x)利用irmo算法进行优化搜索，
以dx＝1m为计算精度，获得轨迹的最佳拟合区间；
[0126]
s2.1根据所确定的六个轨迹参数，确定每个轨迹参数的初始取值范围x
jmax
和x
jmin，
其中，j＝nod＝6；
[0127]
s2.2在步骤s2.1的取值范围内随机生成nop个初始位置点，由nop个初始位置点建立初始种群，所述初始位置点的数值信息通过计算式得出；
[0128]
x(i，j)＝x
jmin
+rand(0，1)(x
jmax-x
jmin
)
[0129]
其中，x(i，j)为生成的第i个初始位置点的第j个轨迹参数；rand(0，1)为在0到1之间的随机数；
[0130]
s2.3确定生成新的预位置点的更新条件，具体的如下，
[0131][0132][0133]
其中，y
i，jk
是指第k代新生成的第i个位置点的第j个参数变量；centre
jk
是指第k代中心位置的第j个参数变量；wk为惯性权值，随迭代次数递减，用于决定算法的收敛速度；k是指当前迭代次数；g是指最大迭代次数；
[0134]
中心位置随着当代最优位置的rbestxk和全局最优位置gbestxk的移动而移动，具体的，
[0135]
centre
k+1
＝centrek+0.4(gbestx
k-centrek)+0.5(rbestx
k-centrek),
[0136]
其中，centrek为第k代的中心位置；centre
k+1
为第k-1代的中心位置。
[0137]
s2.4采用更新条件生成新的预位置点y
ik
，并通过适应度函数f1(x)计算出各预位置点y
ik
对应的适应度值；
[0138]
s2.5更新位置信息，具体的，将第k代预位置点的适应度值fitness(y
ik
)与第k-1代位置点的适应度值fitness(x
ik-1
)进行比较，若fitness(y
ik
)《fitness(x
ik-1
)，则需要更新位置点信息，令fitness(x
ik
)＝fitness(y
ik
)，x
ik
＝y
ik
；否则，令fitness(x
ik
)＝fitness(x
ik-1
)，x
ik
＝x
ik-1
；其中，x
ik
为第k代第i个位置点；x
ik-1
为第k-1代第i个位置点；y
ik
为第k代生成的第i个预位置点。
[0139]
s2.6确定当代最优位置rbestxk和全局最优位置gbestxk，具体的，将更新后的第k代中的fitness(x
ik
)中的最小值的位置点作为当代最优位置rbestxk；将更新后的前k代中的fitness(x
ik
)中的最小值的位置点作为全局最优位置gbestxk；
[0140]
s2.7更新全局最优位置，具体的，比较当代最优位置的rbestxk的适应度值fitness(rbestxk)与全局最优位置gbestxk的适应度值fitness(gbestxk)；若fitness(rbestxk)《fitness(gbestxk)，则需要更新全局最优位置，令fitness(gbestxk)＝fitness(rbestxk)，gbestxk＝rbestxk；否则，令fitness(gbestxk)＝fitness(gbestx
k-1
)，gbestxk＝gbestx
k-1
。
[0141]
s3，将最佳拟合区间确定为轨迹长度优化的最佳优化区间，结合所述最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件并采用irmo算法在最佳优化区间内搜索求解轨迹长度多维非线性目标函数的最优值，得到最佳优化区间内的
钻进轨迹长度的最优钻进轨迹长度；
[0142]
在确定轨迹的最佳拟合区间的基础上，获得六个轨迹参数的最佳优化区间x
bjmax
和x
bjmin
，选择适应度函数f2(x)利用irmo算法对最佳优化区间内轨迹的最短长度进行优化搜索，其算法搜索步骤同s2.1-s2.7。
[0143]
本实施例使用所述最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件来验证钻孔轨迹设计的可行性。
[0144]
优化结果实例：
[0145]
采用irmo算法对“五段法”钻孔轨迹和悬链线轨迹拟合实例：
[0146]
表1悬链线轨迹与“五段法”轨迹参数设置
[0147][0148][0149]
不同算法参数对拟合结果的影响如下：
[0150]
表2不同算法参数计算所得拟合结果对比
[0151][0152]
注：r2值越逼近1表明轨迹的拟合程度越接近。
[0153]
连续20次的搜索拟合，整体sos的优化结果较为稳定(图6)，得到的拟合轨迹图像亦重叠度高(图7)，表明通过irmo算法计算的拟合结果不仅搜索精度高，且非常好的稳定性，证明了irmo算法的优越性。
[0154]
1.拟合基础上确定的参数最佳优化区间如下：
[0155]
表3调整后的轨迹参数约束范围
[0156][0157]
图8，展示了以钻进距离为优化目标的irmo算法连续20次搜索的最优钻孔轨迹。在满足钻孔稳定性的检验的条件下，钻进轨迹长度可以得到了显著的优化。
[0158]
选取了20次优化前后的钻孔轨迹长度对比，可以看出，优化结果非常稳定，标准差仅为0.05813，最大优化长度16.92m，最小优化长度10.62m，平均优化长度14.53m，优化效果十分显著，图9。
[0159]
优化前后的轨迹参数进行对比如下，表4：
[0160]
表4优化前后的轨迹参数
[0161][0162]
本实施例根据水平定向钻进轨迹的悬链线轨迹模型和“直线-圆弧-直线-圆弧-直线”五段式轨迹模型，在掌握工程要求和区域地质条件的情况下，分别建立拟合度和轨迹长度的求解方程，然后利用irmo算法将拟合度求解方程作为算法的适应度函数获得钻孔轨迹的最佳优化区间，进一步将钻孔轨迹长度作为算法的适应度函数在最佳优化区间内搜索求解轨迹长度最优值，同时加入了以泥浆压力理论为基础的轨迹可行性检验条件保证轨迹优化结果的合理性和可行性。该方法能够有效优化水平定向钻进轨迹的长度，精确计算出最优结果下的轨迹参数。
[0163]
实施例二
[0164]
参考图10，本实施例公开了一种二维水平定向钻进轨迹计算装置，其包括如下单元：
[0165]
模型建立单元，用于根据获取的水平定向钻进轨迹起点与终点的水平穿越长度a，纵向高差δz，利用悬链线钻孔轨迹模型建立所述水平定向钻进轨迹的第一轨迹方程，利用“五段法”钻孔轨迹模型建立所述水平定向钻进轨迹的第二轨迹方程并建立所述第二轨迹方程的约束条件；根据所述第一轨迹方程的纵坐标和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数；根据最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件；
[0166]
最佳拟合区间计算单元，用于将所述第二轨迹方程的轨迹参数转化为改进的径向移动算法irmo的基础变量矩阵，根据所述第一轨迹方程的纵坐标和所述第二轨迹方程的纵坐标zi作为拟合程度评估函数并结合irmo算法实现“五段法”钻孔轨迹与悬链线轨迹的最佳拟合，确定两者的最佳拟合区间；
[0167]
轨迹长度计算单元，用于将最佳拟合区间确定为轨迹长度优化的最佳优化区间，结合所述最大泥浆压力map和最小泥浆压力mdp建立所述第二轨迹方程的泥浆压力约束条件并采用irmo算法在最佳优化区间内搜索求解轨迹长度多维非线性目标函数的最优值，得到最佳优化区间内的钻进轨迹长度的最优钻进轨迹长度。
[0168]
具体的，所述“五段法”钻孔轨迹模型的六个轨迹参数如下：入口段直线段长度l1；入口段直线段倾角θ1；入口段过渡圆弧半径r2；中心轨迹直线段长度l3；出口段过渡圆弧半径r4；出口段直线段倾角θ2。
[0169]
具体的，所述最佳拟合区间计算单元还包括：
[0170]
第一计算单元，用于根据所确定的六个轨迹参数，确定每个轨迹参数的初始取值范围x
jmax
和x
jmin
，其中，j＝nod＝6；
[0171]
第二计算单元，用于在第一计算单元获取的取值范围内随机生成nop个初始位置点，由nop个初始位置点建立初始种群，所述初始位置点的数值信息通过计算式得出；
[0172]
x(i，j)＝x
jmin
+rand(0，1)(x
jmax-x
jmin
)
[0173]
其中，x(i，j)为生成的第i个初始位置点的第j个轨迹参数；rand(0，1)为在0到1之间的随机数；
[0174]
第三计算单元，用于确定生成新的预位置点的更新条件，具体的如下，
[0175][0176][0177]
其中，y
i，jk
是指第k代新生成的第i个位置点的第j个参数变量；centre
jk
是指第k代中心位置的第j个参数变量；wk为惯性权值，随迭代次数递减，用于决定算法的收敛速度；k是指当前迭代次数；g是指最大迭代次数；
[0178]
中心位置随着当代最优位置的rbestxk和全局最优位置gbestxk的移动而移动，具体的，
[0179]
centre
k+1
＝centrek+0.4(gbestx
k-centrek)+0.5(rbestx
k-centrek),
[0180]
其中，centrek为第k代的中心位置；centre
k+1
为第k-1代的中心位置；
[0181]
第四计算单元，用于采用更新条件生成新的预位置点y
ik
，并通过拟合程度评估函数f1(x)计算出各预位置点y
ik
对应的适应度值；
[0182]
第五计算单元，用于更新位置信息，具体的，将第k代预位置点的适应度值fitness(y
ik
)与第k-1代位置点的适应度值fitness(x
ik-1
)进行比较，若fitness(y
ik
)《fitness(x
ik-1
)，则需要更新位置点信息，令fitness(x
ik
)＝fitness(y
ik
)，x
ik
＝y
ik
；否则，令fitness(x
ik
)＝fitness(x
ik-1
)，x
ik
＝x
ik-1
；其中，x
ik
为第k代第i个位置点；x
ik-1
为第k-1代第i个位置点；y
ik
为第k代生成的第i个预位置点；
[0183]
第六计算单元，用于确定当代最优位置rbestxk和全局最优位置gbestxk，具体的，将更新后的第k代中的fitness(x
ik
)中的最小值的位置点作为当代最优位置rbestxk；将更新后的前k代中的fitness(x
ik
)中的最小值的位置点作为全局最优位置gbestxk；
[0184]
第七计算单元，用于更新全局最优位置，具体的，比较当代最优位置的rbestxk的适应度值fitness(rbestxk)与全局最优位置gbestxk的适应度值fitness(gbestxk)；若
fitness(rbestxk)《fitness(gbestxk)，则需要更新全局最优位置，令fitness(gbestxk)＝fitness(rbestxk)，gbestxk＝rbestxk；否则，令fitness(gbestxk)＝fitness(gbestx
k-1
)，gbestxk＝gbestx
k-1
。
[0185]
具体的，钻孔轨迹长度适应度函数f2(x)＝l1+r2·
θ1+l3+r4·
θ2+l5；
[0186][0187]
其中，l1为入口段直线段长度、l3为中心轨迹直线段长度、l5为出口段直线段长度，r2为入土圆弧段的半径、r4为出口段过渡圆弧半径，θ1为入口段直线段倾角、θ2为出口段直线段倾角。
[0188]
具体的，所述轨迹长度计算单元具体包括：
[0189]
获得六个轨迹参数的最佳优化区间x
bjmax
和x
bjmin
，选择钻孔轨迹长度适应度函数f2(x)代替第四计算单元中的拟合程度评估函数f1(x)，利用irmo算法对最佳优化区间内轨迹的最短长度进行优化搜索，其算法搜索步骤同第一计算单元至第七计算单元。
[0190]
本实施例根据水平定向钻进轨迹的悬链线轨迹模型和“直线-圆弧-直线-圆弧-直线”五段式轨迹模型，在掌握工程要求和区域地质条件的情况下，分别建立拟合度和轨迹长度的求解方程，然后利用irmo算法将拟合度求解方程作为算法的适应度函数获得钻孔轨迹的最佳优化区间，进一步将钻孔轨迹长度作为算法的适应度函数在最佳优化区间内搜索求解轨迹长度最优值，同时加入了以泥浆压力理论为基础的轨迹可行性检验条件保证轨迹优化结果的合理性和可行性。该方法能够有效优化水平定向钻进轨迹的长度，精确计算出最优结果下的轨迹参数。
[0191]
实施例三
[0192]
参考图11，图11是本实施例的一种二维水平定向钻进轨迹计算设备的结构示意图。该实施例的二维水平定向钻进轨迹计算设备20包括处理器21、存储器22以及存储在所述存储器22中并可在所述处理器21上运行的计算机程序。所述处理器21执行所述计算机程序时实现上述方法实施例中的步骤。或者，所述处理器21执行所述计算机程序时实现上述各装置实施例中各模块/单元的功能。
[0193]
示例性的，所述计算机程序可以被分割成一个或多个模块/单元，所述一个或者多个模块/单元被存储在所述存储器22中，并由所述处理器21执行，以完成本发明。所述一个或多个模块/单元可以是能够完成特定功能的一系列计算机程序指令段，该指令段用于描述所述计算机程序在所述二维水平定向钻进轨迹计算设备20中的执行过程。例如，所述计算机程序可以被分割成实施例二中的各个模块,各模块具体功能请参考上述实施例所述的装置的工作过程，在此不再赘述。
[0194]
所述二维水平定向钻进轨迹计算设备20可包括，但不仅限于，处理器21、存储器22。本领域技术人员可以理解，所述示意图仅仅是二维水平定向钻进轨迹计算设备20的示例，并不构成对二维水平定向钻进轨迹计算设备20的限定，可以包括比图示更多或更少的部件，或者组合某些部件，或者不同的部件，例如所述二维水平定向钻进轨迹计算设备20还可以包括输入输出设备、网络接入设备、总线等。
[0195]
所述处理器21可以是中央处理单元(central processing unit，cpu)，还可以是
其他通用处理器、数字信号处理器(digital signal processor，dsp)、专用集成电路(application specific integrated circuit，asic)、现成可编程门阵列(field-programmable gate array，fpga)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等，所述处理器21是所述二维水平定向钻进轨迹计算设备20的控制中心，利用各种接口和线路连接整个二维水平定向钻进轨迹计算设备20的各个部分。
[0196]
所述存储器22可用于存储所述计算机程序和/或模块，所述处理器21通过运行或执行存储在所述存储器22内的计算机程序和/或模块，以及调用存储在存储器22内的数据，实现所述二维水平定向钻进轨迹计算设备20的各种功能。所述存储器22可主要包括存储程序区和存储数据区，其中，存储程序区可存储操作系统、至少一个功能所需的应用程序(比如声音播放功能、图像播放功能等)等；存储数据区可存储根据手机的使用所创建的数据(比如音频数据、电话本等)等。此外，存储器22可以包括高速随机存取存储器，还可以包括非易失性存储器，例如硬盘、内存、插接式硬盘，智能存储卡(smart media card,smc)，安全数字(secure digital,sd)卡，闪存卡(flash card)、至少一个磁盘存储器件、闪存器件、或其他易失性固态存储器件。
[0197]
其中，所述二维水平定向钻进轨迹计算设备20集成的模块/单元如果以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明实现上述实施例方法中的全部或部分流程，也可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一计算机可读存储介质中，该计算机程序在被处理器21执行时，可实现上述各个方法实施例的步骤。其中，所述计算机程序包括计算机程序代码，所述计算机程序代码可以为源代码形式、对象代码形式、可执行文件或某些中间形式等。所述计算机可读介质可以包括：能够携带所述计算机程序代码的任何实体或装置、记录介质、u盘、移动硬盘、磁碟、光盘、计算机存储器、只读存储器(rom，read-only memory)、随机存取存储器(ram，random access memory)、电载波信号、电信信号以及软件分发介质等。需要说明的是，所述计算机可读介质包含的内容可以根据司法管辖区内立法和专利实践的要求进行适当的增减，例如在某些司法管辖区，根据立法和专利实践，计算机可读介质不包括电载波信号和电信信号。
[0198]
需说明的是，以上所描述的装置实施例仅仅是示意性的，其中所述作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。另外，本发明提供的装置实施例附图中，模块之间的连接关系表示它们之间具有通信连接，具体可以实现为一条或多条通信总线或信号线。本领域普通技术人员在不付出创造性劳动的情况下，即可以理解并实施。
[0199]
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。