git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@11155 f3b2605a-c512-4ea7-a41b-209d697bcdaa

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;
+}

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;

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;

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
 

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<int,int> *hash;