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

This commit is contained in:
sjplimp
2011-04-21 20:05:45 +00:00
parent 0f41ccb943
commit 428bdd0b09
8 changed files with 48 additions and 15 deletions

View File

@ -385,15 +385,15 @@ void PairLubricate::init_style()
int irequest = neighbor->request(this);
// require that atom radii are identical within each type require
// monodisperse system with same radii for all types
// require that atom radii are identical within each type
// require monodisperse system with same radii for all types
double rad,radtype;
for (int i = 1; i <= atom->ntypes; i++) {
if (!atom->radius_consistency(i,radtype))
error->all("Pair colloid requires atoms with same type have same radius");
error->all("Pair lubricate requires monodisperse particles");
if (i > 1 && radtype != rad)
error->all("Pair colloid requires mono-disperse particles");
error->all("Pair lubricate requires monodisperse particles");
rad = radtype;
}
}

View File

@ -12,7 +12,7 @@
------------------------------------------------------------------------- */
/* ----------------------------------------------------------------------
Contributing authors: Randy Schunk (SNL)
Contributing authors: Randy Schunk (Sandia)
------------------------------------------------------------------------- */
#include "math.h"

View File

@ -767,7 +767,7 @@ void AtomVecEllipsoid::unpack_border(int n, int first, double *buf)
type[i] = static_cast<int> (buf[m++]);
mask[i] = static_cast<int> (buf[m++]);
ellipsoid[i] = static_cast<int> (buf[m++]);
if (ellipsoid[i] < 0) ellipsoid[i] = 0;
if (ellipsoid[i] == 0) ellipsoid[i] = -1;
else {
j = nlocal_bonus + nghost_bonus;
if (j == nmax_bonus) grow_bonus();
@ -805,9 +805,8 @@ void AtomVecEllipsoid::unpack_border_vel(int n, int first, double *buf)
type[i] = static_cast<int> (buf[m++]);
mask[i] = static_cast<int> (buf[m++]);
ellipsoid[i] = static_cast<int> (buf[m++]);
if (ellipsoid[i] < 0) ellipsoid[i] = 0;
if (ellipsoid[i] == 0) ellipsoid[i] = -1;
else {
j = nlocal_bonus + nghost_bonus;
if (j == nmax_bonus) grow_bonus();
shape = bonus[j].shape;
quat = bonus[j].quat;
@ -842,7 +841,7 @@ int AtomVecEllipsoid::unpack_border_hybrid(int n, int first, double *buf)
last = first + n;
for (i = first; i < last; i++) {
ellipsoid[i] = static_cast<int> (buf[m++]);
if (ellipsoid[i] < 0) ellipsoid[i] = 0;
if (ellipsoid[i] == 0) ellipsoid[i] = -1;
else {
j = nlocal_bonus + nghost_bonus;
if (j == nmax_bonus) grow_bonus();
@ -933,7 +932,8 @@ int AtomVecEllipsoid::unpack_exchange(double *buf)
angmom[nlocal][2] = buf[m++];
ellipsoid[nlocal] = static_cast<int> (buf[m++]);
if (ellipsoid[nlocal]) {
if (ellipsoid[nlocal] == 0) ellipsoid[nlocal] = -1;
else {
if (nlocal_bonus == nmax_bonus) grow_bonus();
double *shape = bonus[nlocal_bonus].shape;
double *quat = bonus[nlocal_bonus].quat;
@ -1057,7 +1057,8 @@ int AtomVecEllipsoid::unpack_restart(double *buf)
angmom[nlocal][2] = buf[m++];
ellipsoid[nlocal] = static_cast<int> (buf[m++]);
if (ellipsoid[nlocal]) {
if (ellipsoid[nlocal] == 0) ellipsoid[nlocal] = -1;
else {
if (nlocal_bonus == nmax_bonus) grow_bonus();
double *shape = bonus[nlocal_bonus].shape;
double *quat = bonus[nlocal_bonus].quat;

View File

@ -19,6 +19,7 @@
#include "stdlib.h"
#include "math.h"
#include "fix_nh.h"
#include "math_extra.h"
#include "atom.h"
#include "force.h"
#include "comm.h"
@ -31,7 +32,6 @@
#include "domain.h"
#include "memory.h"
#include "error.h"
#include "math_extra.h"
using namespace LAMMPS_NS;

View File

@ -468,7 +468,7 @@ double FixStoreState::memory_usage()
void FixStoreState::grow_arrays(int nmax)
{
memory->grow(values,nmax,nvalues,"fix_store:values");
memory->grow(values,nmax,nvalues,"store/state:values");
if (nvalues == 1) {
if (nmax) vector_atom = &values[0][0];
else vector_atom = NULL;

View File

@ -428,8 +428,8 @@ void quat_to_mat_trans(const double *quat, double mat[3][3])
/* ----------------------------------------------------------------------
compute space-frame inertia tensor of an ellipsoid
quat = orientiation quaternion of ellipsoid
radii = 3 radii of ellipsoid
quat = orientiation quaternion of ellipsoid
return symmetric inertia tensor as 6-vector in Voigt notation
------------------------------------------------------------------------- */
@ -454,6 +454,36 @@ void inertia_ellipsoid(double *radii, double *quat, double mass,
inertia[5] = tensor[0][1];
}
/* ----------------------------------------------------------------------
compute space-frame inertia tensor of a line segment in 2d
length = length of line
theta = orientiation of line
return symmetric inertia tensor as 6-vector in Voigt notation
------------------------------------------------------------------------- */
void inertia_line(double length, double theta, double mass, double *inertia)
{
double p[3][3],ptrans[3][3],itemp[3][3],tensor[3][3];
double q[4],idiag[3];
q[0] = cos(0.5*theta);
q[1] = q[2] = 0.0;
q[3] = sin(0.5*theta);
MathExtra::quat_to_mat(q,p);
MathExtra::quat_to_mat_trans(q,ptrans);
idiag[0] = 0.0;
idiag[1] = 1.0/12.0 * mass * length*length;
idiag[2] = 1.0/12.0 * mass * length*length;
MathExtra::diag_times3(idiag,ptrans,itemp);
MathExtra::times3(p,itemp,tensor);
inertia[0] = tensor[0][0];
inertia[1] = tensor[1][1];
inertia[2] = tensor[2][2];
inertia[3] = tensor[1][2];
inertia[4] = tensor[0][2];
inertia[5] = tensor[0][1];
}
/* ----------------------------------------------------------------------
compute space-frame inertia tensor of a triangle
v0,v1,v2 = 3 vertices of triangle

View File

@ -114,6 +114,8 @@ namespace MathExtra {
void inertia_ellipsoid(double *shape, double *quat, double mass,
double *inertia);
void inertia_line(double length, double theta, double mass,
double *inertia);
void inertia_triangle(double *v0, double *v1, double *v2,
double mass, double *inertia);
}

View File

@ -41,7 +41,7 @@ using namespace LAMMPS_NS;
#define DELTA 4 // must be 2 or larger
// customize for new sections
#define NSECTIONS 23 // change when add to header::section_keywords
#define NSECTIONS 21 // change when add to header::section_keywords
#define MIN(a,b) ((a) < (b) ? (a) : (b))
#define MAX(a,b) ((a) > (b) ? (a) : (b))