ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

first step to flip sign of FFTs

Browse Source
This commit is contained in:
Steve Plimpton
2020-12-03 10:28:03 -07:00
parent 5ea9d97024
commit 42e8bfbcaf
3 changed files with 184 additions and 206 deletions
Download Patch File Download Diff File Expand all files Collapse all files

94
src/KSPACE/pppm.cpp
View File

@ -51,6 +51,7 @@ using namespace MathSpecial;
enum{REVERSE_RHO};
enum{FORWARD_IK,FORWARD_AD,FORWARD_IK_PERATOM,FORWARD_AD_PERATOM};
enum{FORWARD=-1,BACKWARD=1};
#ifdef FFT_SINGLE
#define ZEROF 0.0f
@ -2010,7 +2011,7 @@ void PPPM::poisson_ik()
work1[n++] = ZEROF;
}
fft1->compute(work1,work1,1);
fft1->compute(work1,work1,FORWARD);
// global energy and virial contribution
@ -2056,7 +2057,7 @@ void PPPM::poisson_ik()
return;
}
// compute gradients of V(r) in each of 3 dims by transformimg -ik*V(k)
// compute gradients of V(r) in each of 3 dims by transformimg ik*V(k)
// FFT leaves data in 3d brick decomposition
// copy it into inner portion of vdx,vdy,vdz arrays
@ -2066,12 +2067,12 @@ void PPPM::poisson_ik()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work2[n] = fkx[i]*work1[n+1];
work2[n+1] = -fkx[i]*work1[n];
work2[n] = -fkx[i]*work1[n+1];
work2[n+1] = fkx[i]*work1[n];
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2087,12 +2088,12 @@ void PPPM::poisson_ik()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work2[n] = fky[j]*work1[n+1];
work2[n+1] = -fky[j]*work1[n];
work2[n] = -fky[j]*work1[n+1];
work2[n+1] = fky[j]*work1[n];
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2108,12 +2109,12 @@ void PPPM::poisson_ik()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work2[n] = fkz[k]*work1[n+1];
work2[n+1] = -fkz[k]*work1[n];
work2[n] = -fkz[k]*work1[n+1];
work2[n+1] = fkz[k]*work1[n];
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2132,7 +2133,7 @@ void PPPM::poisson_ik_triclinic()
{
int i,j,k,n;
// compute gradients of V(r) in each of 3 dims by transformimg -ik*V(k)
// compute gradients of V(r) in each of 3 dims by transformimg ik*V(k)
// FFT leaves data in 3d brick decomposition
// copy it into inner portion of vdx,vdy,vdz arrays
@ -2140,12 +2141,12 @@ void PPPM::poisson_ik_triclinic()
n = 0;
for (i = 0; i < nfft; i++) {
work2[n] = fkx[i]*work1[n+1];
work2[n+1] = -fkx[i]*work1[n];
work2[n] = -fkx[i]*work1[n+1];
work2[n+1] = fkx[i]*work1[n];
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2159,12 +2160,12 @@ void PPPM::poisson_ik_triclinic()
n = 0;
for (i = 0; i < nfft; i++) {
work2[n] = fky[i]*work1[n+1];
work2[n+1] = -fky[i]*work1[n];
work2[n] = -fky[i]*work1[n+1];
work2[n+1] = fky[i]*work1[n];
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2178,12 +2179,12 @@ void PPPM::poisson_ik_triclinic()
n = 0;
for (i = 0; i < nfft; i++) {
work2[n] = fkz[i]*work1[n+1];
work2[n+1] = -fkz[i]*work1[n];
work2[n] = -fkz[i]*work1[n+1];
work2[n+1] = fkz[i]*work1[n];
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2211,7 +2212,7 @@ void PPPM::poisson_ad()
work1[n++] = ZEROF;
}
fft1->compute(work1,work1,1);
fft1->compute(work1,work1,FORWARD);
// global energy and virial contribution
@ -2257,7 +2258,7 @@ void PPPM::poisson_ad()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2286,7 +2287,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2308,7 +2309,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2325,7 +2326,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2342,7 +2343,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2359,7 +2360,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2376,7 +2377,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2393,7 +2394,7 @@ void PPPM::poisson_peratom()
n += 2;
}
fft2->compute(work2,work2,-1);
fft2->compute(work2,work2,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -3002,11 +3003,11 @@ int PPPM::timing_1d(int n, double &time1d)
time1 = MPI_Wtime();
for (int i = 0; i < n; i++) {
fft1->timing1d(work1,nfft_both,1);
fft2->timing1d(work1,nfft_both,-1);
fft1->timing1d(work1,nfft_both,FORWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
if (differentiation_flag != 1) {
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
}
}
@ -3032,11 +3033,11 @@ int PPPM::timing_3d(int n, double &time3d)
time1 = MPI_Wtime();
for (int i = 0; i < n; i++) {
fft1->compute(work1,work1,1);
fft2->compute(work1,work1,-1);
fft1->compute(work1,work1,FORWARD);
fft2->compute(work1,work1,BACKWARD);
if (differentiation_flag != 1) {
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
}
}
@ -3313,7 +3314,7 @@ void PPPM::poisson_groups(int AA_flag)
work_A[n++] = ZEROF;
}
fft1->compute(work_A,work_A,1);
fft1->compute(work_A,work_A,FORWARD);
// group B
@ -3323,7 +3324,7 @@ void PPPM::poisson_groups(int AA_flag)
work_B[n++] = ZEROF;
}
fft1->compute(work_B,work_B,1);
fft1->compute(work_B,work_B,FORWARD);
// group-group energy and force contribution,
// keep everything in reciprocal space so
@ -3343,7 +3344,6 @@ void PPPM::poisson_groups(int AA_flag)
if (AA_flag) return;
// multiply by Green's function and s2
// (only for work_A so it is not squared below)
@ -3368,7 +3368,7 @@ void PPPM::poisson_groups(int AA_flag)
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
partial_group = work_A[n+1]*work_B[n] - work_A[n]*work_B[n+1];
partial_group = work_A[n]*work_B[n+1] - work_A[n+1]*work_B[n];
f2group[0] += fkx[i] * partial_group;
n += 2;
}
@ -3379,7 +3379,7 @@ void PPPM::poisson_groups(int AA_flag)
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
partial_group = work_A[n+1]*work_B[n] - work_A[n]*work_B[n+1];
partial_group = work_A[n]*work_B[n+1] - work_A[n+1]*work_B[n];
f2group[1] += fky[j] * partial_group;
n += 2;
}
@ -3390,7 +3390,7 @@ void PPPM::poisson_groups(int AA_flag)
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
partial_group = work_A[n+1]*work_B[n] - work_A[n]*work_B[n+1];
partial_group = work_A[n]*work_B[n+1] - work_A[n+1]*work_B[n];
f2group[2] += fkz[k] * partial_group;
n += 2;
}
@ -3416,7 +3416,7 @@ void PPPM::poisson_groups_triclinic()
n = 0;
for (i = 0; i < nfft; i++) {
partial_group = work_A[n+1]*work_B[n] - work_A[n]*work_B[n+1];
partial_group = work_A[n]*work_B[n+1] - work_A[n+1]*work_B[n];
f2group[0] += fkx[i] * partial_group;
n += 2;
}
@ -3425,7 +3425,7 @@ void PPPM::poisson_groups_triclinic()
n = 0;
for (i = 0; i < nfft; i++) {
partial_group = work_A[n+1]*work_B[n] - work_A[n]*work_B[n+1];
partial_group = work_A[n]*work_B[n+1] - work_A[n+1]*work_B[n];
f2group[1] += fky[i] * partial_group;
n += 2;
}
@ -3434,7 +3434,7 @@ void PPPM::poisson_groups_triclinic()
n = 0;
for (i = 0; i < nfft; i++) {
partial_group = work_A[n+1]*work_B[n] - work_A[n]*work_B[n+1];
partial_group = work_A[n]*work_B[n+1] - work_A[n+1]*work_B[n];
f2group[2] += fkz[i] * partial_group;
n += 2;
}

