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

src/BODY/body_nparticle.cpp

@@ -16,8 +16,9 @@
 #include "math_extra.h"
 #include "atom_vec_body.h"
 #include "atom.h"
-#include "error.h"
 #include "force.h"
+#include "memory.h"
+#include "error.h"
 
 using namespace LAMMPS_NS;
 
@@ -42,6 +43,9 @@ BodyNparticle::BodyNparticle(LAMMPS *lmp, int narg, char **arg) :
 
   icp = new MyPoolChunk<int>(1,1);
   dcp = new MyPoolChunk<double>(3*nmin,3*nmax);
+
+  memory->create(imflag,nmax,"body/nparticle:imflag");
+  memory->create(imdata,nmax,3,"body/nparticle:imdata");
 }
 
 /* ---------------------------------------------------------------------- */
@@ -50,6 +54,8 @@ BodyNparticle::~BodyNparticle()
 {
   delete icp;
   delete dcp;
+  memory->destroy(imflag);
+  memory->destroy(imdata);
 }
 
 /* ---------------------------------------------------------------------- */
@@ -220,3 +226,30 @@ void BodyNparticle::output(int ibonus, int m, double *values)
   values[1] += x[1];
   values[2] += x[2];
 }
+
+/* ---------------------------------------------------------------------- */
+
+int BodyNparticle::image(int ibonus, int *&ivec, double **&darray)
+{
+  double p[3][3];
+  double *x;
+
+  AtomVecBody::Bonus *bonus = &avec->bonus[ibonus];
+  int n = bonus->ivalue[0];
+
+  for (int i = 0; i < n; i++) {
+    imflag[i] = 0;
+    //imflag[i] = SPHERE;
+    MathExtra::quat_to_mat(bonus->quat,p);
+    MathExtra::matvec(p,&bonus->dvalue[3*i],imdata[i]);
+
+    x = atom->x[bonus->ilocal];
+    imdata[i][0] += x[0];
+    imdata[i][1] += x[1];
+    imdata[i][2] += x[2];
+  }
+
+  ivec = imflag;
+  darray = imdata;
+  return n;
+}

src/BODY/body_nparticle.h

@@ -39,6 +39,12 @@ class BodyNparticle : public Body {
   int noutrow(int);
   int noutcol();
   void output(int, int, double *);
+  int image(int, int *&, double **&);
+
+ private:
+  int *imflag;
+  double **imdata;
+
 };
 
 }

src/atom.cpp

@@ -45,7 +45,6 @@ using namespace MathConst;
 #define DELTA_MEMSTR 1024
 #define EPSILON 1.0e-6
 #define CUDA_CHUNK 3000
-#define MAXBODY 20       // max # of lines in one body, also in ReadData class
 
 enum{LAYOUT_UNIFORM,LAYOUT_NONUNIFORM,LAYOUT_TILED};    // several files
 
@@ -713,6 +712,29 @@ int Atom::count_words(const char *line)
   return n;
 }
 
+/* ----------------------------------------------------------------------
+   count and return words in a single line using provided copy buf
+   make copy of line before using strtok so as not to change line
+   trim anything from '#' onward
+------------------------------------------------------------------------- */
+
+int Atom::count_words(const char *line, char *copy)
+{
+  strcpy(copy,line);
+
+  char *ptr;
+  if ((ptr = strchr(copy,'#'))) *ptr = '\0';
+
+  if (strtok(copy," \t\n\r\f") == NULL) {
+    memory->destroy(copy);
+    return 0;
+  }
+  int n = 1;
+  while (strtok(NULL," \t\n\r\f")) n++;
+
+  return n;
+}
+
 /* ----------------------------------------------------------------------
    deallocate molecular topology arrays
    done before realloc with (possibly) new 2nd dimension set to
@@ -756,7 +778,7 @@ void Atom::deallocate_topology()
 }
 
 /* ----------------------------------------------------------------------
-   unpack n lines from Atom section of data file
+   unpack N lines from Atom section of data file
    call style-specific routine to parse line
 ------------------------------------------------------------------------- */
 
