new option *add keep* added to read_dump

doc/src/read_dump.txt

@@ -27,7 +27,7 @@ keyword = {box} or {replace} or {purge} or {trim} or {add} or {label} or {scaled
   {replace} value = {yes} or {no} = overwrite atoms with dump atoms
   {purge} value = {yes} or {no} = delete all atoms before adding dump atoms
   {trim} value = {yes} or {no} = trim atoms not in dump snapshot
-  {add} value = {yes} or {no} = add new dump atoms to system
+  {add} value = {yes} or {keep} or {no} = add new dump atoms to system
   {label} value = field column
     field = one of the listed fields or {id} or {type}
     column = label on corresponding column in dump file
@@ -169,9 +169,9 @@ file, in order, until a match is found.
 
 The dump file must also contain atom IDs, with a column label of "id".
 
-If the {add} keyword is specified with a value of {yes}, as discussed
+If the {add} keyword is specified with a value of {yes} or {keep}, as
-below, the dump file must contain atom types, with a column label of
+discussed below, the dump file must contain atom types, with a column
-"type".
+label of "type".
 
 If a column label you want to read from the dump file is not a match
 to a specified field, the {label} keyword can be used to specify the
@@ -231,26 +231,31 @@ IDs that are in the current system but not in the dump snapshot are
 deleted.  These atoms are unaffected if the {trim} keyword is
 specified with a {no} value.
 
-If the {add} keyword is specified with a {yes} value, then atoms with
+If the {add} keyword is specified with a {no} value (default), then
-IDs that are in the dump snapshot, but not in the current system are
+dump file atoms with IDs that are not in the current system are not
-added to the system.  These dump atoms are ignored if the {add}
+added to the system.  They are simply ignored.
-keyword is specified with a {no} value.
 
-Note that atoms added via the {add} keyword will have only the
+If a {yes} value is specified, the atoms with new IDs are added to the
-attributes read from the dump file due to the {field} arguments.  If
+system but their atom IDs are not preserved.  Instead, after all the
-{x} or {y} or {z} is not specified as a field, a value of 0.0 is used
+atoms are added, new IDs are assigned to them in the same manner as is
-for added atoms.  Added atoms must have an atom type, so this value
+described for the "create_atoms"_create_atoms.html command.  Basically
-must appear in the dump file.
+the largest existing atom ID in the system is identified, and all the
+added atoms are assigned IDs that consecutively follow the largest ID.
+
+If a {keep} value is specified, the atoms with new IDs are added to
+the system and their atom IDs are preserved.  This may lead to
+non-contiguous IDs for the combined system.
+
+Note that atoms added via the {add} keyword will only have the
+attributes read from the dump file due to the {field} arguments.  For
+example, if {x} or {y} or {z} or {q} is not specified as a field, a
+value of 0.0 is used for added atoms.  Added atoms must have an atom
+type, so this value must appear in the dump file.
 
 Any other attributes (e.g. charge or particle diameter for spherical
 particles) will be set to default values, the same as if the
 "create_atoms"_create_atoms.html command were used.
 
-Note that atom IDs are not preserved for new dump snapshot atoms added
-via the {add} keyword.  The procedure for assigning new atom IDS to
-added atoms is the same as is described for the
-"create_atoms"_create_atoms.html command.
-
 :line
 
 Atom coordinates read from the dump file are first converted into

src/read_dump.cpp

@@ -50,6 +50,7 @@ using namespace LAMMPS_NS;
 
 enum{ID,TYPE,X,Y,Z,VX,VY,VZ,Q,IX,IY,IZ,FX,FY,FZ};
 enum{UNSET,NOSCALE_NOWRAP,NOSCALE_WRAP,SCALE_NOWRAP,SCALE_WRAP};
+enum{NOADD,YESADD,KEEPADD};
 
 /* ---------------------------------------------------------------------- */
 
@@ -366,7 +367,7 @@ void ReadDump::header(int fieldinfo)
   // triclinic_snap < 0 means no box info in file
 
   if (triclinic_snap < 0 && boxflag > 0)
-    error->all(FLERR,"No box information in dump. You have to use 'box no'");
+    error->all(FLERR,"No box information in dump, must use 'box no'");
   if (triclinic_snap >= 0) {
     if ((triclinic_snap && !triclinic) ||
         (!triclinic_snap && triclinic))
@@ -479,9 +480,9 @@ void ReadDump::atoms()
     nread += nchunk;
   }
 
-  // if addflag set, add tags to new atoms if possible
+  // if addflag = YESADD, assign IDs to new snapshot atoms
 
-  if (addflag) {
+  if (addflag == YESADD) {
     bigint nblocal = atom->nlocal;
     MPI_Allreduce(&nblocal,&atom->natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
     if (atom->natoms < 0 || atom->natoms >= MAXBIGINT)
@@ -623,7 +624,7 @@ int ReadDump::fields_and_keywords(int narg, char **arg)
   replaceflag = 1;
   purgeflag = 0;
   trimflag = 0;
-  addflag = 0;
+  addflag = NOADD;
   for (int i = 0; i < nfield; i++) fieldlabel[i] = NULL;
   scaleflag = 0;
   wrapflag = 1;
@@ -655,8 +656,9 @@ int ReadDump::fields_and_keywords(int narg, char **arg)
       iarg += 2;
     } else if (strcmp(arg[iarg],"add") == 0) {
       if (iarg+2 > narg) error->all(FLERR,"Illegal read_dump command");
-      if (strcmp(arg[iarg+1],"yes") == 0) addflag = 1;
+      if (strcmp(arg[iarg+1],"yes") == 0) addflag = YESADD;
-      else if (strcmp(arg[iarg+1],"no") == 0) addflag = 0;
+      else if (strcmp(arg[iarg+1],"no") == 0) addflag = NOADD;
+      else if (strcmp(arg[iarg+1],"keep") == 0) addflag = KEEPADD;
       else error->all(FLERR,"Illegal read_dump command");
       iarg += 2;
     } else if (strcmp(arg[iarg],"label") == 0) {
@@ -695,6 +697,8 @@ int ReadDump::fields_and_keywords(int narg, char **arg)
 
   if (purgeflag && (replaceflag || trimflag))
     error->all(FLERR,"If read_dump purges it cannot replace or trim");
+  if (addflag == KEEPADD && atom->tag_enable == 0)
+    error->all(FLERR,"Read_dump cannot add keep with atom IDs disabled");
 
   return narg-iarg;
 }
@@ -740,12 +744,13 @@ void ReadDump::process_atoms(int n)
 {
   int i,m,ifield,itype;
   int xbox,ybox,zbox;
-  tagint tag;
+  tagint mtag;
 
   double **x = atom->x;
   double **v = atom->v;
   double *q = atom->q;
   double **f = atom->f;
+  tagint *tag = atom->tag;
   imageint *image = atom->image;
   int nlocal = atom->nlocal;
   tagint map_tag_max = atom->map_tag_max;
@@ -758,8 +763,8 @@ void ReadDump::process_atoms(int n)
     // NOTE: atom ID in fields is stored as double, not as ubuf
     //       so can only cast it to tagint, thus cannot be full 64-bit ID
 
-    tag = static_cast<tagint> (fields[i][0]);
+    mtag = static_cast<tagint> (fields[i][0]);
-    if (tag <= map_tag_max) m = atom->map(tag);
+    if (mtag <= map_tag_max) m = atom->map(mtag);
     else m = -1;
     if (m < 0 || m >= nlocal) continue;
 
@@ -835,8 +840,15 @@ void ReadDump::process_atoms(int n)
   // create any atoms in chunk that no processor owned
   // add atoms in round-robin sequence on processors
   // cannot do it geometrically b/c dump coords may not be in simulation box
+  // check that dump file snapshot has atom type field
 
-  if (!addflag) return;
+  if (addflag == NOADD) return;
+
+  int tflag = 0;
+  for (ifield = 0; ifield < nfield; ifield++)
+    if (fieldtype[ifield] == TYPE) tflag = 1;
+  if (!tflag) 
+    error->all(FLERR,"Cannot add atoms if dump file does not store atom type");
 
   MPI_Allreduce(ucflag,ucflag_all,n,MPI_INT,MPI_SUM,world);
 
@@ -881,6 +893,7 @@ void ReadDump::process_atoms(int n)
     atom->avec->create_atom(itype,one);
     nadd++;
 
+    tag = atom->tag;
     v = atom->v;
     q = atom->q;
     image = atom->image;
@@ -889,8 +902,12 @@ void ReadDump::process_atoms(int n)
 
     xbox = ybox = zbox = 0;
 
-    for (ifield = 1; ifield < nfield; ifield++) {
+    for (ifield = 0; ifield < nfield; ifield++) {
       switch (fieldtype[ifield]) {
+      case ID:
+        if (addflag == KEEPADD) 
+          tag[m] = static_cast<tagint> (fields[i][ifield]);
+        break;
       case VX:
         v[m][0] = fields[i][ifield];
         break;