159
src/KSPACE/pppm_dipole.cpp
View File

@ -46,6 +46,7 @@ using namespace MathSpecial;
enum{REVERSE_MU};
enum{FORWARD_MU,FORWARD_MU_PERATOM};
enum{FORWARD=-1,BACKWARD=1};
#ifdef FFT_SINGLE
#define ZEROF 0.0f
@ -462,8 +463,8 @@ void PPPMDipole::compute(int eflag, int vflag)
// extra per-atom energy/virial communication
if (evflag_atom)
gc->forward_comm_kspace(this,18,sizeof(FFT_SCALAR),FORWARD_MU_PERATOM,
gc_buf1,gc_buf2,MPI_FFT_SCALAR);
gc_dipole->forward_comm_kspace(this,18,sizeof(FFT_SCALAR),FORWARD_MU_PERATOM,
gc_buf1,gc_buf2,MPI_FFT_SCALAR);
// calculate the force on my particles
@ -610,7 +611,7 @@ void PPPMDipole::allocate()
nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in,
nxlo_out,nxhi_out,nylo_out,nyhi_out,nzlo_out,nzhi_out);
gc->setup(ngc_buf1,ngc_buf2);
gc_dipole->setup(ngc_buf1,ngc_buf2);
npergrid = 9;
@ -643,29 +644,10 @@ void PPPMDipole::deallocate()
memory->destroy(densityy_fft_dipole);
memory->destroy(densityz_fft_dipole);
memory->destroy(greensfn);
memory->destroy(work1);
memory->destroy(work2);
memory->destroy(work3);
memory->destroy(work4);
memory->destroy(vg);
memory->destroy1d_offset(fkx,nxlo_fft);
memory->destroy1d_offset(fky,nylo_fft);
memory->destroy1d_offset(fkz,nzlo_fft);
memory->destroy(gf_b);
memory->destroy2d_offset(rho1d,-order_allocated/2);
memory->destroy2d_offset(drho1d,-order_allocated/2);
memory->destroy2d_offset(rho_coeff,(1-order_allocated)/2);
memory->destroy2d_offset(drho_coeff,(1-order_allocated)/2);
delete fft1;
delete fft2;
delete remap;
delete gc_dipole;
memory->destroy(gc_buf1);
memory->destroy(gc_buf2);
}
/* ----------------------------------------------------------------------
@ -1312,9 +1294,9 @@ void PPPMDipole::poisson_ik_dipole()
n += 2;
}
fft1->compute(work1,work1,1);
fft1->compute(work2,work2,1);
fft1->compute(work3,work3,1);
fft1->compute(work1,work1,FORWARD);
fft1->compute(work2,work2,FORWARD);
fft1->compute(work3,work3,FORWARD);
// global energy and virial contribution
@ -1390,7 +1372,7 @@ void PPPMDipole::poisson_ik_dipole()
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1411,7 +1393,7 @@ void PPPMDipole::poisson_ik_dipole()
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1432,7 +1414,7 @@ void PPPMDipole::poisson_ik_dipole()
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1448,12 +1430,12 @@ void PPPMDipole::poisson_ik_dipole()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work4[n] = fkx[i]*fkx[i]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = -fkx[i]*fkx[i]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
work4[n] = -fkx[i]*fkx[i]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = +fkx[i]*fkx[i]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1469,12 +1451,12 @@ void PPPMDipole::poisson_ik_dipole()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work4[n] = fky[j]*fky[j]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = -fky[j]*fky[j]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
work4[n] = -fky[j]*fky[j]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = fky[j]*fky[j]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1490,12 +1472,12 @@ void PPPMDipole::poisson_ik_dipole()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work4[n] = fkz[k]*fkz[k]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = -fkz[k]*fkz[k]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
work4[n] = -fkz[k]*fkz[k]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = fkz[k]*fkz[k]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1511,12 +1493,12 @@ void PPPMDipole::poisson_ik_dipole()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work4[n] = fkx[i]*fky[j]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = -fkx[i]*fky[j]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
work4[n] = -fkx[i]*fky[j]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = fkx[i]*fky[j]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1532,12 +1514,12 @@ void PPPMDipole::poisson_ik_dipole()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work4[n] = fkx[i]*fkz[k]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = -fkx[i]*fkz[k]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
work4[n] = -fkx[i]*fkz[k]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = fkx[i]*fkz[k]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1553,12 +1535,12 @@ void PPPMDipole::poisson_ik_dipole()
for (k = nzlo_fft; k <= nzhi_fft; k++)
for (j = nylo_fft; j <= nyhi_fft; j++)
for (i = nxlo_fft; i <= nxhi_fft; i++) {
work4[n] = fky[j]*fkz[k]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = -fky[j]*fkz[k]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
work4[n] = -fky[j]*fkz[k]*(work1[n+1]*fkx[i] + work2[n+1]*fky[j] + work3[n+1]*fkz[k]);
work4[n+1] = fky[j]*fkz[k]*(work1[n]*fkx[i] + work2[n]*fky[j] + work3[n]*fkz[k]);
n += 2;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1594,7 +1576,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1617,7 +1599,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1640,7 +1622,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1663,7 +1645,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1686,7 +1668,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1709,7 +1691,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1732,7 +1714,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1755,7 +1737,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1778,7 +1760,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1801,7 +1783,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1824,7 +1806,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1847,7 +1829,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1870,7 +1852,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1893,7 +1875,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1916,7 +1898,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1939,7 +1921,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1962,7 +1944,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -1985,7 +1967,7 @@ void PPPMDipole::poisson_peratom_dipole()
ii++;
}
fft2->compute(work4,work4,-1);
fft2->compute(work4,work4,BACKWARD);
n = 0;
for (k = nzlo_in; k <= nzhi_in; k++)
@ -2013,7 +1995,6 @@ void PPPMDipole::fieldforce_ik_dipole()
// (dx,dy,dz) = distance to "lower left" grid pt
// (mx,my,mz) = global coords of moving stencil pt
double **mu = atom->mu;
double **x = atom->x;
double **f = atom->f;
@ -2417,18 +2398,18 @@ int PPPMDipole::timing_1d(int n, double &time1d)
time1 = MPI_Wtime();
for (int i = 0; i < n; i++) {
fft1->timing1d(work1,nfft_both,1);
fft1->timing1d(work1,nfft_both,1);
fft1->timing1d(work1,nfft_both,1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft1->timing1d(work1,nfft_both,FORWARD);
fft1->timing1d(work1,nfft_both,FORWARD);
fft1->timing1d(work1,nfft_both,FORWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
}
MPI_Barrier(world);
@ -2452,18 +2433,18 @@ int PPPMDipole::timing_3d(int n, double &time3d)
time1 = MPI_Wtime();
for (int i = 0; i < n; i++) {
fft1->compute(work1,work1,1);
fft1->compute(work1,work1,1);
fft1->compute(work1,work1,1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft1->compute(work1,work1,FORWARD);
fft1->compute(work1,work1,FORWARD);
fft1->compute(work1,work1,FORWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
}
MPI_Barrier(world);

137
src/KSPACE/pppm_disp.cpp
View File

@ -51,6 +51,7 @@ enum{FORWARD_IK, FORWARD_AD, FORWARD_IK_PERATOM, FORWARD_AD_PERATOM,
FORWARD_IK_A, FORWARD_AD_A, FORWARD_IK_PERATOM_A, FORWARD_AD_PERATOM_A,
FORWARD_IK_NONE, FORWARD_AD_NONE, FORWARD_IK_PERATOM_NONE,
FORWARD_AD_PERATOM_NONE};
enum{FORWARD=-1,BACKWARD=1};
#ifdef FFT_SINGLE
#define ZEROF 0.0f
@ -4555,7 +4556,7 @@ void PPPMDisp::poisson_ik(FFT_SCALAR* wk1, FFT_SCALAR* wk2,
wk1[n++] = ZEROF;
}
ft1->compute(wk1,wk1,1);
ft1->compute(wk1,wk1,FORWARD);
// if requested, compute energy and virial contribution
@ -4600,12 +4601,12 @@ void PPPMDisp::poisson_ik(FFT_SCALAR* wk1, FFT_SCALAR* wk2,
for (k = nzlo_ft; k <= nzhi_ft; k++)
for (j = nylo_ft; j <= nyhi_ft; j++)
for (i = nxlo_ft; i <= nxhi_ft; i++) {
wk2[n] = 0.5*(kx[i]-kx2[i])*wk1[n+1] + 0.5*(ky[j]-ky2[j])*wk1[n];
wk2[n+1] = -0.5*(kx[i]-kx2[i])*wk1[n] + 0.5*(ky[j]-ky2[j])*wk1[n+1];
wk2[n] = -0.5*(kx[i]-kx2[i])*wk1[n+1] + 0.5*(ky[j]-ky2[j])*wk1[n];
wk2[n+1] = 0.5*(kx[i]-kx2[i])*wk1[n] + 0.5*(ky[j]-ky2[j])*wk1[n+1];
n += 2;
}
ft2->compute(wk2,wk2,-1);
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4622,13 +4623,12 @@ void PPPMDisp::poisson_ik(FFT_SCALAR* wk1, FFT_SCALAR* wk2,
for (k = nzlo_ft; k <= nzhi_ft; k++)
for (j = nylo_ft; j <= nyhi_ft; j++)
for (i = nxlo_ft; i <= nxhi_ft; i++) {
wk2[n] = kz[k]*wk1[n+1];
wk2[n+1] = -kz[k]*wk1[n];
wk2[n] = -kz[k]*wk1[n+1];
wk2[n+1] = kz[k]*wk1[n];
n += 2;
}
ft2->compute(wk2,wk2,-1);
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4647,12 +4647,12 @@ void PPPMDisp::poisson_ik(FFT_SCALAR* wk1, FFT_SCALAR* wk2,
for (k = nzlo_ft; k <= nzhi_ft; k++)
for (j = nylo_ft; j <= nyhi_ft; j++)
for (i = nxlo_ft; i <= nxhi_ft; i++) {
wk2[n] = 0.5*(kz[k]-kz2[k])*wk1[n+1] - wk1[n+1];
wk2[n+1] = -0.5*(kz[k]-kz2[k])*wk1[n] + wk1[n];
wk2[n] = -0.5*(kz[k]-kz2[k])*wk1[n+1] - wk1[n+1];
wk2[n+1] = 0.5*(kz[k]-kz2[k])*wk1[n] + wk1[n];
n += 2;
}
ft2->compute(wk2,wk2,-1);
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4696,7 +4696,7 @@ void PPPMDisp::poisson_ad(FFT_SCALAR* wk1, FFT_SCALAR* wk2,
wk1[n++] = ZEROF;
}
ft1->compute(wk1,wk1,1);
ft1->compute(wk1,wk1,FORWARD);
// if requested, compute energy and virial contribution
@ -4739,10 +4739,9 @@ void PPPMDisp::poisson_ad(FFT_SCALAR* wk1, FFT_SCALAR* wk2,
wk2[n] = wk1[n];
wk2[n+1] = wk1[n+1];
n += 2;
}
ft2->compute(wk2,wk2,-1);
}
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4781,7 +4780,7 @@ void PPPMDisp:: poisson_peratom(FFT_SCALAR* wk1, FFT_SCALAR* wk2, LAMMPS_NS::FFT
n += 2;
}
ft2->compute(wk2,wk2,-1);
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4800,7 +4799,7 @@ void PPPMDisp:: poisson_peratom(FFT_SCALAR* wk1, FFT_SCALAR* wk2, LAMMPS_NS::FFT
n += 2;
}
ft2->compute(wk2,wk2,-1);
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4819,7 +4818,7 @@ void PPPMDisp:: poisson_peratom(FFT_SCALAR* wk1, FFT_SCALAR* wk2, LAMMPS_NS::FFT
n += 2;
}
ft2->compute(wk2,wk2,-1);
ft2->compute(wk2,wk2,BACKWARD);
n = 0;
for (k = nzlo_i; k <= nzhi_i; k++)
@ -4863,7 +4862,7 @@ poisson_2s_ik(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work1_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
}
// two transforms are required when energies and pressures are
// calculated
@ -4876,8 +4875,8 @@ poisson_2s_ik(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work2_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work2_6,work2_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
fft1_6->compute(work2_6,work2_6,FORWARD);
double s2 = scaleinv*scaleinv;
@ -4927,7 +4926,7 @@ poisson_2s_ik(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -4948,7 +4947,7 @@ poisson_2s_ik(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -4969,7 +4968,7 @@ poisson_2s_ik(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -4989,7 +4988,7 @@ poisson_2s_ik(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5037,7 +5036,7 @@ poisson_none_ik(int n1, int n2,FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work1_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
}
// two transforms are required when energies and pressures are calculated
@ -5052,8 +5051,8 @@ poisson_none_ik(int n1, int n2,FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work2_6,work2_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
fft1_6->compute(work2_6,work2_6,FORWARD);
double s2 = scaleinv*scaleinv;
@ -5106,7 +5105,7 @@ poisson_none_ik(int n1, int n2,FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5127,7 +5126,7 @@ poisson_none_ik(int n1, int n2,FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5148,7 +5147,7 @@ poisson_none_ik(int n1, int n2,FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5168,7 +5167,7 @@ poisson_none_ik(int n1, int n2,FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5213,7 +5212,7 @@ poisson_2s_ad(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work1_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
}
// two transforms are required when energies and pressures are calculated
@ -5227,8 +5226,8 @@ poisson_2s_ad(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work2_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work2_6,work2_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
fft1_6->compute(work2_6,work2_6,FORWARD);
double s2 = scaleinv*scaleinv;
@ -5270,7 +5269,7 @@ poisson_2s_ad(FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5312,7 +5311,7 @@ poisson_none_ad(int n1, int n2, FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work1_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
}
// two transforms are required when energies and pressures are calculated
@ -5326,8 +5325,8 @@ poisson_none_ad(int n1, int n2, FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
work2_6[n++] = dfft_2[i];
}
fft1_6->compute(work1_6,work1_6,1);
fft1_6->compute(work2_6,work2_6,1);
fft1_6->compute(work1_6,work1_6,FORWARD);
fft1_6->compute(work2_6,work2_6,FORWARD);
double s2 = scaleinv*scaleinv;
@ -5374,7 +5373,7 @@ poisson_none_ad(int n1, int n2, FFT_SCALAR* dfft_1, FFT_SCALAR* dfft_2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5410,7 +5409,7 @@ poisson_2s_peratom(FFT_SCALAR*** v0_pa_1, FFT_SCALAR*** v1_pa_1, FFT_SCALAR*** v
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5429,7 +5428,7 @@ poisson_2s_peratom(FFT_SCALAR*** v0_pa_1, FFT_SCALAR*** v1_pa_1, FFT_SCALAR*** v
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5448,7 +5447,7 @@ poisson_2s_peratom(FFT_SCALAR*** v0_pa_1, FFT_SCALAR*** v1_pa_1, FFT_SCALAR*** v
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5467,7 +5466,7 @@ poisson_2s_peratom(FFT_SCALAR*** v0_pa_1, FFT_SCALAR*** v1_pa_1, FFT_SCALAR*** v
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5486,7 +5485,7 @@ poisson_2s_peratom(FFT_SCALAR*** v0_pa_1, FFT_SCALAR*** v1_pa_1, FFT_SCALAR*** v
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5505,7 +5504,7 @@ poisson_2s_peratom(FFT_SCALAR*** v0_pa_1, FFT_SCALAR*** v1_pa_1, FFT_SCALAR*** v
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5541,7 +5540,7 @@ poisson_none_peratom(int n1, int n2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5560,7 +5559,7 @@ poisson_none_peratom(int n1, int n2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5579,7 +5578,7 @@ poisson_none_peratom(int n1, int n2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5598,7 +5597,7 @@ poisson_none_peratom(int n1, int n2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5617,7 +5616,7 @@ poisson_none_peratom(int n1, int n2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -5636,7 +5635,7 @@ poisson_none_peratom(int n1, int n2,
n += 2;
}
fft2_6->compute(work2_6,work2_6,-1);
fft2_6->compute(work2_6,work2_6,BACKWARD);
n = 0;
for (k = nzlo_in_6; k <= nzhi_in_6; k++)
@ -8153,11 +8152,11 @@ int PPPMDisp::timing_1d(int n, double &time1d)
if (function[0]) {
for (int i = 0; i < n; i++) {
fft1->timing1d(work1,nfft_both,1);
fft2->timing1d(work1,nfft_both,-1);
fft1->timing1d(work1,nfft_both,FORWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
if (differentiation_flag != 1){
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,-1);
fft2->timing1d(work1,nfft_both,BACKWARD);
fft2->timing1d(work1,nfft_both,BACKWARD);
}
}
}
@ -8171,11 +8170,11 @@ int PPPMDisp::timing_1d(int n, double &time1d)
if (function[1] + function[2] + function[3]) {
for (int i = 0; i < n; i++) {
fft1_6->timing1d(work1_6,nfft_both_6,1);
fft2_6->timing1d(work1_6,nfft_both_6,-1);
fft1_6->timing1d(work1_6,nfft_both_6,FORWARD);
fft2_6->timing1d(work1_6,nfft_both_6,BACKWARD);
if (differentiation_flag != 1){
fft2_6->timing1d(work1_6,nfft_both_6,-1);
fft2_6->timing1d(work1_6,nfft_both_6,-1);
fft2_6->timing1d(work1_6,nfft_both_6,BACKWARD);
fft2_6->timing1d(work1_6,nfft_both_6,BACKWARD);
}
}
}
@ -8203,18 +8202,16 @@ int PPPMDisp::timing_3d(int n, double &time3d)
if (function[1] + function[2] + function[3])
for (int i = 0; i < 2*nfft_both_6; i++) work1_6[i] = ZEROF;
MPI_Barrier(world);
time1 = MPI_Wtime();
if (function[0]) {
for (int i = 0; i < n; i++) {
fft1->compute(work1,work1,1);
fft2->compute(work1,work1,-1);
fft1->compute(work1,work1,FORWARD);
fft2->compute(work1,work1,BACKWARD);
if (differentiation_flag != 1) {
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,-1);
fft2->compute(work1,work1,BACKWARD);
fft2->compute(work1,work1,BACKWARD);
}
}
}
@ -8228,11 +8225,11 @@ int PPPMDisp::timing_3d(int n, double &time3d)
if (function[1] + function[2] + function[3]) {
for (int i = 0; i < n; i++) {
fft1_6->compute(work1_6,work1_6,1);
fft2_6->compute(work1_6,work1_6,-1);
fft1_6->compute(work1_6,work1_6,FORWARD);
fft2_6->compute(work1_6,work1_6,BACKWARD);
if (differentiation_flag != 1) {
fft2_6->compute(work1_6,work1_6,-1);
fft2_6->compute(work1_6,work1_6,-1);
fft2_6->compute(work1_6,work1_6,BACKWARD);
fft2_6->compute(work1_6,work1_6,BACKWARD);
}
}
}