@@ -900,7 +922,7 @@ void Atom::data_atoms(int n, char *buf, tagint id_offset, int type_offset,
 }
 
 /* ----------------------------------------------------------------------
-   unpack n lines from Velocity section of data file
+   unpack N lines from Velocity section of data file
    check that atom IDs are > 0 and <= map_tag_max
    call style-specific routine to parse line
 ------------------------------------------------------------------------- */
@@ -1240,7 +1262,7 @@ void Atom::data_impropers(int n, char *buf, int *count, tagint id_offset,
 }
 
 /* ----------------------------------------------------------------------
-   unpack n lines from atom-style specific section of data file
+   unpack N lines from atom-style specific bonus section of data file
    check that atom IDs are > 0 and <= map_tag_max
    call style-specific routine to parse line
 ------------------------------------------------------------------------- */
@@ -1287,7 +1309,8 @@ void Atom::data_bonus(int n, char *buf, AtomVec *avec_bonus, tagint id_offset)
 }
 
 /* ----------------------------------------------------------------------
-   unpack n lines from atom-style specific section of data file
+   unpack N bodies from Bodies section of data file
+   each body spans multiple lines
    check that atom IDs are > 0 and <= map_tag_max
    call style-specific routine to parse line
 ------------------------------------------------------------------------- */
