答---钻井布局问题研究1

钻井布局问题研究摘要本题是研究勘探部门钻井布局问题，我们通过建立最优化模型来制定钻井的计划，使在最大利用旧井的前提下，尽可能的少打新井，来节约勘探费用。

在进行第一问求解时，我们利用旧井尽可能大为目标，建立目标规划的最优化模型。

通过对对假定网格的横向距离及纵向距离的最大值不起过给定误差值W(0.05单位)，在平面上平等移动网格N使可利用的旧井数尽可能大，利用matlab 编程求解得到，最大可利用旧井数为4个，它们的位置及坐标分别为之p2(1.41,3.50)p4(3.37,3.51)p5(3.40,5.50)p10(8.38,4.50)。

用验证求解的正确LINGO性。

在进行第二问求解中，要使旧井尽可能大为目标，建立目标规划的最优化模型。

因为网格是旋转移动，则用LINGO软件无法计算求解,只能用MATLAB编程求解。

关键词：目标规划非线性规划模型共同约束条件0—1变量一、问题的提出问题背景：勘探部门在某地区找矿。

初步勘探时期已零散地在若干位置上钻井，取得了地质资料。

进入系统勘探时期后，要在一个区域内按纵横等距的网格点来布置井位，进行“撒网式”全面钻探。

由于钻一口井的费用很高，如果新设计的井位与原有井位重合（或相当接近），便可利用旧井的地质资料，不必打这口新井。

因此，应该尽量利用旧井，少打新井，以节约钻探费用。

比如钻一口新井的费用为500万元，利用旧井资料的费用为10万元，则利用一口旧井就节约费用490万元。

设平面上有n个点pi其坐标为（ai，bi）,i=1，2…n表示已经有的n个井位。

新布置的井位是一个正方形网格N的所有结点（所谓“正方形网格”是指每个格子都是正方形的网格；结点是指纵线和横线的交叉点）。

假定每个格子的边长（井位的纵横间距）都是一单位（比如100米）整个网络是可以在平面上任意移动的。

若一个以知点Pi与某个网格结点Xi的距离不超过给定误差W（0.05单位），则认为Pi处的旧井资料可以利用，不必在结点Xi处打新井。

为进行辅助决策，勘探部门要求我们研究如下问题：1）假定网格的横向和纵向是固定的（比如东西向和南北向）并规定两点间的距离为横向距离（横坐标之差绝对值）及纵向距离（纵坐标之差绝对值）的最大值。

在平面上平行移动网格N使可利用的旧井数尽可能大。

试提供数值计算方法，并对下面的数值例子用计算机进行计算。

2）在欧式距离的误差意义下，考虑网格的横向和纵向是固定的情形，给出算法二、问题分析问题一分析题目要求我们在假定网格的横向和纵向是固定的（比如东西向和南北向）并规定两点间的距离为横向距离（横坐标之差绝对值）及纵向距离（纵坐标之差绝对值）的最大值。

在平面上平行移动网格N使可利用的旧井数尽可能大。

通过对目标规划建立非线性规划模型。

由于网格是固定的，整个网络是可以在平面上任意移动。

当一个以知点Pi 与某个网格结点X i的距离不超过给定给定误差W（0.05单位），则认为Pi处的旧井资料可以利用，不必在结点X i处打新井。

在用MATLAB编程软件画出给出的数值例子网格，通过循环结构，取整函数，求绝对值等利用编程找出可利用旧井的最大数问题二分析在欧式距离的误差意义下，考虑网格的横向和纵向是固定的，在这个特定的区域中，可以将网格平移也可以网格任意旋转，总之不管怎样移动或旋转，都要使网格在平面上移动网格N使可利用的旧井数尽可能大。

利用编程软件对其编程求解，要使旧井尽可能大为目标，建立目标规划的最优化模型。

因为网格是旋转移动，则用LINGO软件无法计算求解,只能用MATLAB编程求解。

三、模型假设根据问题分析对模型进行如下假设： 1.不考虑地理环境的影响。

