(整理)二维波动方程第一类吸收边界条件c++实现代码.

合集下载

相关主题

1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。

精品文档

#include "stdafx.h"

#include

using namespace std;

const double pi=4*atan(1.0);

double freq=45;

double sb=7.45;

double t1=2*pi/(sb*4);

double source(double t)

{

//double t2=0.0;

if(t<=t1) return (sin(sb*4*t-pi/2)+1)/10;

else{ double tep=0.0; return tep;}

//return ((1-2*pi*pi*freq*freq*t*t)*exp(-pi*pi*freq*freq*t*t)+1);//Ricker子波}

void update_Vn(double upt,double lowt,double upx1,double lowx1)

{

int i,j,m;

const int Csize=300;

double deg=0;

double stepx1=abs(upx1-lowx1)/(Csize-1);

//double te=sqrt(static_cast(3.0/8.0));

double stept=sqrt(static_cast(1.0/2.0))*stepx1/2.0;//

int tn=static_cast(upt/stept);

double r=stept/stepx1;

double **u_current,**u_old,**u_past;

u_current=new double *[Csize];

u_old=new double*[Csize];

u_past=new double*[Csize];

for(i=0;i

{

u_current[i]=new double [Csize];

u_old[i]=new double[Csize];

u_past[i]=new double[Csize];

}

for(i=0;i

for(j=0;j

{

u_current[i][j]=0;

u_old[i][j]=0;

精品文档

u_past[i][j]=0;

}

double ck[Csize][Csize];

int flag=0;

for(j=0;j

{

for(i=0;i

{

if(j

else ck[i][j]=1;

}

string str;

cout<<"\n 输入保存文件名:";

cin>>str;

ofstream fout(str);

if(!fout)

{

cout<<"\n 不能打开文件"<

}

m=0;

double f0=2.0/(stept*30);

double t0=4.0/f0;

while(m<1500&&((m++)*stept+lowt)

{

for(i=0;i

for(j=0;j

{

if((i!=0&&(i!=Csize-1))&&(j!=0&&j!=(Csize-1)))//cope with the internal points of the interesting domain

{

//if(i==(Csize/2)&&j==(Csize/2))

u_current[i][j]=(r*r/ck[i][j])*(u_old[i+1][j]+u_old[i][j+1]+u_old[i-1][j]+u_old[i][j-1])-2 *u_old[i][j]*(2*r*r/ck[i][j]-1.0)-u_past[i][j]

//+stept*stept*source(m*stept+lowt);//stept*stept*exp(-f0*f0*(m*stept-t0)*(m*stept-t0 ));//stept*stept*source(m*stept+lowt);

//elseu_current[i][j]=(r*r/ck[i][j])*(u_old[i+1][j]+u_old[i][j+1]+u_old[i-1][j]+u_old[i][j -1])-2*u_old[i][j]*(2*r*r/ck[i][j]-1.0)-u_past[i][j];//u[m][i][j]=r*r*ck[i][j]*ck[i][j]*(u[m -1][i+1][j]+u[m-1][i][j+1]+u[m-1][i][j-1]+u[m-1][i-1][j])-2*(2*r*r*ck[i][j]*ck[i][j]-1)* u[m-1][i][j]

//-u[m-2][i][j];