Commit 672a5b8a authored by Florian Wittkamp's avatar Florian Wittkamp

Adaptions for better comparison to WAVE code

parent 2c7d3af5
...@@ -29,7 +29,6 @@ float ** srcpos_loc, float ** signals, int nsrc, int * stype){ ...@@ -29,7 +29,6 @@ float ** srcpos_loc, float ** signals, int nsrc, int * stype){
extern float DX, DY, DZ; extern float DX, DY, DZ;
//extern float DT; //extern float DT;
extern int NT;
int i, j, k, l; int i, j, k, l;
float amp=0.0; float amp=0.0;
...@@ -51,9 +50,7 @@ float ** srcpos_loc, float ** signals, int nsrc, int * stype){ ...@@ -51,9 +50,7 @@ float ** srcpos_loc, float ** signals, int nsrc, int * stype){
//scaled explosive source with respect to spatial and temporal discretization, seismic Moment = 1 Nm //scaled explosive source with respect to spatial and temporal discretization, seismic Moment = 1 Nm
// -> additional temporal derivation // -> additional temporal derivation
if(nt==1){amp=signals[l][nt+1]/(2.0*DX*DY*DZ);} amp=signals[l][nt]/(2.0*DX*DY*DZ);
if((nt>1)&&(nt<NT)){amp=(signals[l][nt+1]-signals[l][nt-1])/(2.0*DX*DY*DZ);}
if(nt==NT){amp=-signals[l][nt-1]/(2.0*DX*DY*DZ);}
sxx[j][i][k]+=amp; sxx[j][i][k]+=amp;
syy[j][i][k]+=amp; syy[j][i][k]+=amp;
......
...@@ -945,6 +945,16 @@ int main(int argc, char **argv){ ...@@ -945,6 +945,16 @@ int main(int argc, char **argv){
* eqsource is a implementation of moment tensor points sources. */ * eqsource is a implementation of moment tensor points sources. */
eqsource(nt,sxx,syy,szz,sxy, syz, sxz, srcpos_loc,signals,nsrc_loc,stype_loc, amon, str, dip, rake); eqsource(nt,sxx,syy,szz,sxy, syz, sxz, srcpos_loc,signals,nsrc_loc,stype_loc, amon, str, dip, rake);
} }
/* store amplitudes at receivers in e.g. sectionvx, sectionvz, sectiondiv, ...*/
if ((SEISMO) && (ntr>0) && (nt==lsamp)){
seismo(nlsamp,ntr,recpos_loc,sectionvx,sectionvy,sectionvz,
sectiondiv,sectioncurl,sectionp,vx,vy,vz,sxx,syy,szz,pi,u);
nlsamp++;
lsamp+=NDT;
}
/* stress free surface ? */ /* stress free surface ? */
if ((FREE_SURF) && (POS[2]==0)){ if ((FREE_SURF) && (POS[2]==0)){
...@@ -958,14 +968,6 @@ int main(int argc, char **argv){ ...@@ -958,14 +968,6 @@ int main(int argc, char **argv){
time_s_exchange[nt]=exchange_s(nt,sxx,syy,szz,sxy,syz,sxz,sbufferlef_to_rig,sbufferrig_to_lef,sbuffertop_to_bot,sbufferbot_to_top,sbufferfro_to_bac,sbufferbac_to_fro, sreq_send, sreq_rec); time_s_exchange[nt]=exchange_s(nt,sxx,syy,szz,sxy,syz,sxz,sbufferlef_to_rig,sbufferrig_to_lef,sbuffertop_to_bot,sbufferbot_to_top,sbufferfro_to_bac,sbufferbac_to_fro, sreq_send, sreq_rec);
/* store amplitudes at receivers in e.g. sectionvx, sectionvz, sectiondiv, ...*/
if ((SEISMO) && (ntr>0) && (nt==lsamp)){
seismo(nlsamp,ntr,recpos_loc,sectionvx,sectionvy,sectionvz,
sectiondiv,sectioncurl,sectionp,vx,vy,vz,sxx,syy,szz,pi,u);
nlsamp++;
lsamp+=NDT;
}
/* save snapshot in file */ /* save snapshot in file */
if ((SNAP) && (nt==lsnap) && (nt<=TSNAP2/DT)){ if ((SNAP) && (nt==lsnap) && (nt<=TSNAP2/DT)){
......
...@@ -53,7 +53,7 @@ float ** wavelet(float ** srcpos_loc, int nsrc){ ...@@ -53,7 +53,7 @@ float ** wavelet(float ** srcpos_loc, int nsrc){
signals=fmatrix(1,nsrc,1,NT); signals=fmatrix(1,nsrc,1,NT);
for (nt=1;nt<=NT;nt++){ for (nt=1;nt<=NT;nt++){
t=(float)nt*DT; t=(float)(nt-1)*DT;
for (k=1;k<=nsrc;k++) { for (k=1;k<=nsrc;k++) {
tshift=srcpos_loc[4][k]; tshift=srcpos_loc[4][k];
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment