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

src/GRANULAR/fix_pour.cpp

@@ -463,9 +463,6 @@ void FixPour::pre_exchange()
   double radtmp,delx,dely,delz,rsq,radsum,rn,h;
   double coord[3];
 
-  int nfix = modify->nfix;
-  Fix **fix = modify->fix;
-
   double denstmp;
   double *sublo = domain->sublo;
   double *subhi = domain->subhi;
@@ -642,8 +639,7 @@ void FixPour::pre_exchange()
           atom->radius[n] = radtmp;
           atom->rmass[n] = 4.0*MY_PI/3.0 * radtmp*radtmp*radtmp * denstmp;
         } else atom->add_molecule_atom(onemols[imol],m,n,maxtag_all);
-        for (j = 0; j < nfix; j++)
-          if (fix[j]->create_attribute) fix[j]->set_arrays(n);
+        modify->create_attribute(n);
       }
     }
 

src/MC/fix_gcmc.cpp

@@ -755,11 +755,7 @@ void FixGCMC::attempt_atomic_insertion()
       atom->v[m][0] = random_unequal->gaussian()*sigma;
       atom->v[m][1] = random_unequal->gaussian()*sigma;
       atom->v[m][2] = random_unequal->gaussian()*sigma;
-
-      int nfix = modify->nfix;
-      Fix **fix = modify->fix;
-      for (int j = 0; j < nfix; j++)
-        if (fix[j]->create_attribute) fix[j]->set_arrays(m);
+      modify->create_attribute(m);
 
       success = 1;
     }
@@ -1105,9 +1101,7 @@ void FixGCMC::attempt_molecule_insertion()
         atom->v[m][2] = vnew[2];
         
         atom->add_molecule_atom(onemols[imol],i,m,maxtag_all);
-        for (int j = 0; j < nfix; j++)
-          if (fix[j]->create_attribute) fix[j]->set_arrays(m);
-
+        modify->create_attribute(m);
       }
     }
 
@@ -1341,11 +1335,7 @@ void FixGCMC::attempt_atomic_insertion_full()
     atom->v[m][1] = random_unequal->gaussian()*sigma;
     atom->v[m][2] = random_unequal->gaussian()*sigma;
     if (charge_flag) atom->q[m] = charge;
-    
-    int nfix = modify->nfix;
-    Fix **fix = modify->fix;
-    for (int j = 0; j < nfix; j++)
-      if (fix[j]->create_attribute) fix[j]->set_arrays(m);
+    modify->create_attribute(m);
   }
 
   atom->natoms++;
@@ -1698,9 +1688,7 @@ void FixGCMC::attempt_molecule_insertion_full()
       atom->v[m][2] = vnew[2];
 
       atom->add_molecule_atom(onemols[imol],i,m,maxtag_all);
-
-      for (int j = 0; j < nfix; j++)
-        if (fix[j]->create_attribute) fix[j]->set_arrays(m);
+      modify->create_attribute(m);
     }
   }
 

src/MISC/fix_deposit.cpp

@@ -287,8 +287,6 @@ void FixDeposit::pre_exchange()
   // attempt an insertion until successful
 
   int dimension = domain->dimension;
-  int nfix = modify->nfix;
-  Fix **fix = modify->fix;
 
   int success = 0;
   int attempt = 0;
@@ -492,8 +490,7 @@ void FixDeposit::pre_exchange()
         atom->v[n][2] = vnew[2];
         if (mode == MOLECULE) 
           atom->add_molecule_atom(onemols[imol],m,n,maxtag_all);
-        for (j = 0; j < nfix; j++)
-          if (fix[j]->create_attribute) fix[j]->set_arrays(n);
+        modify->create_attribute(n);
       }
     }
 

src/USER-FEP/fix_adapt_fep.cpp

@@ -52,6 +52,7 @@ FixAdaptFEP::FixAdaptFEP(LAMMPS *lmp, int narg, char **arg) :
   if (nevery < 0) error->all(FLERR,"Illegal fix adapt/fep command");
 
   dynamic_group_allow = 1;
+  create_attribute = 1;
 
   // count # of adaptations
 
@@ -583,3 +584,13 @@ void FixAdaptFEP::restore_settings()
   if (anypair) force->pair->reinit();
   if (chgflag && force->kspace) force->kspace->qsum_qsq();
 }
+
+/* ----------------------------------------------------------------------
+   initialize one atom's storage values, called when atom is created
+------------------------------------------------------------------------- */
+
+void FixAdaptFEP::set_arrays(int i)
+{
+  if (fix_diam) fix_diam->vstore[i] = atom->radius[i];
+  if (fix_chg) fix_chg->vstore[i] = atom->q[i];
+}

src/USER-FEP/fix_adapt_fep.h

@@ -39,6 +39,7 @@ class FixAdaptFEP : public Fix {
   void post_run();
   void setup_pre_force_respa(int,int);
   void pre_force_respa(int,int,int);
+  void set_arrays(int);
 
  private:
   int nadapt,resetflag,scaleflag,afterflag;

src/compute.cpp

@@ -71,6 +71,7 @@ Compute::Compute(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp)
 
   tempflag = pressflag = peflag = 0;
   pressatomflag = peatomflag = 0;
+  create_attribute = 0;
   tempbias = 0;
 
   timeflag = 0;

src/compute.h

@@ -59,6 +59,8 @@ class Compute : protected Pointers {
   int pressatomflag;  // 1 if Compute calculates per-atom virial
   int peflag;         // 1 if Compute calculates PE (uses Force energies)
   int peatomflag;     // 1 if Compute calculates per-atom PE
+  int create_attribute;    // 1 if compute stores attributes that need
+                           // setting when a new atom is created
 
   int tempbias;       // 0/1 if Compute temp includes self/extra bias
 
@@ -99,6 +101,7 @@ class Compute : protected Pointers {
   virtual void compute_array() {}
   virtual void compute_peratom() {}
   virtual void compute_local() {}
+  virtual void set_arrays(int) {}
 
   virtual int pack_forward_comm(int, int *, double *, int, int *) {return 0;}
   virtual void unpack_forward_comm(int, int, double *) {}

src/compute_displace_atom.cpp

@@ -35,6 +35,7 @@ ComputeDisplaceAtom::ComputeDisplaceAtom(LAMMPS *lmp, int narg, char **arg) :
 
   peratom_flag = 1;
   size_peratom_cols = 4;
+  create_attribute = 1;
 
   // create a new fix STORE style
   // id = compute-ID + COMPUTE_STORE, fix group = compute group
@@ -168,6 +169,19 @@ void ComputeDisplaceAtom::compute_peratom()
   }
 }
 
+/* ----------------------------------------------------------------------
+   initialize one atom's storage values, called when atom is created
+------------------------------------------------------------------------- */
+
+void ComputeDisplaceAtom::set_arrays(int i)
+{
+  double **xoriginal = fix->astore;
+  double **x = atom->x;
+  xoriginal[i][0] = x[i][0];
+  xoriginal[i][1] = x[i][1];
+  xoriginal[i][2] = x[i][2];
+}
+
 /* ----------------------------------------------------------------------
    memory usage of local atom-based array
 ------------------------------------------------------------------------- */

src/compute_displace_atom.h

@@ -30,6 +30,7 @@ class ComputeDisplaceAtom : public Compute {
   ~ComputeDisplaceAtom();
   void init();
   void compute_peratom();
+  void set_arrays(int);
   double memory_usage();
 
  private:

src/compute_msd.cpp

@@ -33,6 +33,7 @@ ComputeMSD::ComputeMSD(LAMMPS *lmp, int narg, char **arg) :
   vector_flag = 1;
   size_vector = 4;
   extvector = 0;
+  create_attribute = 1;
 
   // optional args
 
@@ -206,3 +207,16 @@ void ComputeMSD::compute_vector()
     vector[3] /= nmsd;
   }
 }
+
+/* ----------------------------------------------------------------------
+   initialize one atom's storage values, called when atom is created
+------------------------------------------------------------------------- */
+
+void ComputeMSD::set_arrays(int i)
+{
+  double **xoriginal = fix->astore;
+  double **x = atom->x;
+  xoriginal[i][0] = x[i][0];
+  xoriginal[i][1] = x[i][1];
+  xoriginal[i][2] = x[i][2];
+}

src/compute_msd.h

