地下球体高阻体电测深曲线计算FORTRAN实现程序及绘图结果.docxVIP

高阻球体的电测深曲线计算计算公式采用三级装置的计算公式：将MN尽量取小，与AB相比可以看成是无穷小。计算完、之后，对称四极装置的视电阻率按照下面公式计算： 计算完毕后，将生成的数据文件分别用Grapher和Surfer软件绘图，结果如下。球体上的电测深曲线，ab沿x方向移动，球体为高阻体等视电阻率断面图，ab沿x方向移动球体的电测深曲线，mn中点沿y方向移动，球体为高阻体等视电阻率断面图，mn中点沿y方向移动。计算所用的fortran程序programdianceshenreal::mn=0.1 ,pp1=10,pp2=100 ,r0=1 ,h0=2!mn长度,介质1电阻率，介质2电阻率，球体半径,球心埋深real psa,psb,psx,psy,u12,k,da,db,rm,rn,r!参数分别代表：A极供电时视电阻率，B极供电时视电阻率，x,y方向四极装置视电阻率，!介质1、2的视电阻率比，三级装置的装置系数，a、b到球心的距离，m、n到球心的距离，am长度real cos1,cos2,cos3,cos4!四个角度的余弦值real u,p1,p2,p3,p4!公式中的替代项real xx,xy,ax,ay,bx,by,mx,my,nx,ny!mn中点x,y坐标，a点x,y坐标，b点x,y坐标,m点x,y坐标，n点x,y坐标integeri !循环变量character*3 filenamex(7),filenamey(7)data filenamex/xa,xb,xc,xd,xe,xf,xg/xy=0;xx=0;lus=0!求视电阻率曲线数据（mn中点沿x方向移动）doi=1,7,1xx=xx+0.5open(i,file=filenamex(i)//jieguo.dat)do r=6,0.1,-0.01 call psxy(u12,u,xx,xy,r,mn,pp1,pp2,r0,h0,psx)write(i,(1x,f10.5,1x,f10.5))psx/pp1,renddoclose(i)end do!求等视电阻率断面图数据（mn中点沿x方向移动）open(10,file=xduanmian.dat)do xx=-3,3,0.01do r=3,0.1,-0.01 call psxy(u12,u,xx,xy,r,mn,pp1,pp2,r0,h0,psx) write(10,(1x,f10.5,1x,f10.5,1x,f10.5))psx/pp1,r,xxend do end doclose(10) !求视电阻率曲线数据（mn中点沿y方向移动）datafilenamey/ya,yb,yc,yd,ye,yf,yg/xx=0;xy=0doi=1,7,1xy=xy+0.5open(i,file=filenamey(i)//jieguo.dat)do r=6,0.1,-0.01 call psxy(u12,u,xx,xy,r,mn,pp1,pp2,r0,h0,psy)write(i,(1x,f10.5,1x,f10.5))psy/pp1,rend doclose(i)end do!求等视电阻率断面图数据（mn中点沿y方向移动）open(20,file=yduanmian.dat)doxy=-3,3,0.01do r=3,0.1,-0.01 call psxy(u12,u,xx,xy,r,mn,pp1,pp2,r0,h0,psy) write(20,(1x,f10.5,1x,f10.5,1x,f10.5))psy/pp1,r,xyend do end doclose(20) end program dianceshensubroutine psxy(u12,u,xx,xy,r,mn,pp1,pp2,r0,h0,ps)real psa,psb,ps,u12,k,da,db,rm,rn,rreal xx,xy,ax,ay,bx,by,mx,my,nx,nyreal mn,pp1,pp2,r0,h0real u,p1,p2,p3,p4u12=pp2/pp1u=(u12-1)/(2*u12+1)k=(r+mn/2)*(r+mn/2)/mnax=xx-(r+mn/2);ay=xybx=xx+(r+mn/2);by=xymx=xx-mn/2;my=xynx=xx+mn/2;ny=xyda=sqrt(h0*h0+ax*ax+ay*ay)db=sqrt(h0*h0+bx*bx+by*by)rm=sqrt(mx*mx+my*my+h0*h0)rn=sqrt(nx*nx+ny*ny+h0*h0)call p(da,rm,r,cos1,p1)call p(da,rn,r+mn,cos2,p2)call p(db,rm,r+mn,cos