diff --git a/src/GRANULAR/fix_pour.cpp b/src/GRANULAR/fix_pour.cpp index 1f55e7bc63..c4003f5301 100644 --- a/src/GRANULAR/fix_pour.cpp +++ b/src/GRANULAR/fix_pour.cpp @@ -23,8 +23,6 @@ #include "molecule.h" #include "modify.h" #include "fix_gravity.h" -#include "fix_rigid_small.h" -#include "fix_shake.h" #include "domain.h" #include "region.h" #include "region_block.h" @@ -325,8 +323,9 @@ void FixPour::init() if (rigidflag) { int ifix = modify->find_fix(idrigid); if (ifix < 0) error->all(FLERR,"Fix pour rigid fix does not exist"); - fixrigid = (FixRigidSmall *) modify->fix[ifix]; - if (onemol != fixrigid->onemol) + fixrigid = modify->fix[ifix]; + int tmp; + if (onemol != (Molecule *) fixrigid->extract("onemol",tmp)) error->all(FLERR, "Fix pour and fix rigid/small not using same molecule ID"); } @@ -338,10 +337,10 @@ void FixPour::init() if (shakeflag) { int ifix = modify->find_fix(idshake); if (ifix < 0) error->all(FLERR,"Fix pour shake fix does not exist"); - fixshake = (FixShake *) modify->fix[ifix]; - if (onemol != fixshake->onemol) - error->all(FLERR, - "Fix pour and fix shake not using same molecule ID"); + fixshake = modify->fix[ifix]; + int tmp; + if (onemol != (Molecule *) fixshake->extract("onemol",tmp)) + error->all(FLERR,"Fix pour and fix shake not using same molecule ID"); } } @@ -493,8 +492,8 @@ void FixPour::pre_exchange() coords[i][2] += coord[2]; // coords[3] = particle radius - // default to 0.5, if not defined in Molecule - // same as atom->avec->create_atom() when invoked below + // default to 0.5, if radii not defined in Molecule + // same as atom->avec->create_atom(), invoked below if (onemol->radiusflag) coords[i][3] = onemol->radius[i]; else coords[i][3] = 0.5; @@ -942,3 +941,34 @@ void FixPour::reset_dt() error->all(FLERR,"Cannot change timestep with fix pour"); } +/* ---------------------------------------------------------------------- + extract particle radius for atom type = itype +------------------------------------------------------------------------- */ + +void *FixPour::extract(const char *str, int &itype) +{ + if (strcmp(str,"radius") == 0) { + if (mode == ATOM) { + if (itype == ntype) oneradius = radius_max; + else oneradius = 0.0; + } else { + double *radius = onemol->radius; + int *type = onemol->type; + int natoms = onemol->natoms; + + // check radii of matching types in Molecule + // default to 0.5, if radii not defined in Molecule + // same as atom->avec->create_atom(), invoked in pre_exchange() + + oneradius = 0.0; + for (int i = 0; i < natoms; i++) + if (type[i] == itype-ntype) { + if (radius) oneradius = MAX(oneradius,radius[i]); + else oneradius = 0.5; + } + } + itype = 0; + return &oneradius; + } + return NULL; +} diff --git a/src/GRANULAR/fix_pour.h b/src/GRANULAR/fix_pour.h index 8c11978d28..5f2b853976 100644 --- a/src/GRANULAR/fix_pour.h +++ b/src/GRANULAR/fix_pour.h @@ -25,15 +25,6 @@ FixStyle(pour,FixPour) namespace LAMMPS_NS { class FixPour : public Fix { - friend class PairGranHertzHistory; - friend class PairGranHertzHistoryOMP; - friend class PairGranHooke; - friend class PairGranHookeOMP; - friend class PairGranHookeHistory; - friend class PairGranHookeHistory2; - friend class PairGranHookeHistoryOMP; - friend class PairGranHookeCuda; - public: FixPour(class LAMMPS *, int, char **); ~FixPour(); @@ -41,6 +32,7 @@ class FixPour : public Fix { void init(); void pre_exchange(); void reset_dt(); + void *extract(const char *, int &); private: int ninsert,ntype,seed; @@ -63,8 +55,8 @@ class FixPour : public Fix { int natom; double **coords; int *imageflags; - class FixRigidSmall *fixrigid; - class FixShake *fixshake; + class Fix *fixrigid,*fixshake; + double oneradius; int me,nprocs; int *recvcounts,*displs; diff --git a/src/GRANULAR/pair_gran_hooke_history.cpp b/src/GRANULAR/pair_gran_hooke_history.cpp index 638d89222c..ccebe88f0e 100644 --- a/src/GRANULAR/pair_gran_hooke_history.cpp +++ b/src/GRANULAR/pair_gran_hooke_history.cpp @@ -27,8 +27,6 @@ #include "update.h" #include "modify.h" #include "fix.h" -#include "fix_pour.h" -#include "fix_deposit.h" #include "fix_shear_history.h" #include "comm.h" #include "neighbor.h" @@ -451,44 +449,39 @@ void PairGranHookeHistory::init_style() if (i < modify->nfix) freeze_group_bit = modify->fix[i]->groupbit; else freeze_group_bit = 0; - // check for FixPour and set pour_type and pour_maxrad - - int pour_type = 0; - double pour_maxrad = 0.0; - for (i = 0; i < modify->nfix; i++) - if (strcmp(modify->fix[i]->style,"pour") == 0) break; - if (i < modify->nfix) { - pour_type = ((FixPour *) modify->fix[i])->ntype; - pour_maxrad = ((FixPour *) modify->fix[i])->radius_max; - } - - // check for FixDeposit and set deposit_type and deposit_maxrad - - int deposit_type = 0; - double deposit_maxrad = 0.0; - for (i = 0; i < modify->nfix; i++) - if (strcmp(modify->fix[i]->style,"deposit") == 0) break; - if (i < modify->nfix) { - deposit_type = ((FixDeposit *) modify->fix[i])->ntype; - deposit_maxrad = 0.5; - } - - // check for FixRigid + // check for FixRigid so can extract rigid body masses fix_rigid = NULL; for (i = 0; i < modify->nfix; i++) if (modify->fix[i]->rigid_flag) break; if (i < modify->nfix) fix_rigid = modify->fix[i]; + // check for FixPour and FixDeposit so can extract particle radii + + int ipour; + for (ipour = 0; ipour < modify->nfix; ipour++) + if (strcmp(modify->fix[ipour]->style,"pour") == 0) break; + if (ipour == modify->nfix) ipour = -1; + + int idep; + for (idep = 0; idep < modify->nfix; idep++) + if (strcmp(modify->fix[idep]->style,"deposit") == 0) break; + if (idep == modify->nfix) idep = -1; + // set maxrad_dynamic and maxrad_frozen for each type // include future FixPour and FixDeposit particles as dynamic - for (i = 1; i <= atom->ntypes; i++) + int itype; + for (i = 1; i <= atom->ntypes; i++) { onerad_dynamic[i] = onerad_frozen[i] = 0.0; - if (pour_type) onerad_dynamic[pour_type] = pour_maxrad; - if (deposit_type) onerad_dynamic[deposit_type] = - MAX(onerad_dynamic[deposit_type],deposit_maxrad); - + if (ipour >= 0) + onerad_dynamic[i] = + *((double *) modify->fix[ipour]->extract("radius",itype)); + if (idep >= 0) + onerad_dynamic[i] = + *((double *) modify->fix[idep]->extract("radius",itype)); + } + double *radius = atom->radius; int *mask = atom->mask; int *type = atom->type; diff --git a/src/RIGID/fix_rigid_small.cpp b/src/RIGID/fix_rigid_small.cpp index c00dd22c57..0fe1b607a8 100644 --- a/src/RIGID/fix_rigid_small.cpp +++ b/src/RIGID/fix_rigid_small.cpp @@ -3019,6 +3019,11 @@ void *FixRigidSmall::extract(const char *str, int &dim) return atom2body; } + if (strcmp(str,"onemol") == 0) { + dim = 0; + return onemol; + } + // return vector of rigid body masses, for owned+ghost bodies // used by granular pair styles, indexed by atom2body diff --git a/src/RIGID/fix_rigid_small.h b/src/RIGID/fix_rigid_small.h index fc00da65b0..dcf3ff6af7 100644 --- a/src/RIGID/fix_rigid_small.h +++ b/src/RIGID/fix_rigid_small.h @@ -33,10 +33,6 @@ class FixRigidSmall : public Fix { static FixRigidSmall *frsptr; - // molecules being added on-the-fly as rigid bodies - - class Molecule *onemol; - FixRigidSmall(class LAMMPS *, int, char **); virtual ~FixRigidSmall(); virtual int setmask(); @@ -52,6 +48,8 @@ class FixRigidSmall : public Fix { void grow_arrays(int); void copy_arrays(int, int, int); void set_arrays(int); + void set_molecule(int, int, double *, double *, double *); + int pack_exchange(int, double *); int unpack_exchange(int, double *); int pack_comm(int, int *, double *, int, int *); @@ -72,8 +70,6 @@ class FixRigidSmall : public Fix { double compute_scalar(); double memory_usage(); - void set_molecule(int, int, double *, double *, double *); - protected: int me,nprocs; double dtv,dtf,dtq; @@ -155,6 +151,10 @@ class FixRigidSmall : public Fix { int maxlang; // max size of langextra class RanMars *random; // RNG + // molecules added on-the-fly as rigid bodies + + class Molecule *onemol; + // class data used by ring communication callbacks std::map *hash;