KSpace bug fixes

This commit is contained in:
Steve Plimpton
2022-11-08 11:22:25 -07:00
parent 90b54300e9
commit 65ce9aa791
3 changed files with 95 additions and 67 deletions

View File

@ -868,7 +868,7 @@ void PPPM::deallocate()
memory->destroy(gc_buf2);
memory->destroy3d_offset(density_brick,nzlo_out,nylo_out,nxlo_out);
if (differentiation_flag == 1) {
memory->destroy3d_offset(u_brick,nzlo_out,nylo_out,nxlo_out);
memory->destroy(sf_precoeff1);
@ -882,7 +882,7 @@ void PPPM::deallocate()
memory->destroy3d_offset(vdy_brick,nzlo_out,nylo_out,nxlo_out);
memory->destroy3d_offset(vdz_brick,nzlo_out,nylo_out,nxlo_out);
}
memory->destroy(density_fft);
memory->destroy(greensfn);
memory->destroy(work1);

View File

@ -205,7 +205,7 @@ void PPPMDipole::init()
gc_dipole->set_stencil_atom(-nlower,nupper);
gc_dipole->set_shift_atom(shiftatom_lo,shiftatom_hi);
gc_dipole->set_zfactor(slab_volfactor);
gc_dipole->setup_grid(nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in,
nxlo_out,nxhi_out,nylo_out,nyhi_out,nzlo_out,nzhi_out);
@ -231,6 +231,17 @@ void PPPMDipole::init()
double estimated_accuracy = final_accuracy_dipole();
// allocate K-space dependent memory
// don't invoke allocate peratom(), will be allocated when needed
allocate();
// pre-compute Green's function denomiator expansion
// pre-compute 1d charge distribution coefficients
compute_gf_denom();
compute_rho_coeff();
// print stats
int ngrid_max,nfft_both_max;
@ -250,17 +261,6 @@ void PPPMDipole::init()
ngrid_max,nfft_both_max);
utils::logmesg(lmp,mesg);
}
// allocate K-space dependent memory
// don't invoke allocate peratom(), will be allocated when needed
allocate();
// pre-compute Green's function denomiator expansion
// pre-compute 1d charge distribution coefficients
compute_gf_denom();
compute_rho_coeff();
}
/* ----------------------------------------------------------------------
@ -541,7 +541,7 @@ void PPPMDipole::allocate()
gc_dipole->set_stencil_atom(-nlower,nupper);
gc_dipole->set_shift_atom(shiftatom_lo,shiftatom_hi);
gc_dipole->set_zfactor(slab_volfactor);
gc_dipole->setup_grid(nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in,
nxlo_out,nxhi_out,nylo_out,nyhi_out,nzlo_out,nzhi_out);
@ -552,6 +552,24 @@ void PPPMDipole::allocate()
memory->create(gc_buf1,npergrid*ngc_buf1,"pppm:gc_buf1");
memory->create(gc_buf2,npergrid*ngc_buf2,"pppm:gc_buf2");
// tally local grid sizes
// ngrid = count of owned+ghost grid cells on this proc
// nfft_brick = FFT points in 3d brick-decomposition on this proc
// same as count of owned grid cells
// nfft = FFT points in x-pencil FFT decomposition on this proc
// nfft_both = greater of nfft and nfft_brick
ngrid = (nxhi_out-nxlo_out+1) * (nyhi_out-nylo_out+1) *
(nzhi_out-nzlo_out+1);
nfft_brick = (nxhi_in-nxlo_in+1) * (nyhi_in-nylo_in+1) *
(nzhi_in-nzlo_in+1);
nfft = (nxhi_fft-nxlo_fft+1) * (nyhi_fft-nylo_fft+1) *
(nzhi_fft-nzlo_fft+1);
nfft_both = MAX(nfft,nfft_brick);
// allocate distributed grid data
memory->create3d_offset(densityx_brick_dipole,nzlo_out,nzhi_out,nylo_out,nyhi_out,
@ -638,6 +656,7 @@ void PPPMDipole::deallocate()
delete gc_dipole;
memory->destroy(gc_buf1);
memory->destroy(gc_buf2);
gc_buf1 = gc_buf2 = nullptr;
memory->destroy3d_offset(densityx_brick_dipole,nzlo_out,nylo_out,nxlo_out);
memory->destroy3d_offset(densityy_brick_dipole,nzlo_out,nylo_out,nxlo_out);
@ -664,6 +683,9 @@ void PPPMDipole::deallocate()
delete fft1;
delete fft2;
delete remap;
fft1 = fft2 = nullptr;
remap = nullptr;
}
/* ----------------------------------------------------------------------

View File

@ -389,8 +389,8 @@ void PPPMDisp::init()
alpha = qdist / (cos(0.5*theta) * blen);
}
//if g_ewald and g_ewald_6 have not been specified, set some initial value
// to avoid problems when calculating the energies!
// if g_ewald and g_ewald_6 have not been specified,
// set some initial value, to avoid problems when calculating the energies!
if (!gewaldflag) g_ewald = 1;
if (!gewaldflag_6) g_ewald_6 = 1;
@ -407,6 +407,9 @@ void PPPMDisp::init()
if (accuracy_absolute >= 0.0) accuracy = accuracy_absolute;
else accuracy = accuracy_relative * two_charge_force;
double acc;
double acc_6,acc_real_6,acc_kspace_6;
int iteration = 0;
if (function[0]) {
@ -459,31 +462,9 @@ void PPPMDisp::init()
if (!gewaldflag) adjust_gewald();
// calculate the final accuracy
// calculate the final Coulomb accuracy
double acc = final_accuracy();
// print stats
int ngrid_max,nfft_both_max;
MPI_Allreduce(&ngrid,&ngrid_max,1,MPI_INT,MPI_MAX,world);
MPI_Allreduce(&nfft_both,&nfft_both_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
std::string mesg = fmt::format(" Coulomb G vector (1/distance)= {:.16g}\n",
g_ewald);
mesg += fmt::format(" Coulomb grid = {} {} {}\n",
nx_pppm,ny_pppm,nz_pppm);
mesg += fmt::format(" Coulomb stencil order = {}\n",order);
mesg += fmt::format(" Coulomb estimated absolute RMS force accuracy "
"= {:.8g}\n",acc);
mesg += fmt::format(" Coulomb estimated relative force accuracy = {:.8g}\n",
acc/two_charge_force);
mesg += " using " LMP_FFT_PREC " precision " LMP_FFT_LIB "\n";
mesg += fmt::format(" 3d grid and FFT values/proc = {} {}\n",
ngrid_max,nfft_both_max);
utils::logmesg(lmp,mesg);
}
acc = final_accuracy();
}
iteration = 0;
@ -537,34 +518,11 @@ void PPPMDisp::init()
if (!gewaldflag_6 && accuracy_kspace_6 == accuracy_real_6) adjust_gewald_6();
// calculate the final accuracy
// calculate the final displerson accuracy
double acc,acc_real,acc_kspace;
final_accuracy_6(acc,acc_real,acc_kspace);
// print stats
int ngrid_6_max,nfft_both_6_max;
MPI_Allreduce(&ngrid_6,&ngrid_6_max,1,MPI_INT,MPI_MAX,world);
MPI_Allreduce(&nfft_both_6,&nfft_both_6_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
std::string mesg = fmt::format(" Dispersion G vector (1/distance)= "
"{:.16}\n",g_ewald_6);
mesg += fmt::format(" Dispersion grid = {} {} {}\n",
nx_pppm_6,ny_pppm_6,nz_pppm_6);
mesg += fmt::format(" Dispersion stencil order = {}\n",order_6);
mesg += fmt::format(" Dispersion estimated absolute RMS force accuracy "
"= {:.8}\n",acc);
mesg += fmt::format(" Dispersion estimated relative force accuracy "
"= {:.8}\n",acc/two_charge_force);
mesg += " using " LMP_FFT_PREC " precision " LMP_FFT_LIB "\n";
mesg += fmt::format(" 3d grid and FFT values/proc = {} {}\n",
ngrid_6_max,nfft_both_6_max);
utils::logmesg(lmp,mesg);
}
final_accuracy_6(acc_6,acc_real_6,acc_kspace_6);
}
// allocate K-space dependent memory
allocate();
@ -592,6 +550,54 @@ void PPPMDisp::init()
sf_precoeff1_6,sf_precoeff2_6,sf_precoeff3_6,
sf_precoeff4_6,sf_precoeff5_6,sf_precoeff6_6);
}
// print Coulomb stats
if (function[0]) {
int ngrid_max,nfft_both_max;
MPI_Allreduce(&ngrid,&ngrid_max,1,MPI_INT,MPI_MAX,world);
MPI_Allreduce(&nfft_both,&nfft_both_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
std::string mesg = fmt::format(" Coulomb G vector (1/distance)= {:.16g}\n",
g_ewald);
mesg += fmt::format(" Coulomb grid = {} {} {}\n",
nx_pppm,ny_pppm,nz_pppm);
mesg += fmt::format(" Coulomb stencil order = {}\n",order);
mesg += fmt::format(" Coulomb estimated absolute RMS force accuracy "
"= {:.8g}\n",acc);
mesg += fmt::format(" Coulomb estimated relative force accuracy = {:.8g}\n",
acc/two_charge_force);
mesg += " using " LMP_FFT_PREC " precision " LMP_FFT_LIB "\n";
mesg += fmt::format(" 3d grid and FFT values/proc = {} {}\n",
ngrid_max,nfft_both_max);
utils::logmesg(lmp,mesg);
}
}
// print dipserion stats
if (function[1] + function[2] + function[3]) {
int ngrid_6_max,nfft_both_6_max;
MPI_Allreduce(&ngrid_6,&ngrid_6_max,1,MPI_INT,MPI_MAX,world);
MPI_Allreduce(&nfft_both_6,&nfft_both_6_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
std::string mesg = fmt::format(" Dispersion G vector (1/distance)= "
"{:.16}\n",g_ewald_6);
mesg += fmt::format(" Dispersion grid = {} {} {}\n",
nx_pppm_6,ny_pppm_6,nz_pppm_6);
mesg += fmt::format(" Dispersion stencil order = {}\n",order_6);
mesg += fmt::format(" Dispersion estimated absolute RMS force accuracy "
"= {:.8}\n",acc_6);
mesg += fmt::format(" Dispersion estimated relative force accuracy "
"= {:.8}\n",acc_6/two_charge_force);
mesg += " using " LMP_FFT_PREC " precision " LMP_FFT_LIB "\n";
mesg += fmt::format(" 3d grid and FFT values/proc = {} {}\n",
ngrid_6_max,nfft_both_6_max);
utils::logmesg(lmp,mesg);
}
}
}
/* ----------------------------------------------------------------------