@@ -30,6 +30,7 @@ class ComputeMSD : public Compute {
   virtual ~ComputeMSD();
   void init();
   virtual void compute_vector();
+  void set_arrays(int);
 
  protected:
   int comflag;

src/compute_vacf.cpp

@@ -32,6 +32,7 @@ ComputeVACF::ComputeVACF(LAMMPS *lmp, int narg, char **arg) :
   vector_flag = 1;
   size_vector = 4;
   extvector = 0;
+  create_attribute = 1;
 
   // create a new fix STORE style
   // id = compute-ID + COMPUTE_STORE, fix group = compute group
@@ -137,3 +138,17 @@ void ComputeVACF::compute_vector()
     vector[3] /= nvacf;
   }
 }
+
+/* ----------------------------------------------------------------------
+   initialize one atom's storage values, called when atom is created
+------------------------------------------------------------------------- */
+
+void ComputeVACF::set_arrays(int i)
+{
+  double **voriginal = fix->astore;
+  double **v = atom->v;
+  voriginal[i][0] = v[i][0];
+  voriginal[i][1] = v[i][1];
+  voriginal[i][2] = v[i][2];
+}
+

src/compute_vacf.h

@@ -30,6 +30,7 @@ class ComputeVACF : public Compute {
   ~ComputeVACF();
   void init();
   virtual void compute_vector();
+  void set_arrays(int);
 
  protected:
   bigint nvacf;

src/create_atoms.cpp

@@ -24,6 +24,7 @@
 #include "force.h"
 #include "special.h"
 #include "fix.h"
+#include "compute.h"
 #include "domain.h"
 #include "lattice.h"
 #include "region.h"
@@ -355,7 +356,10 @@ void CreateAtoms::command(int narg, char **arg)
   else if (style == RANDOM) add_random();
   else add_lattice();
 
-  // invoke set_arrays() for fixes that need initialization of new atoms
+  // invoke set_arrays() for fixes/computes/variables
+  //   that need initialization of attributes of new atoms
+  // don't use modify->create_attributes() since would be inefficient
+  //   for large number of atoms
 
   int nlocal = atom->nlocal;
   for (int m = 0; m < modify->nfix; m++) {
@@ -364,6 +368,14 @@ void CreateAtoms::command(int narg, char **arg)
       for (int i = nlocal_previous; i < nlocal; i++)
         fix->set_arrays(i);
   }
+  for (int m = 0; m < modify->ncompute; m++) {
+    Compute *compute = modify->compute[m];
+    if (compute->create_attribute)
+      for (int i = nlocal_previous; i < nlocal; i++)
+        compute->set_arrays(i);
+  }
+  for (int i = nlocal_previous; i < nlocal; i++)
+    input->variable->set_arrays(i);
 
   // restore each equal variable string previously saved
 

src/fix_adapt.cpp

@@ -47,6 +47,7 @@ FixAdapt::FixAdapt(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
   if (nevery < 0) error->all(FLERR,"Illegal fix adapt command");
 
   dynamic_group_allow = 1;
+  create_attribute = 1;
 
   // count # of adaptations
 
@@ -554,3 +555,13 @@ void FixAdapt::restore_settings()
   if (anypair) force->pair->reinit();
   if (chgflag && force->kspace) force->kspace->qsum_qsq();
 }
+
+/* ----------------------------------------------------------------------
+   initialize one atom's storage values, called when atom is created
+------------------------------------------------------------------------- */
+
+void FixAdapt::set_arrays(int i)
+{
+  if (fix_diam) fix_diam->vstore[i] = atom->radius[i];
+  if (fix_chg) fix_chg->vstore[i] = atom->q[i];
+}

src/fix_adapt.h

@@ -39,6 +39,7 @@ class FixAdapt : public Fix {
   void post_run();
   void setup_pre_force_respa(int,int);
   void pre_force_respa(int,int,int);
+  void set_arrays(int);
 
  private:
   int nadapt,resetflag,scaleflag;

src/fix_store.cpp

@@ -28,8 +28,6 @@ FixStore::FixStore(LAMMPS *lmp, int narg, char **arg) : Fix(lmp, narg, arg)
 {
   if (narg != 5) error->all(FLERR,"Illegal fix store command");
 
-  create_attribute = 1;
-
   // syntax: id group style 0/1 nvalue
   // 0/1 flag = not-store or store values in restart file
 
@@ -191,17 +189,3 @@ int FixStore::size_restart(int nlocal)
 {
   return nvalues+1;
 }
-
-/* ----------------------------------------------------------------------
-   initialize one atom's vector/array values, called when atom is created
-   value of 0.0 may not be desired value, but will be unitialized otherwise
-     e.g. by compute msd on group all
-------------------------------------------------------------------------- */
-
-void FixStore::set_arrays(int i)
-{
-  if (vecflag) 
-    vstore[i] = 0.0;
-  else
-    for (int m = 0; m < nvalues; m++) astore[i][m] = 0.0;
-}

src/fix_store.h

@@ -42,7 +42,6 @@ class FixStore : public Fix {
   void unpack_restart(int, int);
   int size_restart(int);
   int maxsize_restart();
-  void set_arrays(int);
 
  private:
   int nvalues;                  // total # of values per atom

src/modify.cpp

@@ -23,6 +23,8 @@
 #include "group.h"
 #include "update.h"
 #include "domain.h"
+#include "input.h"
+#include "variable.h"
 #include "memory.h"
 #include "error.h"
 
@@ -436,6 +438,24 @@ void Modify::post_run()
   for (int i = 0; i < nfix; i++) fix[i]->post_run();
 }
 
+/* ----------------------------------------------------------------------
+   create_attribute call
+   invoked when an atom is added to system during a run
+   necessary so that fixes and computes that store per-atom
+     state can initialize that state for the new atom N
+   computes can store per-atom state via a fix like fix STORE
+     compute has the create_attribute flag, not fix STORE
+------------------------------------------------------------------------- */
+
+void Modify::create_attribute(int n)
+{
+  for (int i = 0; i < nfix; i++)
+    if (fix[i]->create_attribute) fix[i]->set_arrays(n);
+  for (int i = 0; i < ncompute; i++)
+    if (compute[i]->create_attribute) compute[i]->set_arrays(n);
+  input->variable->set_arrays(n);
+}
+
 /* ----------------------------------------------------------------------
    setup rRESPA pre_force call, only for relevant fixes
 ------------------------------------------------------------------------- */

src/modify.h

@@ -59,6 +59,7 @@ class Modify : protected Pointers {
   virtual void end_of_step();
   virtual double thermo_energy();
   virtual void post_run();
+  virtual void create_attribute(int);
 
   virtual void setup_pre_force_respa(int, int);
   virtual void initial_integrate_respa(int, int, int);

src/read_dump.cpp

@@ -31,10 +31,13 @@
 #include "update.h"
 #include "modify.h"
 #include "fix.h"
+#include "compute.h"
 #include "domain.h"
 #include "comm.h"
 #include "force.h"
 #include "irregular.h"
+#include "input.h"
+#include "variable.h"
 #include "error.h"
 #include "memory.h"
 
@@ -902,16 +905,27 @@ void ReadDump::process_atoms(int n)
     }
   }
 
-  // invoke set_arrays() for fixes that need initialization of new atoms
+  // invoke set_arrays() for fixes/computes/variables
+  //   that need initialization of attributes of new atoms
   // same as in CreateAtoms
+  // don't use modify->create_attributes() since would be inefficient
+  //   for large number of atoms
 
   nlocal = atom->nlocal;
-  for (m = 0; m < modify->nfix; m++) {
+  for (int m = 0; m < modify->nfix; m++) {
     Fix *fix = modify->fix[m];
     if (fix->create_attribute)
       for (i = nlocal_previous; i < nlocal; i++)
         fix->set_arrays(i);
   }
+  for (int m = 0; m < modify->ncompute; m++) {
+    Compute *compute = modify->compute[m];
+    if (compute->create_attribute)
+      for (i = nlocal_previous; i < nlocal; i++)
+        compute->set_arrays(i);
+  }
+  for (int i = nlocal_previous; i < nlocal; i++)
+    input->variable->set_arrays(i);
 }
 
 /* ----------------------------------------------------------------------

src/variable.cpp

@@ -696,7 +696,7 @@ void Variable::compute_atom(int ivar, int igroup,
   if (style[ivar] == ATOM) {
     evaluate(data[ivar][0],&tree);
     collapse_tree(tree);
-  } else vstore = reader[ivar]->fix->vstore;
+  } else vstore = reader[ivar]->fixstore->vstore;
 
   int groupbit = group->bitmask[igroup];
   int *mask = atom->mask;
@@ -754,6 +754,15 @@ int Variable::find(char *name)
   return -1;
 }
 
+/* ----------------------------------------------------------------------
+   initialize one atom's storage values in all VarReaders via fix STORE
+   called when atom is created
+------------------------------------------------------------------------- */
+
+void Variable::set_arrays(int i)
+{
+}
+
 /* ----------------------------------------------------------------------
    return 1 if variable is EQUAL style, 0 if not
 ------------------------------------------------------------------------- */
@@ -1471,7 +1480,7 @@ double Variable::evaluate(char *str, Tree **tree)
                        "equal-style variable formula");
           Tree *newtree = new Tree();
           newtree->type = ATOMARRAY;
-          newtree->array = reader[ivar]->fix->vstore;
+          newtree->array = reader[ivar]->fixstore->vstore;
           newtree->nstride = 1;
           newtree->selfalloc = 0;
           newtree->first = newtree->second = NULL;
@@ -1495,7 +1504,7 @@ double Variable::evaluate(char *str, Tree **tree)
 
         } else if (nbracket && style[ivar] == ATOMFILE) {
 
-          peratom2global(1,NULL,reader[ivar]->fix->vstore,1,index,
+          peratom2global(1,NULL,reader[ivar]->fixstore->vstore,1,index,
                          tree,treestack,ntreestack,argstack,nargstack);
 
         } else error->all(FLERR,"Mismatched variable in variable formula");
@@ -3641,7 +3650,7 @@ int Variable::special_function(char *word, char *contents, Tree **tree,
 
       double *result;
       memory->create(result,atom->nlocal,"variable:result");
-      memcpy(result,reader[ivar]->fix->vstore,atom->nlocal*sizeof(double));
+      memcpy(result,reader[ivar]->fixstore->vstore,atom->nlocal*sizeof(double));
 
       int done = reader[ivar]->read_peratom();
       if (done) remove(ivar);
@@ -4275,7 +4284,7 @@ VarReader::VarReader(LAMMPS *lmp, char *name, char *file, int flag) :
   // allocate a new fix STORE, so they persist
   // id = variable-ID + VARIABLE_STORE, fix group = all
 
-  fix = NULL;
+  fixstore = NULL;
   id_fix = NULL;
   buffer = NULL;
 
@@ -4296,7 +4305,7 @@ VarReader::VarReader(LAMMPS *lmp, char *name, char *file, int flag) :
     newarg[3] = (char *) "0";
     newarg[4] = (char *) "1";
     modify->add_fix(5,newarg);
-    fix = (FixStore *) modify->fix[modify->nfix-1];
+    fixstore = (FixStore *) modify->fix[modify->nfix-1];
     delete [] newarg;
 
     buffer = new char[CHUNK*MAXLINE];
@@ -4311,7 +4320,7 @@ VarReader::~VarReader()
 
   // check modify in case all fixes have already been deleted
 
-  if (fix) {
+  if (fixstore) {
     if (modify) modify->delete_fix(id_fix);
     delete [] id_fix;
     delete [] buffer;
@@ -4367,7 +4376,7 @@ int VarReader::read_peratom()
   // set all per-atom values to 0.0
   // values that appear in file will overwrite this
 
-  double *vstore = fix->vstore;
+  double *vstore = fixstore->vstore;
 
   int nlocal = atom->nlocal;
   for (i = 0; i < nlocal; i++) vstore[i] = 0.0;

src/variable.h

@@ -27,6 +27,7 @@ class Variable : protected Pointers {
   void set(char *, int, char **);
   int next(int, char **);
   int find(char *);
+  void set_arrays(int);
   int equalstyle(int);
   int atomstyle(int);
   char *retrieve(char *);
@@ -103,7 +104,7 @@ class Variable : protected Pointers {
 
 class VarReader : protected Pointers {
  public:
-  class FixStore *fix;
+  class FixStore *fixstore;
   char *id_fix;
 
   VarReader(class LAMMPS *, char *, char *, int);