diff --git a/src/RIGID/fix_rigid.cpp b/src/RIGID/fix_rigid.cpp index dae552da37..1b696d45b9 100644 --- a/src/RIGID/fix_rigid.cpp +++ b/src/RIGID/fix_rigid.cpp @@ -2294,7 +2294,7 @@ double FixRigid::compute_scalar() for (int i = 0; i < nbody; i++) { t += masstotal[i] * (fflag[i][0]*vcm[i][0]*vcm[i][0] + - fflag[i][1]*vcm[i][1]*vcm[i][1] + \ + fflag[i][1]*vcm[i][1]*vcm[i][1] + fflag[i][2]*vcm[i][2]*vcm[i][2]); // wbody = angular velocity in body frame diff --git a/src/RIGID/fix_rigid_small.cpp b/src/RIGID/fix_rigid_small.cpp index 46c5459e57..ec9955c8e4 100644 --- a/src/RIGID/fix_rigid_small.cpp +++ b/src/RIGID/fix_rigid_small.cpp @@ -65,6 +65,9 @@ FixRigidSmall::FixRigidSmall(LAMMPS *lmp, int narg, char **arg) : { int i,ibody; + scalar_flag = 1; + extscalar = 0; + global_freq = 1; time_integrate = 1; rigid_flag = 1; virial_flag = 1; @@ -1982,25 +1985,6 @@ void FixRigidSmall::setup_bodies() memory->destroy(itensor); } -/* ---------------------------------------------------------------------- - memory usage of local atom-based arrays -------------------------------------------------------------------------- */ - -double FixRigidSmall::memory_usage() -{ - int nmax = atom->nmax; - double bytes = 2 * nmax * sizeof(int); - bytes += nmax*3 * sizeof(double); - bytes += maxvatom*6 * sizeof(double); // vatom - if (extended) { - bytes += nmax * sizeof(int); - if (orientflag) bytes = nmax*orientflag * sizeof(double); - if (dorientflag) bytes = nmax*3 * sizeof(double); - } - bytes += nmax_body * sizeof(Body); - return bytes; -} - /* ---------------------------------------------------------------------- allocate local atom-based arrays ------------------------------------------------------------------------- */ @@ -2525,6 +2509,67 @@ void FixRigidSmall::reset_dt() dtq = 0.5 * update->dt; } +/* ---------------------------------------------------------------------- + return temperature of collection of rigid bodies + non-active DOF are removed by fflag/tflag and in tfactor +------------------------------------------------------------------------- */ + +double FixRigidSmall::compute_scalar() +{ + double wbody[3],rot[3][3]; + + double *vcm,*inertia; + + double t = 0.0; + + for (int i = 0; i < nlocal_body; i++) { + vcm = body[i].vcm; + t += body[i].mass * (vcm[0]*vcm[0] + vcm[1]*vcm[1] + vcm[2]*vcm[2]); + + // for Iw^2 rotational term, need wbody = angular velocity in body frame + // not omega = angular velocity in space frame + + inertia = body[i].inertia; + MathExtra::quat_to_mat(body[i].quat,rot); + MathExtra::transpose_matvec(rot,body[i].angmom,wbody); + if (inertia[0] == 0.0) wbody[0] = 0.0; + else wbody[0] /= inertia[0]; + if (inertia[1] == 0.0) wbody[1] = 0.0; + else wbody[1] /= inertia[1]; + if (inertia[2] == 0.0) wbody[2] = 0.0; + else wbody[2] /= inertia[2]; + + t += inertia[0]*wbody[0]*wbody[0] + inertia[1]*wbody[1]*wbody[1] + + inertia[2]*wbody[2]*wbody[2]; + } + + double tall; + MPI_Allreduce(&t,&tall,1,MPI_DOUBLE,MPI_SUM,world); + + double tfactor = force->mvv2e / (6.0*nbody * force->boltz); + tall *= tfactor; + return tall; +} + +/* ---------------------------------------------------------------------- + memory usage of local atom-based arrays +------------------------------------------------------------------------- */ + +double FixRigidSmall::memory_usage() +{ + int nmax = atom->nmax; + double bytes = 2 * nmax * sizeof(int); + bytes += nmax*3 * sizeof(double); + bytes += maxvatom*6 * sizeof(double); // vatom + if (extended) { + bytes += nmax * sizeof(int); + if (orientflag) bytes = nmax*orientflag * sizeof(double); + if (dorientflag) bytes = nmax*3 * sizeof(double); + } + bytes += nmax_body * sizeof(Body); + return bytes; +} + /* ---------------------------------------------------------------------- debug method for sanity checking of atom/body data pointers ------------------------------------------------------------------------- */ diff --git a/src/RIGID/fix_rigid_small.h b/src/RIGID/fix_rigid_small.h index 2e5e9f0bc4..c6bd78b381 100644 --- a/src/RIGID/fix_rigid_small.h +++ b/src/RIGID/fix_rigid_small.h @@ -44,7 +44,6 @@ class FixRigidSmall : public Fix { void initial_integrate_respa(int, int, int); void final_integrate_respa(int, int); - double memory_usage(); void grow_arrays(int); void copy_arrays(int, int, int); void set_arrays(int); @@ -60,6 +59,8 @@ class FixRigidSmall : public Fix { int dof(int); void deform(int); void reset_dt(); + double compute_scalar(); + double memory_usage(); protected: int me,nprocs;