981 |
|
|
982 |
|
/* |
983 |
|
Return intertia tensor for entire system and angular momentum Vector. |
984 |
+ |
|
985 |
+ |
|
986 |
+ |
[ Ixx -Ixy -Ixz ] |
987 |
+ |
J =| -Iyx Iyy -Iyz | |
988 |
+ |
[ -Izx -Iyz Izz ] |
989 |
|
*/ |
990 |
|
|
991 |
|
void SimInfo::getInertiaTensor(Mat3x3d &inertiaTensor, Vector3d &angularMomentum){ |
1037 |
|
inertiaTensor(0,1) = -xy; |
1038 |
|
inertiaTensor(0,2) = -xz; |
1039 |
|
inertiaTensor(1,0) = -xy; |
1040 |
< |
inertiaTensor(2,0) = xx + zz; |
1040 |
> |
inertiaTensor(1,1) = xx + zz; |
1041 |
|
inertiaTensor(1,2) = -yz; |
1042 |
|
inertiaTensor(2,0) = -xz; |
1043 |
|
inertiaTensor(2,1) = -yz; |
1065 |
|
SimInfo::MoleculeIterator i; |
1066 |
|
Molecule* mol; |
1067 |
|
|
1068 |
< |
Vector3d thisq(0.0); |
1069 |
< |
Vector3d thisv(0.0); |
1068 |
> |
Vector3d thisr(0.0); |
1069 |
> |
Vector3d thisp(0.0); |
1070 |
|
|
1071 |
< |
double thisMass = 0.0; |
1071 |
> |
double thisMass; |
1072 |
|
|
1073 |
|
for (mol = beginMolecule(i); mol != NULL; mol = nextMolecule(i)) { |
1074 |
< |
thisq = mol->getCom()-com; |
1075 |
< |
thisv = mol->getComVel()-comVel; |
1076 |
< |
thisMass = mol->getMass(); |
1072 |
< |
angularMomentum += cross( thisq, thisv ) * thisMass; |
1074 |
> |
thisMass = mol->getMass(); |
1075 |
> |
thisr = mol->getCom()-com; |
1076 |
> |
thisp = (mol->getComVel()-comVel)*thisMass; |
1077 |
|
|
1078 |
+ |
angularMomentum += cross( thisr, thisp ); |
1079 |
+ |
|
1080 |
|
} |
1081 |
|
|
1082 |
|
#ifdef IS_MPI |