diff --git a/src/angle_zero.cpp b/src/angle_zero.cpp index a63f787840d47e1fbdf007b9a92e77943ffcb363..375d5ad81375fb3b36e03e9a1899097d55a36544 100644 --- a/src/angle_zero.cpp +++ b/src/angle_zero.cpp @@ -17,17 +17,21 @@ #include <math.h> #include <stdlib.h> +#include <string.h> #include "angle_zero.h" #include "atom.h" #include "force.h" +#include "comm.h" +#include "math_const.h" #include "memory.h" #include "error.h" using namespace LAMMPS_NS; +using namespace MathConst; /* ---------------------------------------------------------------------- */ -AngleZero::AngleZero(LAMMPS *lmp) : Angle(lmp) {} +AngleZero::AngleZero(LAMMPS *lmp) : Angle(lmp), coeffflag(1) {} /* ---------------------------------------------------------------------- */ @@ -35,6 +39,7 @@ AngleZero::~AngleZero() { if (allocated && !copymode) { memory->destroy(setflag); + memory->destroy(theta0); } } @@ -48,11 +53,25 @@ void AngleZero::compute(int eflag, int vflag) /* ---------------------------------------------------------------------- */ +void AngleZero::settings(int narg, char **arg) +{ + if ((narg != 0) && (narg != 1)) + error->all(FLERR,"Illegal angle_style command"); + + if (narg == 1) { + if (strcmp("nocoeff",arg[0]) == 0) coeffflag=0; + else error->all(FLERR,"Illegal angle_style command"); + } +} + +/* ---------------------------------------------------------------------- */ + void AngleZero::allocate() { allocated = 1; int n = atom->nangletypes; + memory->create(theta0,n+1,"angle:theta0"); memory->create(setflag,n+1,"angle:setflag"); for (int i = 1; i <= n; i++) setflag[i] = 0; } @@ -63,15 +82,24 @@ void AngleZero::allocate() void AngleZero::coeff(int narg, char **arg) { - if (narg != 1) error->all(FLERR,"Incorrect args for angle coefficients"); + if ((narg < 1) || (coeffflag && narg > 2)) + error->all(FLERR,"Incorrect args for angle coefficients"); + if (!allocated) allocate(); int ilo,ihi; force->bounds(arg[0],atom->nangletypes,ilo,ihi); + double theta0_one = 0.0; + if (coeffflag && (narg == 2)) + theta0_one = force->numeric(FLERR,arg[1]); + + // convert theta0 from degrees to radians + int count = 0; for (int i = ilo; i <= ihi; i++) { setflag[i] = 1; + theta0[i] = theta0_one/180.0 * MY_PI; count++; } @@ -82,14 +110,16 @@ void AngleZero::coeff(int narg, char **arg) double AngleZero::equilibrium_angle(int i) { - return 0.0; + return theta0[i]; } /* ---------------------------------------------------------------------- proc 0 writes out coeffs to restart file ------------------------------------------------------------------------- */ -void AngleZero::write_restart(FILE *fp) {} +void AngleZero::write_restart(FILE *fp) { + fwrite(&theta0[1],sizeof(double),atom->nangletypes,fp); +} /* ---------------------------------------------------------------------- proc 0 reads coeffs from restart file, bcasts them @@ -98,8 +128,23 @@ void AngleZero::write_restart(FILE *fp) {} void AngleZero::read_restart(FILE *fp) { allocate(); + + if (comm->me == 0) { + fread(&theta0[1],sizeof(double),atom->nangletypes,fp); + } + MPI_Bcast(&theta0[1],atom->nangletypes,MPI_DOUBLE,0,world); + for (int i = 1; i <= atom->nangletypes; i++) setflag[i] = 1; } +/* ---------------------------------------------------------------------- + proc 0 writes to data file +------------------------------------------------------------------------- */ + +void AngleZero::write_data(FILE *fp) +{ + for (int i = 1; i <= atom->nangletypes; i++) + fprintf(fp,"%d %g\n",i,theta0[i]/MY_PI*180.0); +} /* ---------------------------------------------------------------------- */ diff --git a/src/angle_zero.h b/src/angle_zero.h index 8d8d966b553cedde14a34416a6787ce3576e6bc2..59c55679a8ed666403a8bdfa33f56c90353feefd 100644 --- a/src/angle_zero.h +++ b/src/angle_zero.h @@ -31,12 +31,19 @@ class AngleZero : public Angle { virtual ~AngleZero(); virtual void compute(int, int); virtual void coeff(int, char **); + virtual void settings(int, char **); + double equilibrium_angle(int); void write_restart(FILE *); void read_restart(FILE *); + void write_data(FILE *); + double single(int, int, int, int); protected: + double *theta0; + int coeffflag; + void allocate(); }; diff --git a/src/bond_zero.cpp b/src/bond_zero.cpp index e8f50deec20516b8a1eee23d12009e68e70299b1..0a354a758eb18e9bf2fe2dc97fe7152c18d31466 100644 --- a/src/bond_zero.cpp +++ b/src/bond_zero.cpp @@ -17,9 +17,11 @@ #include <math.h> #include <stdlib.h> +#include <string.h> #include "bond_zero.h" #include "atom.h" #include "force.h" +#include "comm.h" #include "memory.h" #include "error.h" @@ -27,7 +29,7 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -BondZero::BondZero(LAMMPS *lmp) : Bond(lmp) {} +BondZero::BondZero(LAMMPS *lmp) : Bond(lmp), coeffflag(1) {} /* ---------------------------------------------------------------------- */ @@ -35,6 +37,7 @@ BondZero::~BondZero() { if (allocated && !copymode) { memory->destroy(setflag); + memory->destroy(r0); } } @@ -48,11 +51,25 @@ void BondZero::compute(int eflag, int vflag) /* ---------------------------------------------------------------------- */ +void BondZero::settings(int narg, char **arg) +{ + if ((narg != 0) && (narg != 1)) + error->all(FLERR,"Illegal bond_style command"); + + if (narg == 1) { + if (strcmp("nocoeff",arg[0]) == 0) coeffflag=0; + else error->all(FLERR,"Illegal bond_style command"); + } +} + +/* ---------------------------------------------------------------------- */ + void BondZero::allocate() { allocated = 1; int n = atom->nbondtypes; + memory->create(r0,n+1,"bond:r0"); memory->create(setflag,n+1,"bond:setflag"); for (int i = 1; i <= n; i++) setflag[i] = 0; } @@ -63,15 +80,22 @@ void BondZero::allocate() void BondZero::coeff(int narg, char **arg) { - if (narg != 1) error->all(FLERR,"Incorrect args for bond coefficients"); + if ((narg < 1) || (coeffflag && narg > 2)) + error->all(FLERR,"Incorrect args for bond coefficients"); + if (!allocated) allocate(); int ilo,ihi; force->bounds(arg[0],atom->nbondtypes,ilo,ihi); + double r0_one = 0.0; + if (coeffflag && (narg == 2)) + r0_one = force->numeric(FLERR,arg[1]); + int count = 0; for (int i = ilo; i <= ihi; i++) { setflag[i] = 1; + r0[i] = r0_one; count++; } @@ -84,14 +108,16 @@ void BondZero::coeff(int narg, char **arg) double BondZero::equilibrium_distance(int i) { - return 0.0; + return r0[i]; } /* ---------------------------------------------------------------------- proc 0 writes out coeffs to restart file ------------------------------------------------------------------------- */ -void BondZero::write_restart(FILE *fp) {} +void BondZero::write_restart(FILE *fp) { + fwrite(&r0[1],sizeof(double),atom->nbondtypes,fp); +} /* ---------------------------------------------------------------------- proc 0 reads coeffs from restart file, bcasts them @@ -100,9 +126,27 @@ void BondZero::write_restart(FILE *fp) {} void BondZero::read_restart(FILE *fp) { allocate(); + + if (comm->me == 0) { + fread(&r0[1],sizeof(double),atom->nbondtypes,fp); + } + MPI_Bcast(&r0[1],atom->nbondtypes,MPI_DOUBLE,0,world); + for (int i = 1; i <= atom->nbondtypes; i++) setflag[i] = 1; } +/* ---------------------------------------------------------------------- + proc 0 writes to data file +------------------------------------------------------------------------- */ + +void BondZero::write_data(FILE *fp) +{ + for (int i = 1; i <= atom->nbondtypes; i++) + fprintf(fp,"%d %g\n",i,r0[i]); +} + + + /* ---------------------------------------------------------------------- */ double BondZero::single(int type, double rsq, int i, int j, diff --git a/src/bond_zero.h b/src/bond_zero.h index b9c9a06d98516dfcc0f5614a1b116ed46812b513..a169e094e8ca1129e8a82f59511cdb32f5d47d96 100644 --- a/src/bond_zero.h +++ b/src/bond_zero.h @@ -30,13 +30,20 @@ class BondZero : public Bond { BondZero(class LAMMPS *); virtual ~BondZero(); virtual void compute(int, int); + virtual void settings(int, char **); + void coeff(int, char **); double equilibrium_distance(int); void write_restart(FILE *); void read_restart(FILE *); + void write_data(FILE *); + double single(int, double, int, int, double &); protected: + double *r0; + int coeffflag; + virtual void allocate(); }; diff --git a/src/dihedral_zero.cpp b/src/dihedral_zero.cpp index 3dc9508cf43edc83d7857f9747292f12641e334e..1e4a4ab3ae0b3ca53c0fada325411034e7f231d2 100644 --- a/src/dihedral_zero.cpp +++ b/src/dihedral_zero.cpp @@ -15,12 +15,13 @@ Contributing author: Carsten Svaneborg (SDU) ------------------------------------------------------------------------- */ -#include <mpi.h> #include <math.h> #include <stdlib.h> +#include <string.h> #include "dihedral_zero.h" #include "atom.h" #include "force.h" +#include "comm.h" #include "memory.h" #include "error.h" @@ -28,13 +29,13 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -DihedralZero::DihedralZero(LAMMPS *lmp) : Dihedral(lmp) {} +DihedralZero::DihedralZero(LAMMPS *lmp) : Dihedral(lmp), coeffflag(1) {} /* ---------------------------------------------------------------------- */ DihedralZero::~DihedralZero() { - if (allocated) { + if (allocated && !copymode) { memory->destroy(setflag); } } @@ -49,6 +50,19 @@ void DihedralZero::compute(int eflag, int vflag) /* ---------------------------------------------------------------------- */ +void DihedralZero::settings(int narg, char **arg) +{ + if ((narg != 0) && (narg != 1)) + error->all(FLERR,"Illegal dihedral_style command"); + + if (narg == 1) { + if (strcmp("nocoeff",arg[0]) == 0) coeffflag=0; + else error->all(FLERR,"Illegal dihedral_style command"); + } +} + +/* ---------------------------------------------------------------------- */ + void DihedralZero::allocate() { allocated = 1; @@ -59,12 +73,14 @@ void DihedralZero::allocate() } /* ---------------------------------------------------------------------- - set coeffs for one type + set coeffs for one or more types ------------------------------------------------------------------------- */ void DihedralZero::coeff(int narg, char **arg) { - if (narg != 1) error->all(FLERR,"Incorrect args for dihedral coefficients"); + if ((narg < 1) || (coeffflag && narg > 1)) + error->all(FLERR,"Incorrect args for dihedral coefficients"); + if (!allocated) allocate(); int ilo,ihi; @@ -95,3 +111,12 @@ void DihedralZero::read_restart(FILE *fp) for (int i = 1; i <= atom->ndihedraltypes; i++) setflag[i] = 1; } +/* ---------------------------------------------------------------------- + proc 0 writes to data file +------------------------------------------------------------------------- */ + +void DihedralZero::write_data(FILE *fp) { + for (int i = 1; i <= atom->ndihedraltypes; i++) + fprintf(fp,"%d\n",i); +} + diff --git a/src/dihedral_zero.h b/src/dihedral_zero.h index b075f84989806fe804a8420e81fd17c5797e69f0..8b40a4439fb502ddf0f18ba90b2729d26d61335c 100644 --- a/src/dihedral_zero.h +++ b/src/dihedral_zero.h @@ -9,10 +9,10 @@ the GNU General Public License. See the README file in the top-level LAMMPS directory. - + Identical to dihedral harmonic, except if all k's are zero the force loop is skipped. - + ------------------------------------------------------------------------- */ #ifdef DIHEDRAL_CLASS @@ -24,7 +24,7 @@ DihedralStyle(zero,DihedralZero) #ifndef LMP_DIHEDRAL_ZERO_H #define LMP_DIHEDRAL_ZERO_H -#include "stdio.h" +#include <stdio.h> #include "dihedral.h" namespace LAMMPS_NS { @@ -34,13 +34,17 @@ class DihedralZero : public Dihedral { DihedralZero(class LAMMPS *); virtual ~DihedralZero(); virtual void compute(int, int); - void coeff(int, char **); + virtual void coeff(int, char **); + virtual void settings(int, char **); + void write_restart(FILE *); void read_restart(FILE *); + void write_data(FILE *); protected: + int coeffflag; - void allocate(); + virtual void allocate(); }; } diff --git a/src/improper_zero.cpp b/src/improper_zero.cpp index 438388e10ae4e09b8ec3b233bcb44ccfa705c3fa..958f6621624b198e0e2150d130988bb57aef3f54 100644 --- a/src/improper_zero.cpp +++ b/src/improper_zero.cpp @@ -15,12 +15,13 @@ Contributing author: Carsten Svaneborg (SDU) ------------------------------------------------------------------------- */ -#include <mpi.h> #include <math.h> #include <stdlib.h> +#include <string.h> #include "improper_zero.h" #include "atom.h" #include "force.h" +#include "comm.h" #include "memory.h" #include "error.h" @@ -28,7 +29,7 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -ImproperZero::ImproperZero(LAMMPS *lmp) : Improper(lmp) {} +ImproperZero::ImproperZero(LAMMPS *lmp) : Improper(lmp), coeffflag(1) {} /* ---------------------------------------------------------------------- */ @@ -49,6 +50,19 @@ void ImproperZero::compute(int eflag, int vflag) /* ---------------------------------------------------------------------- */ +void ImproperZero::settings(int narg, char **arg) +{ + if ((narg != 0) && (narg != 1)) + error->all(FLERR,"Illegal improper_style command"); + + if (narg == 1) { + if (strcmp("nocoeff",arg[0]) == 0) coeffflag=0; + else error->all(FLERR,"Illegal improper_style command"); + } +} + +/* ---------------------------------------------------------------------- */ + void ImproperZero::allocate() { allocated = 1; @@ -59,12 +73,14 @@ void ImproperZero::allocate() } /* ---------------------------------------------------------------------- - set coeffs for one type + set coeffs for one or more types ------------------------------------------------------------------------- */ void ImproperZero::coeff(int narg, char **arg) { - if (narg != 1) error->all(FLERR,"Incorrect args for improper coefficients"); + if ((narg < 1) || (coeffflag && narg > 1)) + error->all(FLERR,"Incorrect args for improper coefficients"); + if (!allocated) allocate(); int ilo,ihi; @@ -94,3 +110,13 @@ void ImproperZero::read_restart(FILE *fp) allocate(); for (int i = 1; i <= atom->nimpropertypes; i++) setflag[i] = 1; } + +/* ---------------------------------------------------------------------- + proc 0 writes to data file +------------------------------------------------------------------------- */ + +void ImproperZero::write_data(FILE *fp) { + for (int i = 1; i <= atom->nimpropertypes; i++) + fprintf(fp,"%d\n",i); +} + diff --git a/src/improper_zero.h b/src/improper_zero.h index 7283e858fbc6cd8728a43888aa535bc23d073f38..d999bdd8d9d87343f884583149a3a2e14f13321a 100644 --- a/src/improper_zero.h +++ b/src/improper_zero.h @@ -31,10 +31,15 @@ class ImproperZero : public Improper { virtual ~ImproperZero(); virtual void compute(int, int); virtual void coeff(int, char **); + virtual void settings(int, char **); + void write_restart(FILE *); void read_restart(FILE *); + void write_data(FILE *); protected: + int coeffflag; + virtual void allocate(); }; @@ -45,11 +50,6 @@ class ImproperZero : public Improper { /* ERROR/WARNING messages: -W: Improper problem: %d %ld %d %d %d %d - -Conformation of the 4 listed improper atoms is extreme; you may want -to check your simulation geometry. - E: Incorrect args for improper coefficients Self-explanatory. Check the input script or data file. diff --git a/src/pair_zero.cpp b/src/pair_zero.cpp index c3af3f0efe42a49161d16b189b8b59ee6bff4221..e330ee0361da2189208d3b158918a7aa1a805913 100644 --- a/src/pair_zero.cpp +++ b/src/pair_zero.cpp @@ -30,7 +30,7 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -PairZero::PairZero(LAMMPS *lmp) : Pair(lmp) {} +PairZero::PairZero(LAMMPS *lmp) : Pair(lmp), coeffflag(1) {} /* ---------------------------------------------------------------------- */ @@ -49,7 +49,7 @@ void PairZero::compute(int eflag, int vflag) { if (eflag || vflag) ev_setup(eflag,vflag); else evflag = vflag_fdotr = 0; - + if (vflag_fdotr) virial_fdotr_compute(); } @@ -77,13 +77,17 @@ void PairZero::allocate() void PairZero::settings(int narg, char **arg) { - if (narg != 1) error->all(FLERR,"Illegal pair_style command"); + if ((narg != 1) && (narg != 2)) + error->all(FLERR,"Illegal pair_style command"); cut_global = force->numeric(FLERR,arg[0]); + if (narg == 2) { + if (strcmp("nocoeff",arg[1]) == 0) coeffflag=0; + else error->all(FLERR,"Illegal pair_style command"); + } // reset cutoffs that have been explicitly set - allocate(); int i,j; for (i = 1; i <= atom->ntypes; i++) for (j = i+1; j <= atom->ntypes; j++) @@ -96,8 +100,9 @@ void PairZero::settings(int narg, char **arg) void PairZero::coeff(int narg, char **arg) { - if (narg < 2 || narg > 3) + if ((narg < 2) || (coeffflag && narg > 3)) error->all(FLERR,"Incorrect args for pair coefficients"); + if (!allocated) allocate(); int ilo,ihi,jlo,jhi; @@ -105,7 +110,7 @@ void PairZero::coeff(int narg, char **arg) force->bounds(arg[1],atom->ntypes,jlo,jhi); double cut_one = cut_global; - if (narg == 3) cut_one = force->numeric(FLERR,arg[2]); + if (coeffflag && (narg == 3)) cut_one = force->numeric(FLERR,arg[2]); int count = 0; for (int i = ilo; i <= ihi; i++) { @@ -180,6 +185,7 @@ void PairZero::read_restart(FILE *fp) void PairZero::write_restart_settings(FILE *fp) { fwrite(&cut_global,sizeof(double),1,fp); + fwrite(&coeffflag,sizeof(int),1,fp); } /* ---------------------------------------------------------------------- @@ -191,7 +197,9 @@ void PairZero::read_restart_settings(FILE *fp) int me = comm->me; if (me == 0) { fread(&cut_global,sizeof(double),1,fp); + fread(&coeffflag,sizeof(int),1,fp); } MPI_Bcast(&cut_global,1,MPI_DOUBLE,0,world); + MPI_Bcast(&coeffflag,1,MPI_INT,0,world); } diff --git a/src/pair_zero.h b/src/pair_zero.h index e6cfb32579f996fc4539660fad0f4c8e852d3340..b0994deef866f47aa6039b8e35414a70303ba893 100644 --- a/src/pair_zero.h +++ b/src/pair_zero.h @@ -9,16 +9,16 @@ the GNU General Public License. See the README file in the top-level LAMMPS directory. - - Pair zero is a dummy pair interaction useful for requiring a + + Pair zero is a dummy pair interaction useful for requiring a force cutoff distance in the absense of pair-interactions or with hybrid/overlay if a larger force cutoff distance is required. - + This can be used in conjunction with bond/create to create bonds that are longer than the cutoff of a given force field, or to calculate radial distribution functions for models without pair interactions. - + ------------------------------------------------------------------------- */ #ifdef PAIR_CLASS @@ -50,6 +50,7 @@ class PairZero : public Pair { protected: double cut_global; double **cut; + int coeffflag; virtual void allocate(); }; diff --git a/src/read_data.cpp b/src/read_data.cpp index 4b790f343cd4c9ae8456b6135ec132680ceea202..39afcfbcf8944e82cb29d95c374aafb52245ea5f 100644 --- a/src/read_data.cpp +++ b/src/read_data.cpp @@ -123,6 +123,7 @@ void ReadData::command(int narg, char **arg) // optional args addflag = NONE; + coeffflag = 1; id_offset = 0; offsetflag = shiftflag = 0; toffset = boffset = aoffset = doffset = ioffset = 0; @@ -172,7 +173,9 @@ void ReadData::command(int narg, char **arg) if (domain->dimension == 2 && shift[2] != 0.0) error->all(FLERR,"Non-zero read_data shift z value for 2d simulation"); iarg += 4; - + } else if (strcmp(arg[iarg],"nocoeff") == 0) { + coeffflag = 0; + iarg ++; } else if (strcmp(arg[iarg],"extra/atom/types") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal read_data command"); extra_atom_types = force->inumeric(FLERR,arg[iarg+1]); @@ -279,6 +282,49 @@ void ReadData::command(int narg, char **arg) MPI_Allreduce(&max,&id_offset,1,MPI_LMP_TAGINT,MPI_MAX,world); } + // set up pointer to hold original styles while we replace them with "zero" + Pair *saved_pair = NULL; + Bond *saved_bond = NULL; + Angle *saved_angle = NULL; + Dihedral *saved_dihedral = NULL; + Improper *saved_improper = NULL; + KSpace *saved_kspace = NULL; + + if (coeffflag == 0) { + char *coeffs[2]; + coeffs[0] = (char *) "10.0"; + coeffs[1] = (char *) "nocoeff"; + + saved_pair = force->pair; + force->pair = NULL; + force->create_pair("zero",0); + if (force->pair) force->pair->settings(2,coeffs); + + coeffs[0] = coeffs[1]; + saved_bond = force->bond; + force->bond = NULL; + force->create_bond("zero",0); + if (force->bond) force->bond->settings(1,coeffs); + + saved_angle = force->angle; + force->angle = NULL; + force->create_angle("zero",0); + if (force->angle) force->angle->settings(1,coeffs); + + saved_dihedral = force->dihedral; + force->dihedral = NULL; + force->create_dihedral("zero",0); + if (force->dihedral) force->dihedral->settings(1,coeffs); + + saved_improper = force->improper; + force->improper = NULL; + force->create_improper("zero",0); + if (force->improper) force->improper->settings(1,coeffs); + + saved_kspace = force->kspace; + force->kspace = NULL; + } + // ----------------------------------------------------------------- // perform 1-pass read if no molecular topology in file @@ -778,6 +824,26 @@ void ReadData::command(int narg, char **arg) error->all(FLERR, "Read_data shrink wrap did not assign all atoms correctly"); } + + // restore old styles, when reading with nocoeff flag given + if (coeffflag == 0) { + if (force->pair) delete force->pair; + force->pair = saved_pair; + + if (force->bond) delete force->bond; + force->bond = saved_bond; + + if (force->angle) delete force->angle; + force->angle = saved_angle; + + if (force->dihedral) delete force->dihedral; + force->dihedral = saved_dihedral; + + if (force->improper) delete force->improper; + force->improper = saved_improper; + + force->kspace = saved_kspace; + } } /* ---------------------------------------------------------------------- diff --git a/src/read_data.h b/src/read_data.h index 9402bd18e4c63914b87102c989e70a6a4cd3c698..5463c86f08e2fa849c9254f13f1feec008e203d9 100644 --- a/src/read_data.h +++ b/src/read_data.h @@ -64,7 +64,7 @@ class ReadData : protected Pointers { // optional args - int addflag,offsetflag,shiftflag; + int addflag,offsetflag,shiftflag,coeffflag; tagint addvalue; int toffset,boffset,aoffset,doffset,ioffset; double shift[3];