2.假设网格的横向和纵向是固定的且足够的大3.假设网格在平移的过程只在水平方向和竖直方向移动可将旧井覆盖。

4．假设已知点P i 与某个网络结点X i 的距离不超过给定误差W （0.05单位），则认为P i 处的旧井资料可以利用，不必在结点X i 处打新井。

四、符号说明｛Z ｝ 取整P i （a i ,b i ） p i 旧井坐标(a i ,b i )p i ’(a i ’,b i ’) 平移之后p i ’的坐标(a i ’,b i ’)ϕ旋转前点与原点的连线与x 轴构成的角度X i 网格结点 x 移动后i 点新的x 轴坐标 y 移动后i 点新的y 轴坐标可利用的旧井数 β 网格在平面内旋转的角度j d 旋转移动后第j 点的位置坐标的欧式距离 y d 旋转后第j 点的纵坐标与结点纵坐标的差 x d 旋转后第j 点的横坐标与结点横坐标的差五、模型的建立与求解5.1问题一5.1.1解题思路针对我们在平面上平行移动网格N 使可利用的旧井数尽可能大的问题，我们通过建立一维的平面直角坐标系并在直角坐标系中绘制以1为单位的正方形网格，所绘制的网格必须完全能够覆盖所有的旧井，在平行移动网格时每次以0.01为单位进行水平方向或竖直方向移动。

比较移动后的旧井坐标与网格在水平方向和竖直方向上的差值的绝对值和0.05个单位进行比较，符合题意的选取表示旧井可利用，不符合题意的表示旧井不可利用。

建立非线性规划模型，用MATLAB 编程软件编程求解。

模型一5.1.2模型的分析 【1】0-1变量的引入题目中给了12个坐标，要使网格在平移后旧井的可利用数最大，我们引入{}12...1i 1,0==i f0-1变量来选择旧井可利用。

()121i 1i 0 =⎩⎨⎧=i f i 口旧井可利用第旧井不可利用口第【2】目标确立旧井利用数最大，题目给出12个坐标即表示12个旧井分布，建立最优化模型，确立模型目标【3】约束条件1. 在水平方向上平移后的距离到网格结点的距离2. 在竖直方向上平移后的距离到网格结点的距离3. i 口旧井，引入0—1 5.1.3模型的建立根据上述分析建立模型如下：5.1.3模型的求解 ①目标一求解根据上述模型，利用MATLTB 软件编程求解（MATLAB 程序一）。

