Skip to content
Snippets Groups Projects
Commit ede72a71 authored by Holger Obermaier's avatar Holger Obermaier Committed by Holger Obermaier
Browse files

Add nontemporal hints for stencil of higher order

parent a4f6f39d
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
src/update_s_ssg_acoustic.c
+ 6
− 5
View file @ ede72a71
......@@ -140,7 +140,8 @@ double update_s_acoustic(st_boundary *nb, st_velocity *vel,
for (int j = nb->ny1; j <= nb->ny2; j++) {
for (int i = nb->nx1; i <= nb->nx2; i++) {
#pragma omp simd
float *sp_j_i = stress->sp[j][i];
#pragma omp simd nontemporal(sp_j_i)
for (int k = nb->nz1; k <= nb->nz2; k++) {
/* spatial derivatives of the components of the velocities
* are computed */
......@@ -165,7 +166,7 @@ double update_s_acoustic(st_boundary *nb, st_velocity *vel,
vxxyyzz = vxx + vyy + vzz;
stress->sp[j][i][k] += gv->DT * (g * vxxyyzz);
sp_j_i[k] += gv->DT * (g * vxxyyzz);
}
}
......@@ -246,7 +247,8 @@ double update_s_acoustic(st_boundary *nb, st_velocity *vel,
for (int j = nb->ny1; j <= nb->ny2; j++) {
for (int i = nb->nx1; i <= nb->nx2; i++) {
#pragma omp simd
float *sp_j_i = stress->sp[j][i];
#pragma omp simd nontemporal(sp_j_i)
for (int k = nb->nz1; k <= nb->nz2; k++) {
/* spatial derivatives of the components of the velocities
......@@ -277,8 +279,7 @@ double update_s_acoustic(st_boundary *nb, st_velocity *vel,
vxxyyzz = vxx + vyy + vzz;
stress->sp[j][i][k] += gv->DT * (g * vxxyyzz);
sp_j_i[k] += gv->DT * (g * vxxyyzz);
}
}
}
......
src/update_v_ssg_acoustic.c
+ 60
− 54
View file @ ede72a71
......@@ -125,30 +125,33 @@ void update_v_acoustic(const st_boundary *nb, st_velocity *restrict vel, st_stre
b3 = 0.018063;
b4 = -0.0026274;
} /* Holberg coefficients E=0.1 % */
for (int j = nb->ny1; j <= nb->ny2; j++) {
for (int i = nb->nx1; i <= nb->nx2; i++) {
#pragma omp simd
for (int k = nb->nz1; k <= nb->nz2; k++) {
sp_x = dx * (b1 * (stress->sp[j][i + 1][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i + 2][k] - stress->sp[j][i - 1][k]) +
b3 * (stress->sp[j][i + 3][k] - stress->sp[j][i - 2][k]) +
b4 * (stress->sp[j][i + 4][k] - stress->sp[j][i - 3][k]));
/* updating components of particle velocities */
vel->vx[j][i][k] += (sp_x) / mod_av->rip[j][i][k];
sp_y = dy * (b1 * (stress->sp[j + 1][i][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j + 2][i][k] - stress->sp[j - 1][i][k]) +
b3 * (stress->sp[j + 3][i][k] - stress->sp[j - 2][i][k]) +
b4 * (stress->sp[j + 4][i][k] - stress->sp[j - 3][i][k]));
vel->vy[j][i][k] += (sp_y) / mod_av->rjp[j][i][k];
sp_z = dz * (b1 * (stress->sp[j][i][k + 1] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i][k + 2] - stress->sp[j][i][k - 1]) +
b3 * (stress->sp[j][i][k + 3] - stress->sp[j][i][k - 2]) +
b4 * (stress->sp[j][i][k + 4] - stress->sp[j][i][k - 3]));
vel->vz[j][i][k] += (sp_z) / mod_av->rkp[j][i][k];
}
}
}
break;
for (int j = nb->ny1; j <= nb->ny2; j++) {
for (int i = nb->nx1; i <= nb->nx2; i++) {
float *vx_j_i = vel->vx[j][i];
float *vy_j_i = vel->vy[j][i];
float *vz_j_i = vel->vz[j][i];
#pragma omp simd nontemporal(vx_j_i, vy_j_i, vz_j_i)
for (int k = nb->nz1; k <= nb->nz2; k++) {
sp_x = dx * (b1 * (stress->sp[j][i + 1][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i + 2][k] - stress->sp[j][i - 1][k]) +
b3 * (stress->sp[j][i + 3][k] - stress->sp[j][i - 2][k]) +
b4 * (stress->sp[j][i + 4][k] - stress->sp[j][i - 3][k]));
/* updating components of particle velocities */
vx_j_i[k] += (sp_x) / mod_av->rip[j][i][k];
sp_y = dy * (b1 * (stress->sp[j + 1][i][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j + 2][i][k] - stress->sp[j - 1][i][k]) +
b3 * (stress->sp[j + 3][i][k] - stress->sp[j - 2][i][k]) +
b4 * (stress->sp[j + 4][i][k] - stress->sp[j - 3][i][k]));
vy_j_i[k] += (sp_y) / mod_av->rjp[j][i][k];
sp_z = dz * (b1 * (stress->sp[j][i][k + 1] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i][k + 2] - stress->sp[j][i][k - 1]) +
b3 * (stress->sp[j][i][k + 3] - stress->sp[j][i][k - 2]) +
b4 * (stress->sp[j][i][k + 4] - stress->sp[j][i][k - 3]));
vz_j_i[k] += (sp_z) / mod_av->rkp[j][i][k];
}
}
}
break;
case 10:
b1 = 19845.0 / 16384.0;
b2 = -735.0 / 8192.0;
......@@ -206,36 +209,39 @@ void update_v_acoustic(const st_boundary *nb, st_velocity *restrict vel, st_stre
b5 = 0.0029857;
b6 = -0.00066667;
}
for (int j = nb->ny1; j <= nb->ny2; j++) {
for (int i = nb->nx1; i <= nb->nx2; i++) {
#pragma omp simd
for (int k = nb->nz1; k <= nb->nz2; k++) {
sp_x = dx * (b1 * (stress->sp[j][i + 1][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i + 2][k] - stress->sp[j][i - 1][k]) +
b3 * (stress->sp[j][i + 3][k] - stress->sp[j][i - 2][k]) +
b4 * (stress->sp[j][i + 4][k] - stress->sp[j][i - 3][k]) +
b5 * (stress->sp[j][i + 5][k] - stress->sp[j][i - 4][k]) +
b6 * (stress->sp[j][i + 6][k] - stress->sp[j][i - 5][k]));
/* updating components of particle velocities */
vel->vx[j][i][k] += (sp_x) / mod_av->rip[j][i][k];
sp_y = dy * (b1 * (stress->sp[j + 1][i][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j + 2][i][k] - stress->sp[j - 1][i][k]) +
b3 * (stress->sp[j + 3][i][k] - stress->sp[j - 2][i][k]) +
b4 * (stress->sp[j + 4][i][k] - stress->sp[j - 3][i][k]) +
b5 * (stress->sp[j + 5][i][k] - stress->sp[j - 4][i][k]) +
b6 * (stress->sp[j + 6][i][k] - stress->sp[j - 5][i][k]));
vel->vy[j][i][k] += (sp_y) / mod_av->rjp[j][i][k];
sp_z = dz * (b1 * (stress->sp[j][i][k + 1] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i][k + 2] - stress->sp[j][i][k - 1]) +
b3 * (stress->sp[j][i][k + 3] - stress->sp[j][i][k - 2]) +
b4 * (stress->sp[j][i][k + 4] - stress->sp[j][i][k - 3]) +
b5 * (stress->sp[j][i][k + 5] - stress->sp[j][i][k - 4]) +
b6 * (stress->sp[j][i][k + 6] - stress->sp[j][i][k - 5]));
vel->vz[j][i][k] += (sp_z) / mod_av->rkp[j][i][k];
}
}
}
break;
for (int j = nb->ny1; j <= nb->ny2; j++) {
for (int i = nb->nx1; i <= nb->nx2; i++) {
float *vx_j_i = vel->vx[j][i];
float *vy_j_i = vel->vy[j][i];
float *vz_j_i = vel->vz[j][i];
#pragma omp simd nontemporal(vx_j_i, vy_j_i, vz_j_i)
for (int k = nb->nz1; k <= nb->nz2; k++) {
sp_x = dx * (b1 * (stress->sp[j][i + 1][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i + 2][k] - stress->sp[j][i - 1][k]) +
b3 * (stress->sp[j][i + 3][k] - stress->sp[j][i - 2][k]) +
b4 * (stress->sp[j][i + 4][k] - stress->sp[j][i - 3][k]) +
b5 * (stress->sp[j][i + 5][k] - stress->sp[j][i - 4][k]) +
b6 * (stress->sp[j][i + 6][k] - stress->sp[j][i - 5][k]));
/* updating components of particle velocities */
vx_j_i[k] += (sp_x) / mod_av->rip[j][i][k];
sp_y = dy * (b1 * (stress->sp[j + 1][i][k] - stress->sp[j][i][k]) +
b2 * (stress->sp[j + 2][i][k] - stress->sp[j - 1][i][k]) +
b3 * (stress->sp[j + 3][i][k] - stress->sp[j - 2][i][k]) +
b4 * (stress->sp[j + 4][i][k] - stress->sp[j - 3][i][k]) +
b5 * (stress->sp[j + 5][i][k] - stress->sp[j - 4][i][k]) +
b6 * (stress->sp[j + 6][i][k] - stress->sp[j - 5][i][k]));
vy_j_i[k] += (sp_y) / mod_av->rjp[j][i][k];
sp_z = dz * (b1 * (stress->sp[j][i][k + 1] - stress->sp[j][i][k]) +
b2 * (stress->sp[j][i][k + 2] - stress->sp[j][i][k - 1]) +
b3 * (stress->sp[j][i][k + 3] - stress->sp[j][i][k - 2]) +
b4 * (stress->sp[j][i][k + 4] - stress->sp[j][i][k - 3]) +
b5 * (stress->sp[j][i][k + 5] - stress->sp[j][i][k - 4]) +
b6 * (stress->sp[j][i][k + 6] - stress->sp[j][i][k - 5]));
vz_j_i[k] += (sp_z) / mod_av->rkp[j][i][k];
}
}
}
break;
}
/* absorbing boundary condition (exponential damping) */
......
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

Help | Imprint | Privacy policy | Accessibility | Contact