联剖装置球体 异常的计算

《联剖装置球体 异常的计算》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、课程设计报告——联剖装置球体异常的计算专业：地球物理学姓名：学号：指导教师：一、二、实验C++程序：#include8#include#definepi3.1416doublefun(doublen,doublex){doublePn[1000];Pn[0]=1;Pn[1]=x;if((int)n>=2)Pn[(int)n]=(2*n-1)/n*Pn[(int)n-1]-(n-1)/n*Pn[(int)n-2];//Pn(y)计算的递推公式,3—3公式//returnPn[(int)n];}in

2、tmain(){FILE*fp;inti;doublek,x,n1,n2,n3,n4,t1,t2,t3,t4,temp1,temp2,temp3,temp4,temp5,temp6,temp7,temp8,da,db,rm,rn,Ao,oB,Mo,No,AM=1.5,AN=2.5,MN=1.0,r0=1.0,h=1.5,I=1.0,p1=1.0,p2=0.005,pa[51],pb[51],UMa[51],UMb[51],UNa[51],UNb[51],UMNa[51],UMNb[51];fp=fopen("data.txt",

3、"w+");//定义每个参数的初值//k=2*pi*AM*AN/MN;//计算装置系数//for(i=0,x=-25.0;x<26.0;x++,i++)//定义测线长度及步长//{t1=0.0;n1=0.0;t2=0.0;n2=0.0;t3=0.0;n3=0.0;t4=0.0;n4=0.0;No=fabs(x+MN/2);Mo=fabs(x-MN/2);if(x<0.0){Ao=AN+No;oB=AM-No;}if(x>=0.0){Ao=AN-No;oB=AM+No;}//求测线上供电电极、测量电极和球体中心投影的关系//da=

4、pow(pow(Ao,2)+pow(h,2),0.5);db=pow(pow(oB,2)+pow(h,2),0.5);rm=pow(pow(Mo,2)+pow(h,2),0.5);//求低阻体到供电电极和测量电极的距离//rn=pow(pow(No,2)+pow(h,2),0.5);//printf("%lf",fabs(-5.0))//do{temp1=(p2-p1)*n1/(p1*n1+p2*(n1+1))*pow(r0,2*n1+1)/(pow(da,n1+1)*pow(rm,n1+1))*fun(n1,x);temp

5、2=(p2-p1)*(n1+1)/(p1*(n1+1)+p2*((n1+1)+1))*pow(r0,2*(n1+1)+1)/(pow(da,(n1+1)+1)*pow(rm,(n1+1)+1))*fun(n1+1,x);t1+=temp1;n1++;}while((temp1-temp2)>0.01);//由子程序求电位U，第二项开始进行累加，n值以第n项和第n+1项之差小于1%//do{temp3=(p2-p1)*n2/(p1*n2+p2*(n2+1))*pow(r0,2*n2+1)/(pow(db,n2+1)*pow(rm,

6、n2+1))*fun(n2,x);temp4=(p2-p1)*(n2+1)/(p1*(n2+1)+p2*((n2+1)+1))*pow(r0,2*(n2+1)+1)/(pow(db,(n2+1)+1)*pow(rm,(n2+1)+1))*fun(n2+1,x);t2+=temp3;n2++;}while((temp3-temp4)>0.01);do{temp5=(p2-p1)*n3/(p1*n3+p2*(n3+1))*pow(r0,2*n3+1)/(pow(da,n3+1)*pow(rn,n3+1))*fun(n3,x);tem

7、p6=(p2-p1)*(n3+1)/(p1*(n3+1)+p2*((n3+1)+1))*pow(r0,2*(n3+1)+1)/(pow(da,(n3+1)+1)*pow(rn,(n3+1)+1))*fun(n3+1,x);t3+=temp5;n3++;}while((temp5-temp6)>0.01);do{temp7=(p2-p1)*n4/(p1*n4+p2*(n4+1))*pow(r0,2*n4+1)/(pow(db,n4+1)*pow(rn,n4+1))*fun(n4,x);temp8=(p2-p1)*(n4+1)/(p

8、1*(n4+1)+p2*((n4+1)+1))*pow(r0,2*(n4+1)+1)/(pow(db,(n4+1)+1)*pow(rn,(n4+1)+1))*fun(n4+1,x);t4+=temp7;n4++;}while((temp7-temp8)>0.01