得到的最大旧井数：（如表1所示）表1∑==121i i f Z Max 12...105.0]5.0[=≤++-+i f y b y b i i i 12...1i i 0i 1==口不可利用的旧井第口可利用的旧井第i f ∑==121i ifZ Max {}⎪⎩⎪⎨⎧=∈=≤++-+=≤++-+12...1i 1,012...105.0]5.0[12...1i 05.0]5.0[s.t ii i i i i i f i f y b y b f x a x a 12...105.0]5.0[=≤++-+i f x a x a i i i附表1附表25.2问题二5.2.1解题思路在欧式距离的误差意义下，考虑网格的横向和纵向是固定的，在这个特定的区域中，可以将网格平移也可以网格任意旋转，总之不管怎样移动或旋转，都要使网格在平面上移动网格N 使可利用的旧井数尽可能大。

利用 编程软件对其编程求解，要使旧井尽可能大为目标，建立目标规划的最优化模型。

因为网格是旋转移动，则用LINGO 软件无法计算求解,只能用MATLAB 编程求解。

模型二5.2.2模型分析【1】0-1变量的引入题目中给了12个坐标，要使网格在平移后旧井的可利用数最大，我们引入0-1变量来选择旧井可利用。

()121i 1i 0 =⎩⎨⎧=i f i 口旧井可利用第旧井不可利用口第MATLAB【2】目标确立 旧井利用数最大，题目给出12个坐标即表示12个旧井分布，通过建立最优化模型，确立模型目标【3】约束条件1.网格在平面内旋转，以坐标点p i （a i ,b i ）为原点旋转，满足的条件为2.旋转前点与原点的连线与x 轴构成的角度3.网格在平面内旋转的角度4.旋转后第j 点的横坐标与结点横坐标的差5.旋转后第j 点的纵坐标与结点纵坐标的差6. 旋转移动后第j 点的位置坐标的欧式距离7.平移后可利用的第i 口旧井，引入0—1变量5.2.3模型的建立根据上述分析建立模型如下：∑==121i i f Z Max )121i 1i 0=⎩⎨⎧=i f i口旧井可利用第旧井不可利用口第∑==121i ifZ Max 12...105.0]5.0[:=≤++--i f a x a x x i i i i i 12...105.0]5.0[:=≤++--i f b y b y y i i i i i ⎭⎬⎫⎩⎨⎧∈+=2,0πθθϕβ]5.0)([)(+-=j f j f d x ]5.0)([)(+-=j l j l d y 05.0)]5.0)([)(()]5.0)([)((22≤+-++-=j l j l j f j f d j )/h(j))arccos(a(j =ϕ⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎧≤+-++-=+-=+-=⎭⎬⎫⎩⎨⎧∈+===≤++--=≤++--05.0)]5.0)([)(()]5.0)([)((]5.0)([)(]5.0)([)(2,0)/h(j))arccos(a(j 12 (105).0]5.0[12...105.0]5.0[22j y x i i i i i i i i i i j l j l j f j f d j l j l d j f j f d i f b y b y i f a x a x πθθϕβϕ5.2.4模型的求解 ①目标二求解根据上述模型，利用MATLTB 软件编程求解（MATLAB 程序二）。

得到的最大旧井数：（如表2）表2附表2-1附表2-2附录一附录一MATLAB(一)编程如下： x=[0.50 1.41 3.00 3.37 3.40 4.72 4.72 5.43 7.57 8.38 8.98 9.50]; y=[2.00 3.50 1.50 3.515.502.006.244.102.014.50e=0;for i=0:0.01:1for j=0:0.01:1;a=x+i;b=y+j;d=round(a);g=round(b);s=0;for k=1:1:12;if abs(a(k)-d(k))<=0.05 & abs(b(k)-g(k))<=0.05;s=s+1;p=abs(a(k)-d(k));q=abs(b(k)-g(k));c(k)=k;else c(k)=0;endendif s>e;e=s;m=a;n=b;t=c;endendendemntpqplot(m,n,'*');grid附录二利用编程验证答案的正确性。

LINGOsets:demand/1..12/:a,b,k;endsetsdata:a=0.50 1.41 3.00 3.37 3.40 4.72 4.72 5.43 7.57 8.388.98 9.50;b=2.00 3.50 1.50 3.51 5.50 2.00 6.24 4.10 2.01 4.50enddataMax=@sum(demand(h):k(h));@for(demand(h):k(h)*@abs(a(h)+x-@floor(a(h)+x+0.5))<=0.05); @for(demand(h):k(h)*@abs(b(h)+y-@floor(b(h)+y+0.5))<=0.05); @for(demand(h):@bin(k(h)));x<1;y<1;附录三MATLAB(二)编程如下：x=[0.50 1.41 3.00 3.37 3.40 4.72 4.72 5.43 7.57 8.38 8.98 9.50];y=[2.00 3.50 1.50 3.51 5.50 2.00 6.24 4.10 2.01 4.50 3.41 0.80];e=0;for i=1:12;a=x-x(i);b=y-y(i);for k=0:0.01:pi/2;s=0;for j=1:12;if abs(i-j)>0h(j)=sqrt((a(j))^2+(b(j))^2);Q(j)=acos(a(j)/h(j));f(j)=h(j)*cos(k+Q(j));l(j)=h(j)*sin(k+Q(j));elsef(j)=j;l(j)=j;endif sqrt((f(j)-round(f(j)))^2+(l(j)-round(l(j)))^2)<=0.05;s=s+1;c(j)=j;elsec(j)=0;endendif s>e;e=s;t=c;q=f;w=l;endendendteplot(q,w,'*')grid。