Skip to content
Snippets Groups Projects
Commit 08b12eaf authored by Jan Thorbecke's avatar Jan Thorbecke
Browse files

fixed OpenMP bug in viscoacoustic4.c for Threads>1

parent 8624790f
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
fdelmodc/demo/fdelmodc_visco.scr
+ 1
− 0
View file @ 08b12eaf
...@@ -37,6 +37,7 @@ export filero=model_ro.su ...@@ -37,6 +37,7 @@ export filero=model_ro.su
export fileqp=relax_cp.su export fileqp=relax_cp.su
export fileqs=relax_cs.su export fileqs=relax_cs.su
export OMP_NUM_THREADS=4
../fdelmodc \ ../fdelmodc \
file_cp=$filecp file_cs=$filecs file_den=$filero \ file_cp=$filecp file_cs=$filecs file_den=$filero \
......
fdelmodc/fdelmodc.c
+ 1
− 3
View file @ 08b12eaf
...@@ -40,9 +40,7 @@ int acoustic2(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixs ...@@ -40,9 +40,7 @@ int acoustic2(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixs
int acoustic4Block(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixsrc, int izsrc, float **src_nwav, float *vx, int acoustic4Block(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixsrc, int izsrc, float **src_nwav, float *vx,
float *vz, float *p, float *rox, float *roz, float *l2m, int verbose); float *vz, float *p, float *rox, float *roz, float *l2m, int verbose);
int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixsrc, int izsrc, float **src_nwav, float *vx, float *vz, float *p, float *rox, float *roz, float *l2m, float *tss, float *tep, float *q, int verbose);
int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixsrc, int izsrc, float **src_nwav, float *vx, float
*vz, float *p, float *rox, float *roz, float *l2m, float *tss, float *tep, float *q, int verbose);
int elastic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixsrc, int izsrc, float **src_nwav, float *vx, float *vz, float *tzz, float *txx, float *txz, float *rox, float *roz, float *l2m, float *lam, float *mul, int verbose); int elastic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, int ixsrc, int izsrc, float **src_nwav, float *vx, float *vz, float *tzz, float *txx, float *txz, float *rox, float *roz, float *l2m, float *lam, float *mul, int verbose);
......
fdelmodc/viscoacoustic4.c
+ 1
− 7
View file @ 08b12eaf
...@@ -112,14 +112,13 @@ int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, in ...@@ -112,14 +112,13 @@ int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, in
boundariesP(mod, bnd, vx, vz, p, NULL, NULL, rox, roz, l2m, NULL, NULL, itime, verbose); boundariesP(mod, bnd, vx, vz, p, NULL, NULL, rox, roz, l2m, NULL, NULL, itime, verbose);
/* calculate p/tzz for all grid points except on the virtual boundary */ /* calculate p/tzz for all grid points except on the virtual boundary */
#pragma omp for private (iz,ix, Tpp) nowait schedule(guided,1)
for (ix=mod.ioPx; ix<mod.iePx; ix++) { for (ix=mod.ioPx; ix<mod.iePx; ix++) {
#pragma omp for private (iz) nowait schedule(guided,1)
#pragma ivdep #pragma ivdep
for (iz=mod.ioPz; iz<mod.iePz; iz++) { for (iz=mod.ioPz; iz<mod.iePz; iz++) {
dxvx[iz] = c1*(vx[(ix+1)*n1+iz] - vx[ix*n1+iz]) + dxvx[iz] = c1*(vx[(ix+1)*n1+iz] - vx[ix*n1+iz]) +
c2*(vx[(ix+2)*n1+iz] - vx[(ix-1)*n1+iz]); c2*(vx[(ix+2)*n1+iz] - vx[(ix-1)*n1+iz]);
} }
#pragma omp for private (iz) nowait schedule(guided,1)
#pragma ivdep #pragma ivdep
for (iz=mod.ioPz; iz<mod.iePz; iz++) { for (iz=mod.ioPz; iz<mod.iePz; iz++) {
dzvz[iz] = c1*(vz[ix*n1+iz+1] - vz[ix*n1+iz]) + dzvz[iz] = c1*(vz[ix*n1+iz+1] - vz[ix*n1+iz]) +
...@@ -127,14 +126,12 @@ int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, in ...@@ -127,14 +126,12 @@ int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, in
} }
/* help variables to let the compiler vectorize the loops */ /* help variables to let the compiler vectorize the loops */
#pragma omp for private (iz) nowait schedule(guided,1)
#pragma ivdep #pragma ivdep
for (iz=mod.ioPz; iz<mod.iePz; iz++) { for (iz=mod.ioPz; iz<mod.iePz; iz++) {
Tpp = tep[ix*n1+iz]*tss[ix*n1+iz]; Tpp = tep[ix*n1+iz]*tss[ix*n1+iz];
Tlm[iz] = (1.0-Tpp)*tss[ix*n1+iz]*l2m[ix*n1+iz]*0.5; Tlm[iz] = (1.0-Tpp)*tss[ix*n1+iz]*l2m[ix*n1+iz]*0.5;
Tlp[iz] = l2m[ix*n1+iz]*Tpp; Tlp[iz] = l2m[ix*n1+iz]*Tpp;
} }
#pragma omp for private (iz) schedule(guided,1)
#pragma ivdep #pragma ivdep
for (iz=mod.ioPz; iz<mod.iePz; iz++) { for (iz=mod.ioPz; iz<mod.iePz; iz++) {
Tt1[iz] = 1.0/(ddt+0.5*tss[ix*n1+iz]); Tt1[iz] = 1.0/(ddt+0.5*tss[ix*n1+iz]);
...@@ -142,18 +139,15 @@ int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, in ...@@ -142,18 +139,15 @@ int viscoacoustic4(modPar mod, srcPar src, wavPar wav, bndPar bnd, int itime, in
} }
/* the update with the relaxation correction */ /* the update with the relaxation correction */
#pragma omp for private (iz) schedule(guided,1)
#pragma ivdep #pragma ivdep
for (iz=mod.ioPz; iz<mod.iePz; iz++) { for (iz=mod.ioPz; iz<mod.iePz; iz++) {
p[ix*n1+iz] -= Tlp[iz]*(dzvz[iz]+dxvx[iz]) + q[ix*n1+iz]; p[ix*n1+iz] -= Tlp[iz]*(dzvz[iz]+dxvx[iz]) + q[ix*n1+iz];
} }
#pragma omp for private (iz) schedule(guided,1)
#pragma ivdep #pragma ivdep
for (iz=mod.ioPz; iz<mod.iePz; iz++) { for (iz=mod.ioPz; iz<mod.iePz; iz++) {
q[ix*n1+iz] = (Tt2[iz]*q[ix*n1+iz] + Tlm[iz]*(dxvx[iz]+dzvz[iz]))*Tt1[iz]; q[ix*n1+iz] = (Tt2[iz]*q[ix*n1+iz] + Tlm[iz]*(dxvx[iz]+dzvz[iz]))*Tt1[iz];
p[ix*n1+iz] -= q[ix*n1+iz]; p[ix*n1+iz] -= q[ix*n1+iz];
} }
} }
/* Add stress source */ /* Add stress source */
......
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