33,311
社区成员
发帖
与我相关
我的任务
分享
#include <stdio.h>
#include <math.h>
#define pi 3.14159265
#define wv(k) sin(pi*(k*dt/(g*dt))*360/180)
double u[200][100]={0};
double w[200][100]={0};
void main()
{
FILE *p;
p=fopen("huang.txt","w");
int a,b,h,t;
a=200;
b=300;
h=50;
t=2*h/a;
/*grid space*/
/*
dx=1;%m
XM=200;%m
dy=1;%m
YM=100;%m*/
/*wavelet*/
int f=50;
double dt;
/*T=1/f=0.02*/
dt=0.0001;
int g;
/*g=T/dt*/
g=200;
int wfxi,wfyi,j,k,m,n;
/*source position*/
int sx,sy;
sy=0;
sx=100;
/*wave field 1*/
double t1,d,x[200],y[50];
int ix[200],iy[50],i;
for (k=0;k<100;k++)
{
t1=t*double(k)/100;
d=t1*a;
for(i=int(sx-d);i<=sx+d;i++)
{
x[i]=double(i);
ix[i]=int(x[i]);
y[i]=double(sqrt(d*d-(i-sx)*(i-sx)));
iy[i]=int(y[i]);
}
for (i=0;i<(2*d);i++)
{
wfxi=ix[i];
wfyi=iy[i];
for (j=0;j<4;j++)
{
if ((wfyi-j+1)>=1)
u[wfxi][wfyi-j+1]=double(wv(j));
}
}
for (i=0;i<200;i=i++)
{
w[i][k]=u[i][0];
}
}
//wave field 2
double t2,d2;
for (k=50;k<100;k++)
{
t2=t/2+t*double(k)/100;
d2=(t2-t/2)*a;
for(i=int(sx-d);i<=(sx+d);i++)
{
x[i]=double(i);
ix[i]=int(x[i]);
y[i]=double(h-sqrt(d2*d2-(i-sx)*(i-sx)));
iy[i]=int(y[i]);
}
for (i=0;i<(2*d2);i++)
{
wfxi=ix[i];
wfyi=iy[i];
for (j=0;j<4;j++)
{
if ((wfyi-j+1)>=1)
u[wfxi][wfyi-j+1]=-double(wv(j));
}
}
for (i=0;i<200;i=i++)
{
w[i][k]=u[i][0];
}
}
for(m=0;m<100;m++)
{
for(n=0;n<200;n++)
{
fprintf(p,"%f ",w[n][m]);
}
fprintf(p,"\n");
}
fclose(p);
}