《计算地球物理学.docx》由会员分享,可在线阅读,更多相关《计算地球物理学.docx(5页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、计算地球物理学计算地球物理学一、双曲线形方程算法和程序:在如下问题中,对下列给定定值,用程序求解波动方程uxx,t=c2uxx(x,t),其中0xa且0tb,边界条件为:u(0,t)=0且u(a,t)=0,0tbu(x,0)=f(x),0xa(x,0)=g(x),0xaux用surf和contour命令画图得到近似值解。1.设a=1,b=1,c=1,f(x)=sin(x),g(x)=0。为了方便起见,选择h=0.1,k=0.1。2.设a=1,b=1,c=2,f(x)=x-x2,g(x)=0。为了方便起见,选择h=0.1,k=0.05。解:程序代码:functionu,r,x,y=finedif
2、(f,g,a,b,c,h,k)%finedion波动方程的差分方法程序%f:初始条件方程,字符型sring;%g:边界条件方程,字符型sring;%a:位置x的上限0,a;%b:时间t的上限0,b;%c:方程系数;%h:x的剖分步长;%k:t的剖分步长;n=a/h+1;m=b/k+1;r=c*k/h;r2=r2;r22=r2/2;s1=1-r2;s2=2-2*r2;U=zeros(n,m);%赋值边界条件fori=2:n-1U(i,1)=feval(f,h*(i-1);U(i,2)=s1*feval(f,h*(i-1)+k*feval(g,h*(i-1)+r22*(feval(f,h*i)+feval(f,h*(i-2);end%求取个点数值forj=3:mfori=2:(n-1)U(i,j)=s2*U(i,j-1)+r2*(U(i-1,j-1)+U(i+1,j-1)-U(i,j-2);endendu=U%坐标量展示:x=0:h:a;y=0:k:bend问题1:稳定性条件分析与运算结果:r=1结果稳定结果图展示:二、抛物型方程的算法和程序:求解热传导方程:ut(x,t)=c2uxx(x,t),其中0