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

doc/atom_style.txt

@@ -71,14 +71,22 @@ All of the styles define point particles, except the {sphere},
 particles.
 
 For the {sphere} style, the particles are spheres and each stores a
-per-particle diameter and mass.  For the {ellipsoid} style, the
+per-particle diameter and mass.  If the diameter > 0.0, the particle
-particles are ellipsoids and each stores a per-particle shape vector
+is a finite-size sphere.  If the diameter = 0.0, it is a point
-with the 3 diamters of the ellipsoid.  For the {electron} style, the
+particle.
-particles representing electrons are 3d Gaussians with a specified
+
-position and bandwidth or uncertainty in position, which is
+For the {ellipsoid} style, the particles are ellipsoids and each
-represented by the eradius = electron size.  For the {peri} style, the
+stores a flag which indicates whether it is a finite-size ellipsoid or
-particles are spherical and each stores a per-particle mass and
+a point particle.  If it is an ellipsoid, it also stores a shape
-volume.
+vector with the 3 diamters of the ellipsoid and a quaternion 4-vector
+with its orientation.
+
+For the {electron} style, the particles representing electrons are 3d
+Gaussians with a specified position and bandwidth or uncertainty in
+position, which is represented by the eradius = electron size.
+
+For the {peri} style, the particles are spherical and each stores a
+per-particle mass and volume.
 
 :line
 

doc/create_atoms.txt

@@ -139,7 +139,7 @@ to change these values.
 charge = 0.0
 dipole moment magnitude = 0.0
 diameter = 1.0
-shape = 1.0 1.0 1.0
+shape = 0.0 0.0 0.0
 density = 1.0
 volume = 1.0
 velocity = 0.0 0.0 0.0
@@ -148,11 +148,19 @@ angular momentum = 0.0 0.0 0.0
 quaternion = (1,0,0,0)
 bonds, angles, dihedrals, impropers = none :ul
 
-Note that this means the {sphere} and {ellipsoid} atom styles set the
+Note that the {sphere} atom style sets the default particle diameter
-diameter/shape and density to 1.0 and thus calculates a mass for the
+to 1.0 as well as the density.  This means the mass for the particle
-particle, which is PI/6 * diameter^3 = 0.5236.  The {peri} style sets
+is not 1.0, but is PI/6 * diameter^3 = 0.5236.
-the volume and density to 1.0 and thus also set the mass for the
+
-particle to 1.0.
+Note that the {ellipsoid} atom style sets the default particle shape
+to (0.0 0.0 0.0) and the density to 1.0 which means it is a point
+particle, not an ellipsoid, and has a mass of 1.0.
+
+Note that the {peri} style sets the default volume and density to 1.0
+and thus also set the mass for the particle to 1.0.
+
+The "set"_set.html command can be used to override many of these
+default settings.
 
 [Restrictions:]
 

doc/dump.txt

@@ -41,9 +41,7 @@ args = list of arguments for a particular style :l
 			  vx, vy, vz, fx, fy, fz,
                           q, mux, muy, muz, mu,
                           radius, omegax, omegay, omegaz,
-			  angmomx, angmomy, angmomz,
+			  angmomx, angmomy, angmomz, tqx, tqy, tqz,
-			  shapex,shapey, shapez,
-			  quatw, quati, quatj, quatk, tqx, tqy, tqz,
 			  spin, eradius, ervel, erforce,
 			  c_ID, c_ID\[N\], f_ID, f_ID\[N\], v_name :pre
 
@@ -63,8 +61,6 @@ args = list of arguments for a particular style :l
       radius = radius of spherical particle
       omegax,omegay,omegaz = angular velocity of extended particle
       angmomx,angmomy,angmomz = angular momentum of extended particle
-      shapex,shapey,shapez = 3 diameters of ellipsoidal particle
-      quatw,quati,quatj,quatk = quaternion components for aspherical particles
       tqx,tqy,tqz = torque on extended particles
       spin = electron spin
       eradius = electron radius
@@ -415,15 +411,6 @@ The {angmomx}, {angmomy}, and {angmomz} attributes are specific to
 extended aspherical particles that have an angular momentum.  Only
 the {ellipsoid} atom style defines this quantity.
 
