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

simplify output in KSPACE package

Browse Source
This commit is contained in:
Axel Kohlmeyer
2020-06-29 08:17:27 -04:00
parent b1d4ca2b03
commit 7de4655e69
11 changed files with 183 additions and 338 deletions
Download Patch File Download Diff File Expand all files Collapse all files

40
src/KSPACE/ewald.cpp
View File

@ -21,6 +21,7 @@
#include "ewald.h"
#include <mpi.h>
#include <cmath>
#include <string>
#include "atom.h"
#include "comm.h"
#include "force.h"
@ -29,6 +30,8 @@
#include "math_const.h"
#include "memory.h"
#include "error.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace MathConst;
@ -88,10 +91,7 @@ Ewald::~Ewald()
void Ewald::init()
{
if (comm->me == 0) {
if (screen) fprintf(screen,"Ewald initialization ...\n");
if (logfile) fprintf(logfile,"Ewald initialization ...\n");
}
if (comm->me == 0) utils::logmesg(lmp,"Ewald initialization ...\n");
// error check
@ -185,28 +185,16 @@ void Ewald::init()
// stats
if (comm->me == 0) {
if (screen) {
fprintf(screen," G vector (1/distance) = %g\n",g_ewald);
fprintf(screen," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(screen," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(screen," KSpace vectors: actual max1d max3d = %d %d %d\n",
kcount,kmax,kmax3d);
fprintf(screen," kxmax kymax kzmax = %d %d %d\n",
kxmax,kymax,kzmax);
}
if (logfile) {
fprintf(logfile," G vector (1/distance) = %g\n",g_ewald);
fprintf(logfile," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(logfile," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(logfile," KSpace vectors: actual max1d max3d = %d %d %d\n",
kcount,kmax,kmax3d);
fprintf(logfile," kxmax kymax kzmax = %d %d %d\n",
kxmax,kymax,kzmax);
}
std::string mesg = fmt::format(" G vector (1/distance) = {}\n",g_ewald);
mesg += fmt::format(" estimated absolute RMS force accuracy = {}\n",
estimated_accuracy);
mesg += fmt::format(" estimated relative force accuracy = {}\n",
estimated_accuracy/two_charge_force);
mesg += fmt::format(" KSpace vectors: actual max1d max3d = {} {} {}\n",
kcount,kmax,kmax3d);
mesg += fmt::format(" kxmax kymax kzmax = {} {} {}\n",
kxmax,kymax,kzmax);
utils::logmesg(lmp,mesg);
}
}

40
src/KSPACE/ewald_dipole.cpp
View File

@ -18,6 +18,7 @@
#include "ewald_dipole.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include <cmath>
#include "atom.h"
#include "comm.h"
@ -29,6 +30,8 @@
#include "memory.h"
#include "error.h"
#include "update.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace MathConst;
@ -63,10 +66,7 @@ EwaldDipole::~EwaldDipole()
void EwaldDipole::init()
{
if (comm->me == 0) {
if (screen) fprintf(screen,"EwaldDipole initialization ...\n");
if (logfile) fprintf(logfile,"EwaldDipole initialization ...\n");
}
if (comm->me == 0) utils::logmesg(lmp,"EwaldDipole initialization ...\n");
// error check
@ -192,28 +192,16 @@ void EwaldDipole::init()
// stats
if (comm->me == 0) {
if (screen) {
fprintf(screen," G vector (1/distance) = %g\n",g_ewald);
fprintf(screen," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(screen," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(screen," KSpace vectors: actual max1d max3d = %d %d %d\n",
kcount,kmax,kmax3d);
fprintf(screen," kxmax kymax kzmax = %d %d %d\n",
kxmax,kymax,kzmax);
}
if (logfile) {
fprintf(logfile," G vector (1/distance) = %g\n",g_ewald);
fprintf(logfile," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(logfile," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(logfile," KSpace vectors: actual max1d max3d = %d %d %d\n",
kcount,kmax,kmax3d);
fprintf(logfile," kxmax kymax kzmax = %d %d %d\n",
kxmax,kymax,kzmax);
}
std::string mesg = fmt::format(" G vector (1/distance) = {}\n",g_ewald);
mesg += fmt::format(" estimated absolute RMS force accuracy = {}\n",
estimated_accuracy);
mesg += fmt::format(" estimated relative force accuracy = {}\n",
estimated_accuracy/two_charge_force);
mesg += fmt::format(" KSpace vectors: actual max1d max3d = {} {} {}\n",
kcount,kmax,kmax3d);
mesg += fmt::format(" kxmax kymax kzmax = {} {} {}\n",
kxmax,kymax,kzmax);
utils::logmesg(lmp,mesg);
}
}

40
src/KSPACE/ewald_dipole_spin.cpp
View File

@ -18,6 +18,7 @@
#include "ewald_dipole_spin.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include <cmath>
#include "atom.h"
#include "comm.h"
@ -28,6 +29,8 @@
#include "memory.h"
#include "error.h"
#include "update.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace MathConst;
@ -61,10 +64,7 @@ EwaldDipoleSpin::~EwaldDipoleSpin() {}
void EwaldDipoleSpin::init()
{
if (comm->me == 0) {
if (screen) fprintf(screen,"EwaldDipoleSpin initialization ...\n");
if (logfile) fprintf(logfile,"EwaldDipoleSpin initialization ...\n");
}
if (comm->me == 0) utils::logmesg(lmp,"EwaldDipoleSpin initialization ...\n");
// error check
@ -182,28 +182,16 @@ void EwaldDipoleSpin::init()
// stats
if (comm->me == 0) {
if (screen) {
fprintf(screen," G vector (1/distance) = %g\n",g_ewald);
fprintf(screen," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(screen," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(screen," KSpace vectors: actual max1d max3d = %d %d %d\n",
kcount,kmax,kmax3d);
fprintf(screen," kxmax kymax kzmax = %d %d %d\n",
kxmax,kymax,kzmax);
}
if (logfile) {
fprintf(logfile," G vector (1/distance) = %g\n",g_ewald);
fprintf(logfile," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(logfile," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(logfile," KSpace vectors: actual max1d max3d = %d %d %d\n",
kcount,kmax,kmax3d);
fprintf(logfile," kxmax kymax kzmax = %d %d %d\n",
kxmax,kymax,kzmax);
}
std::string mesg = fmt::format(" G vector (1/distance) = {}\n",g_ewald);
mesg += fmt::format(" estimated absolute RMS force accuracy = {}\n",
estimated_accuracy);
mesg += fmt::format(" estimated relative force accuracy = {}\n",
estimated_accuracy/two_charge_force);
mesg += fmt::format(" KSpace vectors: actual max1d max3d = {} {} {}\n",
kcount,kmax,kmax3d);
mesg += fmt::format(" kxmax kymax kzmax = {} {} {}\n",
kxmax,kymax,kzmax);
utils::logmesg(lmp,mesg);
}
}

29
src/KSPACE/ewald_disp.cpp
View File

@ -18,6 +18,7 @@
#include "ewald_disp.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include <cmath>
#include "math_vector.h"
#include "math_const.h"
@ -30,6 +31,8 @@
#include "memory.h"
#include "error.h"
#include "update.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace MathConst;
@ -89,10 +92,7 @@ void EwaldDisp::init()
nbox = -1;
bytes = 0.0;
if (!comm->me) {
if (screen) fprintf(screen,"EwaldDisp initialization ...\n");
if (logfile) fprintf(logfile,"EwaldDisp initialization ...\n");
}
if (!comm->me) utils::logmesg(lmp,"EwaldDisp initialization ...\n");
triclinic_check();
if (domain->dimension == 2)
@ -169,11 +169,9 @@ void EwaldDisp::init()
if (qsqsum == 0.0 && bsbsum == 0.0 && M2 == 0.0)
error->all(FLERR,"Cannot use Ewald/disp solver "
"on system with no charge, dipole, or LJ particles");
if (fabs(qsum) > SMALL && comm->me == 0) {
char str[128];
sprintf(str,"System is not charge neutral, net charge = %g",qsum);
error->warning(FLERR,str);
}
if (fabs(qsum) > SMALL && comm->me == 0)
error->warning(FLERR,fmt::format("System is not charge neutral, "
"net charge = {}",qsum));
if (!function[1] && !function[2]) dispersionflag = 0;
if (!function[3]) dipoleflag = 0;
@ -229,10 +227,9 @@ void EwaldDisp::init()
}
}
if (!comm->me) {
if (screen) fprintf(screen, " G vector = %g, accuracy = %g\n", g_ewald,accuracy);
if (logfile) fprintf(logfile, " G vector = %g accuracy = %g\n", g_ewald,accuracy);
}
if (!comm->me)
utils::logmesg(lmp,fmt::format(" G vector = {}, accuracy = {}\n",
g_ewald,accuracy));
g_ewald_6 = g_ewald;
deallocate_peratom();
@ -294,10 +291,8 @@ void EwaldDisp::setup()
if (!(first_output||comm->me)) {
first_output = 1;
if (screen) fprintf(screen,
" vectors: nbox = %d, nkvec = %d\n", nbox, nkvec);
if (logfile) fprintf(logfile,
" vectors: nbox = %d, nkvec = %d\n", nbox, nkvec);
utils::logmesg(lmp,fmt::format(" vectors: nbox = {}, nkvec = {}\n",
nbox, nkvec));
}
}

33
src/KSPACE/fix_tune_kspace.cpp
View File

@ -19,6 +19,8 @@
#include <cmath>
#include <cstring>
#include <limits>
#include <string>
#include "comm.h"
#include "update.h"
#include "force.h"
#include "kspace.h"
@ -29,6 +31,8 @@
#include "neighbor.h"
#include "modify.h"
#include "compute.h"
#include "utils.h"
#include "fmt/format.h"
#define SWAP(a,b) {temp=(a);(a)=(b);(b)=temp;}
#define SIGN(a,b) ((b) >= 0.0 ? fabs(a) : -fabs(a))
@ -145,11 +149,12 @@ void FixTuneKspace::pre_exchange()
update_pair_style(new_pair_style,pair_cut_coul);
update_kspace_style(new_kspace_style,new_acc_str);
} else if (niter == 4) {
store_old_kspace_settings();
if (screen) fprintf(screen,"ewald_time = %g\npppm_time = %g\nmsm_time = %g",
ewald_time, pppm_time, msm_time);
if (logfile) fprintf(logfile,"ewald_time = %g\npppm_time = %g\nmsm_time = %g",
ewald_time, pppm_time, msm_time);
store_old_kspace_settings();
if (comm->me == 0)
utils::logmesg(lmp,fmt::format("ewald_time = {}\n"
"pppm_time = {}\n"
"msm_time = {}\n",
ewald_time, pppm_time, msm_time));
// switch to fastest one
strcpy(new_kspace_style,"ewald");
snprintf(new_pair_style,64,"%s/long",base_pair_style);
@ -240,8 +245,8 @@ void FixTuneKspace::update_pair_style(char *new_pair_style,
p_pair_settings_file = tmpfile();
force->pair->write_restart(p_pair_settings_file);
rewind(p_pair_settings_file);
if (screen) fprintf(screen,"Creating new pair style: %s\n",new_pair_style);
if (logfile) fprintf(logfile,"Creating new pair style: %s\n",new_pair_style);
if (comm->me == 0)
utils::logmesg(lmp,fmt::format("Creating new pair style: {}\n",new_pair_style));
// delete old pair style and create new one
force->create_pair(new_pair_style,1);
@ -250,8 +255,9 @@ void FixTuneKspace::update_pair_style(char *new_pair_style,
double *pcutoff = (double *) force->pair->extract("cut_coul",itmp);
double current_cutoff = *pcutoff;
if (screen) fprintf(screen,"Coulomb cutoff for real space: %g\n", current_cutoff);
if (logfile) fprintf(logfile,"Coulomb cutoff for real space: %g\n", current_cutoff);
if (comm->me == 0)
utils::logmesg(lmp,fmt::format("Coulomb cutoff for real space: {}\n",
current_cutoff));
// close temporary file
fclose(p_pair_settings_file);
@ -319,8 +325,9 @@ void FixTuneKspace::adjust_rcut(double time)
int itmp;
double *p_cutoff = (double *) force->pair->extract("cut_coul",itmp);
double current_cutoff = *p_cutoff;
if (screen) fprintf(screen,"Old Coulomb cutoff for real space: %g\n",current_cutoff);
if (logfile) fprintf(logfile,"Old Coulomb cutoff for real space: %g\n",current_cutoff);
if (comm->me == 0)
utils::logmesg(lmp,fmt::format("Old Coulomb cutoff for real space: {}\n",
current_cutoff));
// use Brent's method from Numerical Recipes to find optimal real space cutoff
@ -390,8 +397,8 @@ void FixTuneKspace::adjust_rcut(double time)
// report the new cutoff
double *new_cutoff = (double *) force->pair->extract("cut_coul",itmp);
current_cutoff = *new_cutoff;
if (screen) fprintf(screen,"Adjusted Coulomb cutoff for real space: %g\n", current_cutoff);
if (logfile) fprintf(logfile,"Adjusted Coulomb cutoff for real space: %g\n", current_cutoff);
if (comm->me == 0)
utils::logmesg(lmp,fmt::format("Adjusted Coulomb cutoff for real space: {}\n", current_cutoff));
store_old_kspace_settings();
update_pair_style(new_pair_style,pair_cut_coul);

22
src/KSPACE/msm_cg.cpp
View File

@ -26,6 +26,8 @@
#include "force.h"
#include "neighbor.h"
#include "memory.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
@ -140,20 +142,12 @@ void MSMCG::compute(int eflag, int vflag)
charged_frac = static_cast<double>(charged_all) * 100.0
/ static_cast<double>(atom->natoms);
if (me == 0) {
if (screen)
fprintf(screen,
" MSM/cg optimization cutoff: %g\n"
" Total charged atoms: %.1f%%\n"
" Min/max charged atoms/proc: %.1f%% %.1f%%\n",
smallq,charged_frac,charged_fmin,charged_fmax);
if (logfile)
fprintf(logfile,
" MSM/cg optimization cutoff: %g\n"
" Total charged atoms: %.1f%%\n"
" Min/max charged atoms/proc: %.1f%% %.1f%%\n",
smallq,charged_frac,charged_fmin,charged_fmax);
}
if (me == 0)
utils::logmesg(lmp,fmt::format(" MSM/cg optimization cutoff: {}\n"
" Total charged atoms: {:.1f}%\n"
" Min/max charged atoms/proc: {:.1f}%"
" {:.1f}%\n",smallq,
charged_frac,charged_fmin,charged_fmax));
}
// only need to rebuild this list after a neighbor list update

12
src/KSPACE/pppm.cpp
View File

@ -220,11 +220,8 @@ void PPPM::init()
error->all(FLERR,"Incorrect boundaries with slab PPPM");
}
if (order < 2 || order > MAXORDER) {
char str[128];
sprintf(str,"PPPM order cannot be < 2 or > than %d",MAXORDER);
error->all(FLERR,str);
}
if (order < 2 || order > MAXORDER)
error->all(FLERR,fmt::format("PPPM order cannot be < 2 or > {}",MAXORDER));
// compute two charge force
@ -1275,10 +1272,7 @@ void PPPM::adjust_gewald()
g_ewald -= dx;
if (fabs(newton_raphson_f()) < SMALL) return;
}
char str[128];
sprintf(str, "Could not compute g_ewald");
error->all(FLERR, str);
error->all(FLERR, "Could not compute g_ewald");
}
/* ----------------------------------------------------------------------

22
src/KSPACE/pppm_cg.cpp
View File

@ -28,6 +28,8 @@
#include "memory.h"
#include "math_const.h"
#include "remap.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace MathConst;
@ -151,20 +153,12 @@ void PPPMCG::compute(int eflag, int vflag)
charged_frac = static_cast<double>(charged_all) * 100.0
/ static_cast<double>(atom->natoms);
if (me == 0) {
if (screen)
fprintf(screen,
" PPPM/cg optimization cutoff: %g\n"
" Total charged atoms: %.1f%%\n"
" Min/max charged atoms/proc: %.1f%% %.1f%%\n",
smallq,charged_frac,charged_fmin,charged_fmax);
if (logfile)
fprintf(logfile,
" PPPM/cg optimization cutoff: %g\n"
" Total charged atoms: %.1f%%\n"
" Min/max charged atoms/proc: %.1f%% %.1f%%\n",
smallq,charged_frac,charged_fmin,charged_fmax);
}
if (me == 0)
utils::logmesg(lmp,fmt::format(" PPPM/cg optimization cutoff: {}\n"
" Total charged atoms: {:.1f}%\n"
" Min/max charged atoms/proc: {:.1f}%"
" {:.1f}%\n",smallq,
charged_frac,charged_fmin,charged_fmax));
}
// only need to rebuild this list after a neighbor list update

58
src/KSPACE/pppm_dipole.cpp
View File

@ -18,6 +18,7 @@
#include "pppm_dipole.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include <cmath>
#include "atom.h"
#include "comm.h"
@ -30,6 +31,8 @@
#include "memory.h"
#include "error.h"
#include "update.h"
#include "utils.h"
#include "fmt/format.h"
#include "math_const.h"
#include "math_special.h"
@ -102,10 +105,7 @@ PPPMDipole::~PPPMDipole()
void PPPMDipole::init()
{
if (me == 0) {
if (screen) fprintf(screen,"PPPMDipole initialization ...\n");
if (logfile) fprintf(logfile,"PPPMDipole initialization ...\n");
}
if (me == 0) utils::logmesg(lmp,"PPPMDipole initialization ...\n");
// error check
@ -143,11 +143,9 @@ void PPPMDipole::init()
error->all(FLERR,"Incorrect boundaries with slab PPPMDipole");
}
if (order < 2 || order > MAXORDER) {
char str[128];
sprintf(str,"PPPMDipole order cannot be < 2 or > than %d",MAXORDER);
error->all(FLERR,str);
}
if (order < 2 || order > MAXORDER)
error->all(FLERR,fmt::format("PPPMDipole order cannot be < 2 or > {}",
MAXORDER));
// compute two charge force
@ -246,37 +244,17 @@ void PPPMDipole::init()
MPI_Allreduce(&nfft_both,&nfft_both_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
#ifdef FFT_SINGLE
const char fft_prec[] = "single";
#else
const char fft_prec[] = "double";
#endif
if (screen) {
fprintf(screen," G vector (1/distance) = %g\n",g_ewald);
fprintf(screen," grid = %d %d %d\n",nx_pppm,ny_pppm,nz_pppm);
fprintf(screen," stencil order = %d\n",order);
fprintf(screen," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(screen," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(screen," using %s precision FFTs\n",fft_prec);
fprintf(screen," 3d grid and FFT values/proc = %d %d\n",
ngrid_max,nfft_both_max);
}
if (logfile) {
fprintf(logfile," G vector (1/distance) = %g\n",g_ewald);
fprintf(logfile," grid = %d %d %d\n",nx_pppm,ny_pppm,nz_pppm);
fprintf(logfile," stencil order = %d\n",order);
fprintf(logfile," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(logfile," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(logfile," using %s precision FFTs\n",fft_prec);
fprintf(logfile," 3d grid and FFT values/proc = %d %d\n",
ngrid_max,nfft_both_max);
}
std::string mesg = fmt::format(" G vector (1/distance) = {}\n",g_ewald);
mesg += fmt::format(" grid = {} {} {}\n",nx_pppm,ny_pppm,nz_pppm);
mesg += fmt::format(" stencil order = {}\n",order);
mesg += fmt::format(" estimated absolute RMS force accuracy = {}\n",
estimated_accuracy);
mesg += fmt::format(" estimated relative force accuracy = {}\n",
estimated_accuracy/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);
}
// allocate K-space dependent memory

58
src/KSPACE/pppm_dipole_spin.cpp
View File

@ -18,6 +18,7 @@
#include "pppm_dipole_spin.h"
#include <mpi.h>
#include <cstring>
#include <string>
#include "atom.h"
#include "comm.h"
#include "gridcomm.h"
@ -27,6 +28,8 @@
#include "memory.h"
#include "error.h"
#include "update.h"
#include "utils.h"
#include "fmt/format.h"
#include "math_const.h"
@ -87,10 +90,7 @@ PPPMDipoleSpin::~PPPMDipoleSpin()
void PPPMDipoleSpin::init()
{
if (me == 0) {
if (screen) fprintf(screen,"PPPMDipoleSpin initialization ...\n");
if (logfile) fprintf(logfile,"PPPMDipoleSpin initialization ...\n");
}
if (me == 0) utils::logmesg(lmp,"PPPMDipoleSpin initialization ...\n");
// error check
@ -123,11 +123,9 @@ void PPPMDipoleSpin::init()
error->all(FLERR,"Incorrect boundaries with slab PPPMDipoleSpin");
}
if (order < 2 || order > MAXORDER) {
char str[128];
sprintf(str,"PPPMDipoleSpin order cannot be < 2 or > than %d",MAXORDER);
error->all(FLERR,str);
}
if (order < 2 || order > MAXORDER)
error->all(FLERR,fmt::format("PPPMDipoleSpin order cannot be < 2 or > {}",
MAXORDER));
// compute two charge force
@ -226,37 +224,17 @@ void PPPMDipoleSpin::init()
MPI_Allreduce(&nfft_both,&nfft_both_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
#ifdef FFT_SINGLE
const char fft_prec[] = "single";
#else
const char fft_prec[] = "double";
#endif
if (screen) {
fprintf(screen," G vector (1/distance) = %g\n",g_ewald);
fprintf(screen," grid = %d %d %d\n",nx_pppm,ny_pppm,nz_pppm);
fprintf(screen," stencil order = %d\n",order);
fprintf(screen," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(screen," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(screen," using %s precision FFTs\n",fft_prec);
fprintf(screen," 3d grid and FFT values/proc = %d %d\n",
ngrid_max,nfft_both_max);
}
if (logfile) {
fprintf(logfile," G vector (1/distance) = %g\n",g_ewald);
fprintf(logfile," grid = %d %d %d\n",nx_pppm,ny_pppm,nz_pppm);
fprintf(logfile," stencil order = %d\n",order);
fprintf(logfile," estimated absolute RMS force accuracy = %g\n",
estimated_accuracy);
fprintf(logfile," estimated relative force accuracy = %g\n",
estimated_accuracy/two_charge_force);
fprintf(logfile," using %s precision FFTs\n",fft_prec);
fprintf(logfile," 3d grid and FFT values/proc = %d %d\n",
ngrid_max,nfft_both_max);
}
std::string mesg = fmt::format(" G vector (1/distance) = {}\n",g_ewald);
mesg += fmt::format(" grid = {} {} {}\n",nx_pppm,ny_pppm,nz_pppm);
mesg += fmt::format(" stencil order = {}\n",order);
mesg += fmt::format(" estimated absolute RMS force accuracy = {}\n",
estimated_accuracy);
mesg += fmt::format(" estimated relative force accuracy = {}\n",
estimated_accuracy/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);
}
// allocate K-space dependent memory

167
src/KSPACE/pppm_disp.cpp
View File

@ -20,6 +20,7 @@
#include <mpi.h>
#include <cstring>
#include <cmath>
#include <string>
#include "math_const.h"
#include "atom.h"
#include "comm.h"
@ -34,6 +35,8 @@
#include "remap_wrap.h"
#include "memory.h"
#include "error.h"
#include "utils.h"
#include "fmt/format.h"
using namespace LAMMPS_NS;
using namespace MathConst;
@ -254,10 +257,7 @@ PPPMDisp::~PPPMDisp()
void PPPMDisp::init()
{
if (me == 0) {
if (screen) fprintf(screen,"PPPMDisp initialization ...\n");
if (logfile) fprintf(logfile,"PPPMDisp initialization ...\n");
}
if (me == 0) utils::logmesg(lmp,"PPPMDisp initialization ...\n");
// error check
@ -277,11 +277,9 @@ void PPPMDisp::init()
error->all(FLERR,"Incorrect boundaries with slab PPPMDisp");
}
if (order > MAXORDER || order_6 > MAXORDER) {
char str[128];
sprintf(str,"PPPMDisp coulomb order cannot be greater than %d",MAXORDER);
error->all(FLERR,str);
}
if (order > MAXORDER || order_6 > MAXORDER)
error->all(FLERR,fmt::format("PPPMDisp coulomb or dispersion order cannot"
" be greater than {}",MAXORDER));
// compute two charge force
@ -330,9 +328,9 @@ void PPPMDisp::init()
else error->all(FLERR,"Unsupported mixing rule in kspace_style pppm/disp");
break;
default:
sprintf(str, "Unsupported order in kspace_style "
"pppm/disp, pair_style %s", force->pair_style);
error->all(FLERR,str);
error->all(FLERR,std::string("Unsupported order in kspace_style "
"pppm/disp, pair_style ")
+ force->pair_style);
}
function[k] = 1;
}
@ -340,11 +338,8 @@ void PPPMDisp::init()
// warn, if function[0] is not set but charge attribute is set!
if (!function[0] && atom->q_flag && me == 0) {
char str[128];
sprintf(str, "Charges are set, but coulombic solver is not used");
error->warning(FLERR, str);
}
if (!function[0] && atom->q_flag && me == 0)
error->warning(FLERR, "Charges are set, but coulombic solver is not used");
// show error message if pppm/disp is not used correctly
@ -481,31 +476,19 @@ void PPPMDisp::init()
MPI_Allreduce(&nfft_both,&nfft_both_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
if (screen) {
fprintf(screen," Coulomb G vector (1/distance)= %g\n",g_ewald);
fprintf(screen," Coulomb grid = %d %d %d\n",nx_pppm,ny_pppm,nz_pppm);
fprintf(screen," Coulomb stencil order = %d\n",order);
fprintf(screen," Coulomb estimated absolute RMS force accuracy = %g\n",
acc);
fprintf(screen," Coulomb estimated relative force accuracy = %g\n",
acc/two_charge_force);
fprintf(screen," using " LMP_FFT_PREC " precision " LMP_FFT_LIB "\n");
fprintf(screen," 3d grid and FFT values/proc = %d %d\n",
ngrid_max, nfft_both_max);
}
if (logfile) {
fprintf(logfile," Coulomb G vector (1/distance) = %g\n",g_ewald);
fprintf(logfile," Coulomb grid = %d %d %d\n",nx_pppm,ny_pppm,nz_pppm);
fprintf(logfile," Coulomb stencil order = %d\n",order);
fprintf(logfile,
" Coulomb estimated absolute RMS force accuracy = %g\n",
acc);
fprintf(logfile," Coulomb estimated relative force accuracy = %g\n",
acc/two_charge_force);
fprintf(logfile," using " LMP_FFT_PREC " precision " LMP_FFT_LIB "\n");
fprintf(logfile," 3d grid and FFT values/proc = %d %d\n",
ngrid_max, nfft_both_max);
}
std::string mesg = fmt::format(" Coulomb G vector (1/distance)= {}\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 "
"= {}\n",acc);
mesg += fmt::format(" Coulomb estimated relative force accuracy = {}\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);
}
}
@ -573,46 +556,19 @@ void PPPMDisp::init()
MPI_Allreduce(&nfft_both_6,&nfft_both_max,1,MPI_INT,MPI_MAX,world);
if (me == 0) {
#ifdef FFT_SINGLE
const char fft_prec[] = "single";
#else
const char fft_prec[] = "double";
#endif
if (screen) {
fprintf(screen," Dispersion G vector (1/distance)= %g\n",g_ewald_6);
fprintf(screen," Dispersion grid = %d %d %d\n",
nx_pppm_6,ny_pppm_6,nz_pppm_6);
fprintf(screen," Dispersion stencil order = %d\n",order_6);
fprintf(screen," Dispersion estimated absolute "
"RMS force accuracy = %g\n",acc);
fprintf(screen," Dispersion estimated absolute "
"real space RMS force accuracy = %g\n",acc_real);
fprintf(screen," Dispersion estimated absolute "
"kspace RMS force accuracy = %g\n",acc_kspace);
fprintf(screen," Dispersion estimated relative force accuracy = %g\n",
acc/two_charge_force);
fprintf(screen," using %s precision FFTs\n",fft_prec);
fprintf(screen," 3d grid and FFT values/proc dispersion = %d %d\n",
ngrid_max,nfft_both_max);
}
if (logfile) {
fprintf(logfile," Dispersion G vector (1/distance) = %g\n",g_ewald_6);
fprintf(logfile," Dispersion grid = %d %d %d\n",
nx_pppm_6,ny_pppm_6,nz_pppm_6);
fprintf(logfile," Dispersion stencil order = %d\n",order_6);
fprintf(logfile," Dispersion estimated absolute "
"RMS force accuracy = %g\n",acc);
fprintf(logfile," Dispersion estimated absolute "
"real space RMS force accuracy = %g\n",acc_real);
fprintf(logfile," Dispersion estimated absolute "
"kspace RMS force accuracy = %g\n",acc_kspace);
fprintf(logfile," Disperion estimated relative force accuracy = %g\n",
acc/two_charge_force);
fprintf(logfile," using %s precision FFTs\n",fft_prec);
fprintf(logfile," 3d grid and FFT values/proc dispersion = %d %d\n",
ngrid_max,nfft_both_max);
}
std::string mesg = fmt::format(" Dispersion G vector (1/distance)= "
"{}\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 "
"= {}\n",acc);
mesg += fmt::format(" Dispersion estimated relative force accuracy "
"= {}\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);
}
}
@ -1285,10 +1241,8 @@ void PPPMDisp::init_coeffs() // local pair coeffs
delete [] B;
B = NULL;
if (function[3] + function[2]) { // no mixing rule or arithmetic
if (function[2] && me == 0) {
if (screen) fprintf(screen," Optimizing splitting of Dispersion coefficients\n");
if (logfile) fprintf(logfile," Optimizing splitting of Dispersion coefficients\n");
}
if (function[2] && me == 0)
utils::logmesg(lmp," Optimizing splitting of Dispersion coefficients\n");
// allocate data for eigenvalue decomposition
double **A=NULL;
@ -1349,11 +1303,9 @@ void PPPMDisp::init_coeffs() // local pair coeffs
}
err = bmax/amax;
if (err > 1.0e-4) {
char str[128];
sprintf(str,"Estimated error in splitting of dispersion coeffs is %g",err);
error->warning(FLERR, str);
}
if (err > 1.0e-4 && comm->me == 0)
error->warning(FLERR,fmt::format("Estimated error in splitting of "
"dispersion coeffs is {}",err));
// set B
B = new double[nsplit*n+nsplit];
for (int i = 0; i< nsplit; i++) {
@ -1374,31 +1326,24 @@ void PPPMDisp::init_coeffs() // local pair coeffs
function[3] = 0;
function[2] = 0;
function[1] = 1;
if (me == 0) {
if (screen) fprintf(screen," Using geometric mixing for reciprocal space\n");
if (logfile) fprintf(logfile," Using geometric mixing for reciprocal space\n");
}
if (me == 0)
utils::logmesg(lmp," Using geometric mixing for reciprocal space\n");
}
if (function[2] && nsplit <= 6) {
if (me == 0) {
if (screen) fprintf(screen," Using %d instead of 7 structure factors\n",nsplit);
if (logfile) fprintf(logfile," Using %d instead of 7 structure factors\n",nsplit);
}
if (me == 0)
utils::logmesg(lmp,fmt::format(" Using {} instead of 7 structure "
"factors\n",nsplit));
function[3] = 1;
function[2] = 0;
}
if (function[2] && (nsplit > 6)) {
if (me == 0) {
if (screen) fprintf(screen," Using 7 structure factors\n");
if (logfile) fprintf(logfile," Using 7 structure factors\n");
}
if (me == 0) utils::logmesg(lmp," Using 7 structure factors\n");
if ( B ) delete [] B;
}
if (function[3]) {
if (me == 0) {
if (screen) fprintf(screen," Using %d structure factors\n",nsplit);
if (logfile) fprintf(logfile," Using %d structure factors\n",nsplit);
}
if (me == 0)
utils::logmesg(lmp,fmt::format(" Using {} structure factors\n",
nsplit));
if (nsplit > 9) error->warning(FLERR, "Simulations might be very slow because of large number of structure factors");
}
@ -2884,9 +2829,7 @@ void PPPMDisp::adjust_gewald()
// Failed to converge
char str[128];
sprintf(str, "Could not compute g_ewald");
error->all(FLERR, str);
error->all(FLERR, "Could not compute g_ewald");
}
@ -3517,9 +3460,7 @@ void PPPMDisp::adjust_gewald_6()
// Failed to converge
char str[128];
sprintf(str, "Could not adjust g_ewald_6");
error->all(FLERR, str);
error->all(FLERR, "Could not adjust g_ewald_6");
}