ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

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

Browse Source
This commit is contained in:
sjplimp
2016-09-07 16:12:51 +00:00
parent d3694613fd
commit ee4a1f0452
1 changed files with 2 additions and 2 deletions
Download Patch File Download Diff File Expand all files Collapse all files

4
src/group.cpp
View File

@ -1665,7 +1665,6 @@ void Group::inertia(int igroup, double *cm, double itensor[3][3], int iregion)
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
compute angular velocity omega from L and I compute angular velocity omega from L and I
diagonalize I instead of inverting it, to allow for a singular matrix
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
void Group::omega(double *angmom, double inertia[3][3], double *w) void Group::omega(double *angmom, double inertia[3][3], double *w)
@ -1684,7 +1683,7 @@ void Group::omega(double *angmom, double inertia[3][3], double *w)
// non-singular I matrix // non-singular I matrix
// use L = Iw, inverting I to solve for w // use L = Iw, inverting I to solve for w
// this should give exact zeroing of omega by velocity command // this should give exact zeroing of angular momentum by velocity command
if (determinant > EPSILON) { if (determinant > EPSILON) {
@ -1718,6 +1717,7 @@ void Group::omega(double *angmom, double inertia[3][3], double *w)
// handle (nearly) singular I matrix // handle (nearly) singular I matrix
// typically due to 2-atom group or linear molecule // typically due to 2-atom group or linear molecule
// use jacobi() and angmom_to_omega() to calculate valid omega // use jacobi() and angmom_to_omega() to calculate valid omega
// less exact answer than matrix inversion, due to iterative Jacobi method
} else { } else {
int ierror = MathExtra::jacobi(inertia,idiag,evectors); int ierror = MathExtra::jacobi(inertia,idiag,evectors);