-The {shapex}, {shapey}, and {shapez} attributes are specific to
-extended ellipsoidal particles that have a finite size and shape, such
-those defined with an atom style of {ellipsoidal}.
-
-The {quatw}, {quati}, {quatj}, {quatk} attributes are for aspherical
-particles defined with an atom style of {ellipsoid}.  They are the
-components of the quaternion that defines the orientation of the
-particle.
-
 The {tqx}, {tqy}, {tqz} attributes are for extended spherical or
 aspherical particles that can sustain a rotational torque due
 to interactions with other particles.

doc/fix_store_state.txt

@@ -21,8 +21,7 @@ input = one or more atom attributes :l
 		        vx, vy, vz, fx, fy, fz,
                         q, mux, muy, muz,
                         radius, omegax, omegay, omegaz,
-                        angmomx, angmomy, angmomz,
+                        angmomx, angmomy, angmomz, tqx, tqy, tqz
-		        quatw, quati, quatj, quatk, tqx, tqy, tqz
 			c_ID, c_ID\[N\], f_ID, f_ID\[N\], v_name :pre
 
       id = atom ID
@@ -40,7 +39,6 @@ input = one or more atom attributes :l
       radius = radius of extended spherical particle
       omegax,omegay,omegaz = angular velocity of extended particle
       angmomx,angmomy,angmomz = angular momentum of extended particle
-      quatw,quati,quatj,quatk = quaternion components for aspherical particles
       tqx,tqy,tqz = torque on extended particles
       c_ID = per-atom vector calculated by a compute with ID
       c_ID\[I\] = Ith column of per-atom array calculated by a compute with ID

doc/set.txt

@@ -136,23 +136,25 @@ many processors are being used.
 Keyword {quat} uses the specified values to create a quaternion
 (4-vector) that represents the orientation of the selected atoms.
 Note that particles defined by "atom_style ellipsoid"_atom_style.html
-have 3 shape paraeters whicha are used to specify the aspect ratios of
+have 3 shape parameters.  The 3 values must be non-zero for each
-an ellipsoidal particle, which is oriented by default with its x-axis
+particle set by this command.  They are used to specify the aspect
-along the simulation box's x-axis, and similarly for y and z.  If this
+ratios of an ellipsoidal particle, which is oriented by default with
-body is rotated (via the right-hand rule) by an angle theta around a
+its x-axis along the simulation box's x-axis, and similarly for y and
-unit rotation vector (a,b,c), then the quaternion that represents its
+z.  If this body is rotated (via the right-hand rule) by an angle
-new orientation is given by (cos(theta/2), a*sin(theta/2),
+theta around a unit rotation vector (a,b,c), then the quaternion that
-b*sin(theta/2), c*sin(theta/2)).  The theta and a,b,c values are the
+represents its new orientation is given by (cos(theta/2),
-arguments to the {quat} keyword.  LAMMPS normalizes the quaternion in
+a*sin(theta/2), b*sin(theta/2), c*sin(theta/2)).  The theta and a,b,c
-case (a,b,c) was not specified as a unit vector.  For 2d systems, the
+values are the arguments to the {quat} keyword.  LAMMPS normalizes the
-a,b,c values are ignored, since a rotation vector of (0,0,1) is the
+quaternion in case (a,b,c) was not specified as a unit vector.  For 2d
-only valid choice.
+systems, the a,b,c values are ignored, since a rotation vector of
+(0,0,1) is the only valid choice.
 
 Keyword {quat/random} randomizes the orientation of the quaternion of
 the selected atoms.  Random numbers are used in such a way that the
 orientation of a particular atom is the same, regardless of how many
 processors are being used.  For 2d systems, only orientations in the
-xy plane are generated.
+xy plane are generated.  As with keyword {quat}, the 3 shape values
+must be non-zero for each particle set by this command.
 
 Keyword {diameter} sets the size of the selected atoms.  The particles
 must be finite-size spheres as defined by the "atom_style
