diff --git a/src/integrate.cpp b/src/integrate.cpp index 3343337185..beef03c302 100644 --- a/src/integrate.cpp +++ b/src/integrate.cpp @@ -14,6 +14,9 @@ #include "stdlib.h" #include "integrate.h" #include "update.h" +#include "force.h" +#include "pair.h" +#include "kspace.h" #include "modify.h" #include "compute.h" @@ -38,6 +41,18 @@ Integrate::~Integrate() delete [] vlist_atom; } +/* ---------------------------------------------------------------------- */ + +void Integrate::init() +{ + // allow pair and Kspace compute() to be turned off via modify flags + + if (force->pair && force->pair->compute_flag) pair_compute_flag = 1; + else pair_compute_flag = 0; + if (force->kspace && force->kspace->compute_flag) kspace_compute_flag = 1; + else kspace_compute_flag = 0; +} + /* ---------------------------------------------------------------------- setup lists of computes for global and per-atom PE and pressure ------------------------------------------------------------------------- */ diff --git a/src/integrate.h b/src/integrate.h index 853390ae12..bc5a43fd9e 100644 --- a/src/integrate.h +++ b/src/integrate.h @@ -22,7 +22,7 @@ class Integrate : protected Pointers { public: Integrate(class LAMMPS *, int, char **); virtual ~Integrate(); - virtual void init() = 0; + virtual void init(); virtual void setup() = 0; virtual void setup_minimal(int) = 0; virtual void run(int) = 0; @@ -42,6 +42,9 @@ class Integrate : protected Pointers { class Compute **vlist_global; class Compute **vlist_atom; + int pair_compute_flag; // 0 if pair->compute is skipped + int kspace_compute_flag; // 0 if kspace->compute is skipped + void ev_setup(); void ev_set(bigint); }; diff --git a/src/kspace.cpp b/src/kspace.cpp index 9913b44809..4c527f887f 100644 --- a/src/kspace.cpp +++ b/src/kspace.cpp @@ -28,6 +28,9 @@ using namespace LAMMPS_NS; KSpace::KSpace(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp) { energy = 0.0; + virial[0] = virial[1] = virial[2] = virial[3] = virial[4] = virial[5] = 0.0; + + compute_flag = 1; order = 5; gridflag = 0; gewaldflag = 0; @@ -53,6 +56,15 @@ KSpace::~KSpace() memory->destroy(vatom); } +/* ---------------------------------------------------------------------- */ + +void KSpace::compute_dummy(int eflag, int vflag) +{ + if (eflag || vflag) ev_setup(eflag,vflag); + else evflag = evflag_atom = eflag_global = vflag_global = + eflag_atom = vflag_atom = 0; +} + /* ---------------------------------------------------------------------- setup for energy, virial computation see integrate::ev_set() for values of eflag (0-3) and vflag (0-6) @@ -149,6 +161,12 @@ void KSpace::modify_params(int narg, char **arg) error->warning(FLERR,"Kspace_modify slab param < 2.0 may " "cause unphysical behavior"); slabflag = 1; + } else if (strcmp(arg[iarg],"compute") == 0) { + if (iarg+2 > narg) error->all(FLERR,"Illegal kspace_modify command"); + if (strcmp(arg[iarg+1],"yes") == 0) compute_flag = 1; + else if (strcmp(arg[iarg+1],"no") == 0) compute_flag = 0; + else error->all(FLERR,"Illegal kspace_modify command"); + iarg += 2; } else error->all(FLERR,"Illegal kspace_modify command"); } } diff --git a/src/kspace.h b/src/kspace.h index 30605b4193..78c8b477ee 100644 --- a/src/kspace.h +++ b/src/kspace.h @@ -29,10 +29,15 @@ class KSpace : protected Pointers { double g_ewald; int nx_pppm,ny_pppm,nz_pppm; + int compute_flag; // 0 if skip compute() + KSpace(class LAMMPS *, int, char **); virtual ~KSpace(); void modify_params(int, char **); void *extract(const char *); + void compute_dummy(int, int); + + // general child-class methods virtual void init() = 0; virtual void setup() = 0; diff --git a/src/min.cpp b/src/min.cpp index 48e297b3fa..2f9a461ba3 100644 --- a/src/min.cpp +++ b/src/min.cpp @@ -151,6 +151,13 @@ void Min::init() rho_flag = 0; if (atom->rho_flag) rho_flag = 1; + // allow pair and Kspace compute() to be turned off via modify flags + + if (force->pair && force->pair->compute_flag) pair_compute_flag = 1; + else pair_compute_flag = 0; + if (force->kspace && force->kspace->compute_flag) kspace_compute_flag = 1; + else kspace_compute_flag = 0; + // orthogonal vs triclinic simulation box triclinic = domain->triclinic; @@ -172,10 +179,6 @@ void Min::init() neighbor->dist_check = 1; niter = neval = 0; - - // style-specific initialization - - init_style(); } /* ---------------------------------------------------------------------- @@ -250,7 +253,8 @@ void Min::setup() force_clear(); modify->setup_pre_force(vflag); - if (force->pair) force->pair->compute(eflag,vflag); + if (pair_compute_flag) force->pair->compute(eflag,vflag); + else if (force->pair) force->pair->compute_dummy(eflag,vflag); if (atom->molecular) { if (force->bond) force->bond->compute(eflag,vflag); @@ -261,7 +265,8 @@ void Min::setup() if (force->kspace) { force->kspace->setup(); - force->kspace->compute(eflag,vflag); + if (kspace_compute_flag) force->kspace->compute(eflag,vflag); + else force->kspace->compute_dummy(eflag,vflag); } if (force->newton) comm->reverse_comm(); @@ -321,7 +326,8 @@ void Min::setup_minimal(int flag) force_clear(); modify->setup_pre_force(vflag); - if (force->pair) force->pair->compute(eflag,vflag); + if (pair_compute_flag) force->pair->compute(eflag,vflag); + else if (force->pair) force->pair->compute_dummy(eflag,vflag); if (atom->molecular) { if (force->bond) force->bond->compute(eflag,vflag); @@ -332,7 +338,8 @@ void Min::setup_minimal(int flag) if (force->kspace) { force->kspace->setup(); - force->kspace->compute(eflag,vflag); + if (kspace_compute_flag) force->kspace->compute(eflag,vflag); + else force->kspace->compute_dummy(eflag,vflag); } if (force->newton) comm->reverse_comm(); @@ -458,7 +465,7 @@ double Min::energy_force(int resetflag) timer->stamp(); - if (force->pair) { + if (pair_compute_flag) { force->pair->compute(eflag,vflag); timer->stamp(TIME_PAIR); } @@ -471,7 +478,7 @@ double Min::energy_force(int resetflag) timer->stamp(TIME_BOND); } - if (force->kspace) { + if (kspace_compute_flag) { force->kspace->compute(eflag,vflag); timer->stamp(TIME_KSPACE); } diff --git a/src/min.h b/src/min.h index fde35f5c80..bb21df0580 100644 --- a/src/min.h +++ b/src/min.h @@ -30,7 +30,7 @@ class Min : protected Pointers { Min(class LAMMPS *); virtual ~Min(); - void init(); + virtual void init(); void setup(); void setup_minimal(int); void run(int); @@ -66,6 +66,9 @@ class Min : protected Pointers { int torqueflag,erforceflag; int e_flag,rho_flag; + int pair_compute_flag; // 0 if pair->compute is skipped + int kspace_compute_flag; // 0 if kspace->compute is skipped + int narray; // # of arrays stored by fix_minimize class FixMinimize *fix_minimize; // fix that stores auxiliary data diff --git a/src/min_fire.cpp b/src/min_fire.cpp index 18081f17b5..c061eceeda 100644 --- a/src/min_fire.cpp +++ b/src/min_fire.cpp @@ -44,8 +44,10 @@ MinFire::MinFire(LAMMPS *lmp) : Min(lmp) {} /* ---------------------------------------------------------------------- */ -void MinFire::init_style() +void MinFire::init() { + Min::init(); + dt = update->dt; } diff --git a/src/min_fire.h b/src/min_fire.h index 0160ab4afb..926b1e6ecc 100644 --- a/src/min_fire.h +++ b/src/min_fire.h @@ -28,7 +28,7 @@ class MinFire : public Min { public: MinFire(class LAMMPS *); ~MinFire() {} - void init_style(); + void init(); void setup_style(); void reset_vectors(); int iterate(int); diff --git a/src/min_hftn.cpp b/src/min_hftn.cpp index 7c1c851f45..60e4890bd4 100644 --- a/src/min_hftn.cpp +++ b/src/min_hftn.cpp @@ -104,11 +104,13 @@ MinHFTN::~MinHFTN (void) } /* ---------------------------------------------------------------------- - Public method init_style + Public method init ------------------------------------------------------------------------- */ -void MinHFTN::init_style() +void MinHFTN::init() { + Min::init(); + for (int i = 1; i < NUM_HFTN_ATOM_BASED_VECTORS; i++) { if (_daExtraGlobal[i] != NULL) delete [] _daExtraGlobal[i]; diff --git a/src/min_hftn.h b/src/min_hftn.h index 4cfaf40247..42c181930c 100644 --- a/src/min_hftn.h +++ b/src/min_hftn.h @@ -42,8 +42,7 @@ class MinHFTN : public Min MinHFTN (LAMMPS *); ~MinHFTN (void); - - void init_style(); + void init(); void setup_style(); void reset_vectors(); int iterate (int); diff --git a/src/min_linesearch.cpp b/src/min_linesearch.cpp index 73d9ef83f0..b50971a35b 100644 --- a/src/min_linesearch.cpp +++ b/src/min_linesearch.cpp @@ -77,8 +77,10 @@ MinLineSearch::~MinLineSearch() /* ---------------------------------------------------------------------- */ -void MinLineSearch::init_style() +void MinLineSearch::init() { + Min::init(); + if (linestyle == 0) linemin = &MinLineSearch::linemin_backtrack; else if (linestyle == 1) linemin = &MinLineSearch::linemin_quadratic; else if (linestyle == 2) linemin = &MinLineSearch::linemin_forcezero; diff --git a/src/min_linesearch.h b/src/min_linesearch.h index bbe325b5f4..55dbad17a6 100644 --- a/src/min_linesearch.h +++ b/src/min_linesearch.h @@ -22,7 +22,7 @@ class MinLineSearch : public Min { public: MinLineSearch(class LAMMPS *); ~MinLineSearch(); - void init_style(); + void init(); void setup_style(); void reset_vectors(); diff --git a/src/min_quickmin.cpp b/src/min_quickmin.cpp index fe9df160b5..80479827da 100644 --- a/src/min_quickmin.cpp +++ b/src/min_quickmin.cpp @@ -41,8 +41,10 @@ MinQuickMin::MinQuickMin(LAMMPS *lmp) : Min(lmp) {} /* ---------------------------------------------------------------------- */ -void MinQuickMin::init_style() +void MinQuickMin::init() { + Min::init(); + dt = update->dt; } diff --git a/src/min_quickmin.h b/src/min_quickmin.h index dcdc227761..a69c30f616 100644 --- a/src/min_quickmin.h +++ b/src/min_quickmin.h @@ -28,7 +28,7 @@ class MinQuickMin : public Min { public: MinQuickMin(class LAMMPS *); ~MinQuickMin() {} - void init_style(); + void init(); void setup_style(); void reset_vectors(); int iterate(int); diff --git a/src/pair.cpp b/src/pair.cpp index 234235211f..722feec102 100644 --- a/src/pair.cpp +++ b/src/pair.cpp @@ -64,6 +64,7 @@ Pair::Pair(LAMMPS *lmp) : Pointers(lmp) // pair_modify settings + compute_flag = 1; offset_flag = 0; mix_flag = GEOMETRIC; tail_flag = 0; @@ -127,6 +128,12 @@ void Pair::modify_params(int narg, char **arg) else if (strcmp(arg[iarg+1],"no") == 0) tail_flag = 0; else error->all(FLERR,"Illegal pair_modify command"); iarg += 2; + } else if (strcmp(arg[iarg],"compute") == 0) { + if (iarg+2 > narg) error->all(FLERR,"Illegal pair_modify command"); + if (strcmp(arg[iarg+1],"yes") == 0) compute_flag = 1; + else if (strcmp(arg[iarg+1],"no") == 0) compute_flag = 0; + else error->all(FLERR,"Illegal pair_modify command"); + iarg += 2; } else error->all(FLERR,"Illegal pair_modify command"); } } @@ -262,6 +269,14 @@ double Pair::mix_distance(double sig1, double sig2) return value; } +/* ---------------------------------------------------------------------- */ + +void Pair::compute_dummy(int eflag, int vflag) +{ + if (eflag || vflag) ev_setup(eflag,vflag); + else evflag = 0; +} + /* ---------------------------------------------------------------------- setup for energy, virial computation see integrate::ev_set() for values of eflag (0-3) and vflag (0-6) diff --git a/src/pair.h b/src/pair.h index 53fa3e9c43..d86b7005d7 100644 --- a/src/pair.h +++ b/src/pair.h @@ -72,6 +72,8 @@ class Pair : protected Pointers { class NeighList *listmiddle; class NeighList *listouter; + int compute_flag; // 0 if skip compute() + Pair(class LAMMPS *); virtual ~Pair(); @@ -84,6 +86,7 @@ class Pair : protected Pointers { void write_file(int, char **); void init_bitmap(double, double, int, int &, int &, int &, int &); virtual void modify_params(int, char **); + void compute_dummy(int, int); // need to be public, so can be called by pair_style reaxc diff --git a/src/respa.cpp b/src/respa.cpp index 3e54b58424..aedfa03027 100644 --- a/src/respa.cpp +++ b/src/respa.cpp @@ -244,6 +244,8 @@ Respa::~Respa() void Respa::init() { + Integrate::init(); + // warn if no fixes if (modify->nfix == 0 && comm->me == 0) @@ -359,17 +361,17 @@ void Respa::setup() force->dihedral->compute(eflag,vflag); if (level_improper == ilevel && force->improper) force->improper->compute(eflag,vflag); - if (level_pair == ilevel && force->pair) + if (level_pair == ilevel && pair_compute_flag) force->pair->compute(eflag,vflag); - if (level_inner == ilevel && force->pair) + if (level_inner == ilevel && pair_compute_flag) force->pair->compute_inner(); - if (level_middle == ilevel && force->pair) + if (level_middle == ilevel && pair_compute_flag) force->pair->compute_middle(); - if (level_outer == ilevel && force->pair) + if (level_outer == ilevel && pair_compute_flag) force->pair->compute_outer(eflag,vflag); if (level_kspace == ilevel && force->kspace) { force->kspace->setup(); - force->kspace->compute(eflag,vflag); + if (kspace_compute_flag) force->kspace->compute(eflag,vflag); } if (newton[ilevel]) comm->reverse_comm(); copy_f_flevel(ilevel); @@ -423,17 +425,17 @@ void Respa::setup_minimal(int flag) force->dihedral->compute(eflag,vflag); if (level_improper == ilevel && force->improper) force->improper->compute(eflag,vflag); - if (level_pair == ilevel && force->pair) + if (level_pair == ilevel && pair_compute_flag) force->pair->compute(eflag,vflag); - if (level_inner == ilevel && force->pair) + if (level_inner == ilevel && pair_compute_flag) force->pair->compute_inner(); - if (level_middle == ilevel && force->pair) + if (level_middle == ilevel && pair_compute_flag) force->pair->compute_middle(); - if (level_outer == ilevel && force->pair) + if (level_outer == ilevel && pair_compute_flag) force->pair->compute_outer(eflag,vflag); if (level_kspace == ilevel && force->kspace) { force->kspace->setup(); - force->kspace->compute(eflag,vflag); + if (kspace_compute_flag) force->kspace->compute(eflag,vflag); } if (newton[ilevel]) comm->reverse_comm(); copy_f_flevel(ilevel); @@ -557,23 +559,23 @@ void Respa::recurse(int ilevel) force->improper->compute(eflag,vflag); timer->stamp(TIME_BOND); } - if (level_pair == ilevel && force->pair) { + if (level_pair == ilevel && pair_compute_flag) { force->pair->compute(eflag,vflag); timer->stamp(TIME_PAIR); } - if (level_inner == ilevel && force->pair) { + if (level_inner == ilevel && pair_compute_flag) { force->pair->compute_inner(); timer->stamp(TIME_PAIR); } - if (level_middle == ilevel && force->pair) { + if (level_middle == ilevel && pair_compute_flag) { force->pair->compute_middle(); timer->stamp(TIME_PAIR); } - if (level_outer == ilevel && force->pair) { + if (level_outer == ilevel && pair_compute_flag) { force->pair->compute_outer(eflag,vflag); timer->stamp(TIME_PAIR); } - if (level_kspace == ilevel && force->kspace) { + if (level_kspace == ilevel && kspace_compute_flag) { force->kspace->compute(eflag,vflag); timer->stamp(TIME_KSPACE); } diff --git a/src/verlet.cpp b/src/verlet.cpp index 223f02c954..5190281173 100644 --- a/src/verlet.cpp +++ b/src/verlet.cpp @@ -46,6 +46,8 @@ Verlet::Verlet(LAMMPS *lmp, int narg, char **arg) : void Verlet::init() { + Integrate::init(); + // warn if no fixes if (modify->nfix == 0 && comm->me == 0) @@ -118,7 +120,8 @@ void Verlet::setup() force_clear(); modify->setup_pre_force(vflag); - if (force->pair) force->pair->compute(eflag,vflag); + if (pair_compute_flag) force->pair->compute(eflag,vflag); + else if (force->pair) force->pair->compute_dummy(eflag,vflag); if (atom->molecular) { if (force->bond) force->bond->compute(eflag,vflag); @@ -129,7 +132,8 @@ void Verlet::setup() if (force->kspace) { force->kspace->setup(); - force->kspace->compute(eflag,vflag); + if (kspace_compute_flag) force->kspace->compute(eflag,vflag); + else force->kspace->compute_dummy(eflag,vflag); } if (force->newton) comm->reverse_comm(); @@ -172,7 +176,8 @@ void Verlet::setup_minimal(int flag) force_clear(); modify->setup_pre_force(vflag); - if (force->pair) force->pair->compute(eflag,vflag); + if (pair_compute_flag) force->pair->compute(eflag,vflag); + else if (force->pair) force->pair->compute_dummy(eflag,vflag); if (atom->molecular) { if (force->bond) force->bond->compute(eflag,vflag); @@ -183,7 +188,8 @@ void Verlet::setup_minimal(int flag) if (force->kspace) { force->kspace->setup(); - force->kspace->compute(eflag,vflag); + if (kspace_compute_flag) force->kspace->compute(eflag,vflag); + else force->kspace->compute_dummy(eflag,vflag); } if (force->newton) comm->reverse_comm(); @@ -256,7 +262,7 @@ void Verlet::run(int n) timer->stamp(); - if (force->pair) { + if (pair_compute_flag) { force->pair->compute(eflag,vflag); timer->stamp(TIME_PAIR); } @@ -269,7 +275,7 @@ void Verlet::run(int n) timer->stamp(TIME_BOND); } - if (force->kspace) { + if (kspace_compute_flag) { force->kspace->compute(eflag,vflag); timer->stamp(TIME_KSPACE); }