@@ -1297,8 +1320,10 @@ void Atom::data_bodies(int n, char *buf, AtomVecBody *avec_body,
 {
   int j,m,tagdata,ninteger,ndouble;
 
-  char **ivalues = new char*[10*MAXBODY];
+  int maxint = 0;
-  char **dvalues = new char*[10*MAXBODY];
+  int maxdouble = 0;
+  char **ivalues = NULL;
+  char **dvalues = NULL;
 
   // loop over lines of body data
   // tokenize the lines into ivalues and dvalues
@@ -1307,9 +1332,21 @@ void Atom::data_bodies(int n, char *buf, AtomVecBody *avec_body,
   for (int i = 0; i < n; i++) {
     if (i == 0) tagdata = ATOTAGINT(strtok(buf," \t\n\r\f")) + id_offset;
     else tagdata = ATOTAGINT(strtok(NULL," \t\n\r\f")) + id_offset;
+
     ninteger = atoi(strtok(NULL," \t\n\r\f"));
     ndouble = atoi(strtok(NULL," \t\n\r\f"));
 
+    if (ninteger > maxint) {
+      delete [] ivalues;
+      maxint = ninteger;
+      ivalues = new char*[maxint];
+    }
+    if (ndouble > maxdouble) {
+      delete [] dvalues;
+      maxdouble = ndouble;
+      dvalues = new char*[maxdouble];
+    }
+
     for (j = 0; j < ninteger; j++)
       ivalues[j] = strtok(NULL," \t\n\r\f");
     for (j = 0; j < ndouble; j++)

src/atom.h

@@ -207,6 +207,7 @@ class Atom : protected Pointers {
 
   int parse_data(const char *);
   int count_words(const char *);
+  int count_words(const char *, char *);
 
   void deallocate_topology();
 

src/body.h

@@ -45,6 +45,7 @@ class Body : protected Pointers {
   virtual int noutrow(int) = 0;
   virtual int noutcol() = 0;
   virtual void output(int, int, double *) = 0;
+  virtual int image(int, double, double, int *&, double **&) = 0;
 };
 
 }

src/dump_image.cpp

@@ -19,6 +19,8 @@
 #include "image.h"
 #include "atom.h"
 #include "atom_vec.h"
+#include "atom_vec_body.h"
+#include "body.h"
 #include "molecule.h"
 #include "domain.h"
 #include "group.h"
@@ -36,6 +38,7 @@ using namespace MathConst;
 #define BIG 1.0e20
 
 enum{NUMERIC,ATOM,TYPE,ELEMENT,ATTRIBUTE};
+enum{SPHERE,LINE};           // also in Body child classes
 enum{STATIC,DYNAMIC};
 enum{NO,YES};
 
@@ -103,6 +106,7 @@ DumpImage::DumpImage(LAMMPS *lmp, int narg, char **arg) :
   // set defaults for optional args
 
   atomflag = YES;
+  bodyflag = NO;
   if (atom->nbondtypes == 0) bondflag = NO;
   else {
     bondflag = YES;
@@ -142,6 +146,15 @@ DumpImage::DumpImage(LAMMPS *lmp, int narg, char **arg) :
       else error->all(FLERR,"Illegal dump image command");
       iarg += 2;
 
+    } else if (strcmp(arg[iarg],"body") == 0) {
+      if (iarg+4 > narg) error->all(FLERR,"Illegal dump image command");
+      if (strcmp(arg[iarg+1],"yes") == 0) bodyflag = YES;
+      else if (strcmp(arg[iarg+1],"no") == 0) bodyflag = NO;
+      else error->all(FLERR,"Illegal dump image command");
+      bodyflag1 = force->numeric(FLERR,arg[iarg+2]);
+      bodyflag2 = force->numeric(FLERR,arg[iarg+3]);
+      iarg += 4;
+
     } else if (strcmp(arg[iarg],"bond") == 0) {
       if (iarg+3 > narg) error->all(FLERR,"Illegal dump image command");
       if (atom->nbondtypes == 0)
@@ -320,6 +333,14 @@ DumpImage::DumpImage(LAMMPS *lmp, int narg, char **arg) :
     } else error->all(FLERR,"Illegal dump image command");
   }
 
+  // error check for bodyflag
+
+  if (bodyflag) {
+    AtomVecBody *avec = (AtomVecBody *) atom->style_match("body");
+    if (!avec) error->all(FLERR,"Dump image body yes requires atom style body");
+    bptr = avec->bptr;
+  }
+
   // allocate image buffer now that image size is known
 
   image->buffers();
@@ -658,9 +679,11 @@ void DumpImage::view_params()
 
 void DumpImage::create_image()
 {
-  int i,j,m,n,itype,atom1,atom2,imol,iatom,btype;
+  int i,j,k,m,n,itype,atom1,atom2,imol,iatom,btype,ibonus;
   tagint tagprev;
   double diameter,delx,dely,delz;
+  int *bodyvec;
+  double **bodyarray;
   double *color,*color1,*color2;
   double xmid[3];
 
@@ -668,6 +691,7 @@ void DumpImage::create_image()
 
   if (atomflag) {
     double **x = atom->x;
+    int *body = atom->body;
 
     m = 0;
     for (i = 0; i < nchoose; i++) {
@@ -695,7 +719,20 @@ void DumpImage::create_image()
         diameter = buf[m+1];
       }
 
-      image->draw_sphere(x[j],color,diameter);
+      if (!body || !bodyflag || body[j] < 0)
+        image->draw_sphere(x[j],color,diameter);
+      else {
+        ibonus = body[i];
+        n = bptr->image(ibonus,bodyflag1,bodyflag2,bodyvec,bodyarray);
+        for (k = 0; k < n; k++) {
+          if (bodyvec[k] == SPHERE)
+            image->draw_sphere(bodyarray[k],color,bodyarray[k][3]);
+          else if (bodyvec[k] == LINE)
+            image->draw_cylinder(&bodyarray[k][0],&bodyarray[k][3],
+                                 color,bodyarray[k][6],3);
+        }
+      }
+
       m += size_one;
     }
   }

src/dump_image.h

@@ -40,6 +40,8 @@ class DumpImage : public DumpCustom {
   int acolor,adiam;                // what determines color/diam of atoms
   double adiamvalue;               // atom diameter value
   int atomflag,bondflag;           // 0/1 for draw atoms,bonds
+  int bodyflag;                    // 0/1 for draw atoms as bodies
+  double bodyflag1,bodyflag2;      // user params for drawing bodies
   int bcolor,bdiam;                // what determines color/diam of bonds
   double bdiamvalue;               // bond diameter value
   char *thetastr,*phistr;          // variables for view theta,phi
@@ -63,6 +65,7 @@ class DumpImage : public DumpCustom {
   double **colortype,**colorelement,**bcolortype;   // per-type colors
 
   class Image *image;              // class that renders each image
+  class Body *bptr;                // class for Body particles
 
   int *chooseghost;                // extended choose array for comm
   double **bufcopy;                // buffer for communicating bond/atom info

src/read_data.cpp

@@ -52,7 +52,7 @@ using namespace LAMMPS_NS;
 #define LB_FACTOR 1.1
 #define CHUNK 1024
 #define DELTA 4            // must be 2 or larger
-#define MAXBODY 20         // max # of lines in one body, also in Atom class
+#define MAXBODY 32         // max # of lines in one body
 
                            // customize for new sections
 #define NSECTIONS 25       // change when add to header::section_keywords
@@ -73,6 +73,7 @@ ReadData::ReadData(LAMMPS *lmp) : Pointers(lmp)
 {
   MPI_Comm_rank(world,&me);
   line = new char[MAXLINE];
+  copy = new char[MAXLINE];
   keyword = new char[MAXLINE];
   style = new char[MAXLINE];
   buffer = new char[CHUNK*MAXLINE];
@@ -98,6 +99,7 @@ ReadData::ReadData(LAMMPS *lmp) : Pointers(lmp)
 ReadData::~ReadData()
 {
   delete [] line;
+  delete [] copy;
   delete [] keyword;
   delete [] style;
   delete [] buffer;
@@ -1462,12 +1464,12 @@ void ReadData::bonus(bigint nbonus, AtomVec *ptr, const char *type)
    read all body data
    variable amount of info per body, described by ninteger and ndouble
    to find atoms, must build atom map if not a molecular system
-   if not firstpass, just read but no processing of data
+   if not firstpass, just read past data, but no processing of data
 ------------------------------------------------------------------------- */
 
 void ReadData::bodies(int firstpass)
 {
-  int i,m,nchunk,nline,nmax,ninteger,ndouble,tmp,onebody;
+  int i,m,nchunk,nline,nmax,ninteger,ndouble,nword,onebody,tmp;
   char *eof;
 
   int mapflag = 0;
@@ -1498,19 +1500,38 @@ void ReadData::bodies(int firstpass)
         sscanf(&buffer[m],"%d %d %d",&tmp,&ninteger,&ndouble);
         m += strlen(&buffer[m]);
 
+        // read lines one at a time into buffer
+        // make copy of line and count words
+        // count to ninteger and ndouble until have enough lines
+
         onebody = 0;
-        if (ninteger) onebody += (ninteger-1)/10 + 1;
+
-        if (ndouble) onebody += (ndouble-1)/10 + 1;
+        nword = 0;
+        while (nword < ninteger) {
+          eof = fgets(&buffer[m],MAXLINE,fp);
+          if (eof == NULL) error->one(FLERR,"Unexpected end of data file");
+          nword += atom->count_words(&buffer[m],copy);
+          m += strlen(&buffer[m]);
+          onebody++;
+        }
+        if (nword > ninteger) 
+          error->one(FLERR,"Too many value in body lines in data file");
+
+        nword = 0;
+        while (nword < ndouble) {
+          eof = fgets(&buffer[m],MAXLINE,fp);
+          if (eof == NULL) error->one(FLERR,"Unexpected end of data file");
+          nword += atom->count_words(&buffer[m],copy);
+          m += strlen(&buffer[m]);
+          onebody++;
+        }
+        if (nword > ndouble) 
+          error->one(FLERR,"Too many value in body lines in data file");
+
         if (onebody+1 > MAXBODY)
           error->one(FLERR,
                      "Too many lines in one body in data file - boost MAXBODY");
 
-        for (i = 0; i < onebody; i++) {
-          eof = fgets(&buffer[m],MAXLINE,fp);
-          if (eof == NULL) error->one(FLERR,"Unexpected end of data file");
-          m += strlen(&buffer[m]);
-        }
-
         nchunk++;
         nline += onebody+1;
       }

src/read_data.h

@@ -33,7 +33,7 @@ class ReadData : protected Pointers {
 
  private:
   int me,compressed;
-  char *line,*keyword,*buffer,*style;
+  char *line,*copy,*keyword,*buffer,*style;
   FILE *fp;
   char **arg;
   int narg,maxarg;