@@ -167,13 +169,13 @@ particles must be aspherical ellipsoids as defined by the "atom_style
 ellipsoid"_atom_style.html command.  The {Sx}, {Sy}, {Sz} settings are
 the 3 diameters of the ellipsoid in each direction.  All 3 can be set
 to the same value, which means the ellipsoid is effectively a sphere.
-Or then can all be set to 0.0 which means the particle will be treated
+They can also all be set to 0.0 which means the particle will be
-as a point particle.  Note that this command does not adjust the
+treated as a point particle.  Note that this command does not adjust
-particle mass, even if it was defined with a density, e.g. via the
+the particle mass, even if it was defined with a density, e.g. via the
 "read_data"_read_data.html command.
 
-Keyword {mas} sets the mass of all selected particles.  The
+Keyword {mass} sets the mass of all selected particles.  The particles
-particles must have a per-atom mass attribute, as defined by the
+must have a per-atom mass attribute, as defined by the
 "atom_style"_atom_style.html command.  See the "mass" command for how
 to set mass values on a per-type basis.
 
@@ -184,7 +186,7 @@ to set mass values on a per-type basis.  If the atom has a radius
 attribute (see "atom_style sphere"_atom_style.html) and its radius is
 non-zero, its mass is set from the density and particle volume.  The
 same is true if the atom has a shape attribute (see "atom_style
-ellipsoid"_atom_style.html) and its shape parameters are non-zero.
+ellipsoid"_atom_style.html) and its 3 shape parameters are non-zero.
 Otherwise the mass is set to the density value directly.
 
 Keyword {volume} sets the volume of all selected particles.

examples/ellipse/in.ellipse.gayberne

@@ -10,12 +10,11 @@ create_box   2 box
 create_atoms 1 box
 
 set	     group all type/fraction 2 0.1 95392
-set	     group all quat/random 18238
-
 set 	     type 1 mass 1.0
 set 	     type 2 mass 1.5
 set 	     type 1 shape 1 1 1
 set 	     type 2 shape 3 1 1
+set	     group all quat/random 18238
 
 compute	     rot all temp/asphere
 group	     spheroid type 1
@@ -41,11 +40,13 @@ thermo	     100
 
 timestep     0.002
 
+compute	     q all property/atom quatw quati quatj quatk
+
 dump	     1 all custom 100 dump.ellipse.gayberne &
-	     id type x y z quatw quati quatj quatk
+	     id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
 
 #dump	     1 all custom 100 dump.ellipse.resquared &
-#	     id type x y z quatw quati quatj quatk
+#	     id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
 
 fix	     1 all npt/asphere temp 2.0 2.0 0.1 iso 0.0 1.0 1.0 &
 	       mtk no pchain 0 tchain 1

examples/ellipse/in.ellipse.resquared

@@ -10,12 +10,11 @@ create_box   2 box
 create_atoms 1 box
 
 set	     group all type/fraction 2 0.1 95392
-set	     group all quat/random 18238
-
 set 	     type 1 mass 1.0
 set 	     type 2 mass 1.5
 set 	     type 1 shape 1 1 1
 set 	     type 2 shape 3 1 1
+set	     group all quat/random 18238
 
 compute	     rot all temp/asphere
 group	     spheroid type 1
@@ -41,11 +40,13 @@ thermo	     100
 
 timestep     0.002
 
+compute	     q all property/atom quatw quati quatj quatk
+
 #dump	     1 all custom 100 dump.ellipse.gayberne &
-#	     id type x y z quatw quati quatj quatk
+#	     id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
 
 dump	     1 all custom 100 dump.ellipse.resquared &
-	     id type x y z quatw quati quatj quatk
+	     id type x y z c_q[1] c_q[2] c_q[3] c_q[4]
 
 fix	     1 all npt/asphere temp 2.0 2.0 0.1 iso 0.0 1.0 1.0 &
 	       mtk no pchain 0 tchain 1

tools/restart2data.cpp

@@ -92,6 +92,7 @@ class Data {
   int style_hybrid,style_molecular,style_peri,style_sphere;
 
   bigint natoms;
+  bigint nellipsoids;
   bigint nbonds,nangles,ndihedrals,nimpropers;
   int ntypes,nbondtypes,nangletypes,ndihedraltypes,nimpropertypes;
   int bond_per_atom,angle_per_atom,dihedral_per_atom,improper_per_atom;
@@ -202,6 +203,7 @@ class Data {
   double *s0,*x0x,*x0y,*x0z;
   double *shapex,*shapey,*shapez;
   double *quatw,*quati,*quatj,*quatk,*angmomx,*angmomy,*angmomz;
+  int *ellipsoid;
   int *bond_type,*angle_type,*dihedral_type,*improper_type;
   int *bond_atom1,*bond_atom2;
   int *angle_atom1,*angle_atom2,*angle_atom3;
@@ -916,22 +918,25 @@ int atom_ellipsoid(double *buf, Data &data, int iatoms)
   data.vy[iatoms] = buf[m++];
   data.vz[iatoms] = buf[m++];
 
+  data.rmass[iatoms] = buf[m++];
+  data.angmomx[iatoms] = buf[m++];
+  data.angmomy[iatoms] = buf[m++];
+  data.angmomz[iatoms] = buf[m++];
+  data.ellipsoid[iatoms] = static_cast<int> (buf[m++]);
+
+  if (data.ellipsoid[iatoms]) {
+    data.nellipsoids++;
     data.shapex[iatoms] = buf[m++];
     data.shapey[iatoms] = buf[m++];
     data.shapez[iatoms] = buf[m++];
-  data.rmass[iatoms] = buf[m++];
-  if (data.shapex[iatoms] == 0.0) data.density[iatoms] = data.rmass[iatoms];
-  else 
-    data.density[iatoms] = data.rmass[iatoms] / 
-      (4.0*PI/3.0 * 
-       data.shapex[iatoms]*data.shapey[iatoms]*data.shapez[iatoms]);
     data.quatw[iatoms] = buf[m++];
     data.quati[iatoms] = buf[m++];
     data.quatj[iatoms] = buf[m++];
     data.quatk[iatoms] = buf[m++];
-  data.angmomx[iatoms] = buf[m++];
+    data.density[iatoms] = data.rmass[iatoms] / 
-  data.angmomy[iatoms] = buf[m++];
+      (4.0*PI/3.0 * 
-  data.angmomz[iatoms] = buf[m++];
+       data.shapex[iatoms]*data.shapey[iatoms]*data.shapez[iatoms]);
+  } else data.density[iatoms] = data.rmass[iatoms];
 
   return m;
 }
@@ -1218,18 +1223,19 @@ void allocate_full(Data &data)
 
 void allocate_ellipsoid(Data &data)
 {
-  data.shapex = new double[data.natoms];
-  data.shapey = new double[data.natoms];
-  data.shapez = new double[data.natoms];
   data.rmass = new double[data.natoms];
   data.density = new double[data.natoms];
-  data.quatw = new double[data.natoms];
-  data.quati = new double[data.natoms];
-  data.quatj = new double[data.natoms];
-  data.quatk = new double[data.natoms];
   data.angmomx = new double[data.natoms];
   data.angmomy = new double[data.natoms];
   data.angmomz = new double[data.natoms];
+  data.ellipsoid = new int[data.natoms];
+  data.quatw = new double[data.natoms];
+  data.shapex = new double[data.natoms];
+  data.shapey = new double[data.natoms];
+  data.shapez = new double[data.natoms];
+  data.quati = new double[data.natoms];
+  data.quatj = new double[data.natoms];
+  data.quatk = new double[data.natoms];
 }
 
 void allocate_sphere(Data &data)
@@ -2625,7 +2631,10 @@ void improper(FILE *fp, Data &data)
 // initialize Data
 // ---------------------------------------------------------------------
 
-Data::Data() {}
+Data::Data()
+{
+  nellipsoids = 0;
+}
 
 // ---------------------------------------------------------------------
 // print out stats on problem
@@ -2639,6 +2648,9 @@ void Data::stats()
   printf("  Ntimestep = " BIGINT_FORMAT "\n",ntimestep);
   printf("  Nprocs = %d\n",nprocs);
   printf("  Natoms = " BIGINT_FORMAT "\n",natoms);
+
+  if (nellipsoids) printf("  Nellipsoids = " BIGINT_FORMAT "\n",nellipsoids);
+
   printf("  Nbonds = " BIGINT_FORMAT "\n",nbonds);
   printf("  Nangles = " BIGINT_FORMAT "\n",nangles);
   printf("  Ndihedrals = " BIGINT_FORMAT "\n",ndihedrals);
@@ -2651,6 +2663,7 @@ void Data::stats()
   printf("  Angle style = %s\n",angle_style);
   printf("  Dihedral style = %s\n",dihedral_style);
   printf("  Improper style = %s\n",improper_style);
+
   printf("  Xlo xhi = %g %g\n",xlo,xhi);
   printf("  Ylo yhi = %g %g\n",ylo,yhi);
   printf("  Zlo zhi = %g %g\n",zlo,zhi);
@@ -2670,6 +2683,7 @@ void Data::write(FILE *fp, FILE *fp2)
 	  BIGINT_FORMAT ", procs = %d\n\n",ntimestep,nprocs);
 
   fprintf(fp,BIGINT_FORMAT " atoms\n",natoms);
+  if (nellipsoids) fprintf(fp,BIGINT_FORMAT " ellipsoids\n",nellipsoids);
   if (nbonds) fprintf(fp,BIGINT_FORMAT " bonds\n",nbonds);
   if (nangles) fprintf(fp,BIGINT_FORMAT " angles\n",nangles);
   if (ndihedrals) fprintf(fp,BIGINT_FORMAT " dihedrals\n",ndihedrals);
@@ -3240,6 +3254,17 @@ void Data::write(FILE *fp, FILE *fp2)
       }
   }
  
+  if (nellipsoids) {
+    fprintf(fp,"\nEllipsoids\n\n");
+    for (uint64_t i = 0; i < natoms; i++) {
+      if (ellipsoid[i])
+	fprintf(fp,"%d %-1.16e %-1.16e %-1.16e "
+		"%-1.16e %-1.16e %-1.16e %-1.16e \n",
+		tag[i],2.0*shapex[i],2.0*shapey[i],2.0*shapez[i],
+		quatw[i],quati[i],quatj[i],quatk[i]);
+    }
+  }
+
   if (nbonds) {
     fprintf(fp,"\nBonds\n\n");
     for (uint64_t i = 0; i < nbonds; i++)
@@ -3302,16 +3327,16 @@ void Data::write_atom_charge(FILE *fp, int i, int ix, int iy, int iz)
 
 void Data::write_atom_dipole(FILE *fp, int i, int ix, int iy, int iz)
 {
-  fprintf(fp,"%d %d %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %d %d %d",
+  fprintf(fp,"%d %d %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e "
+	  "%-1.16e %-1.16e %d %d %d",
 	  tag[i],type[i],q[i],x[i],y[i],z[i],
 	  mux[i],muy[i],muz[i],ix,iy,iz);
 }
 
 void Data::write_atom_ellipsoid(FILE *fp, int i, int ix, int iy, int iz)
 {
-  fprintf(fp,"%d %d %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %d %d %d",
+  fprintf(fp,"%d %d %d %-1.16e %-1.16e %-1.16e %-1.16e %d %d %d",
-	  tag[i],type[i],shapex[i],shapey[i],shapez[i],density[i],
+	  tag[i],type[i],ellipsoid[i],density[i],x[i],y[i],z[i],ix,iy,iz);
-	  x[i],y[i],z[i],quatw[i],quati[i],quatj[i],quatk[i],ix,iy,iz);
 }
 
 void Data::write_atom_full(FILE *fp, int i, int ix, int iy, int iz)
@@ -3367,9 +3392,7 @@ void Data::write_atom_dipole_extra(FILE *fp, int i)
 
 void Data::write_atom_ellipsoid_extra(FILE *fp, int i)
 {
-  fprintf(fp," %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e %-1.16e",
+  fprintf(fp," %d %-1.16e",ellipsoid[i],density[i]);
-	  shapex[i],shapey[i],shapez[i],density[i],
-	  quatw[i],quati[i],quatj[i],quatk[i]);
 }
 
 void Data::write_atom_full_extra(FILE *fp, int i)