45#include "brains/Thermo.hpp"
55#include "types/FixedChargeAdapter.hpp"
56#include "types/FluctuatingChargeAdapter.hpp"
57#include "types/MultipoleAdapter.hpp"
58#include "utils/Constants.hpp"
59#include "utils/simError.h"
61#include "math/AlphaHull.hpp"
62#include "math/ConvexHull.hpp"
68 RealType Thermo::getTranslationalKinetic() {
71 if (!snap->hasTranslationalKineticEnergy) {
72 SimInfo::MoleculeIterator miter;
73 vector<StuntDouble*>::iterator iiter;
78 RealType kinetic(0.0);
82 for (sd = mol->beginIntegrableObject(iiter); sd != NULL;
83 sd = mol->nextIntegrableObject(iiter)) {
88 mass * (vel[0] * vel[0] + vel[1] * vel[1] + vel[2] * vel[2]);
93 MPI_Allreduce(MPI_IN_PLACE, &kinetic, 1, MPI_REALTYPE, MPI_SUM,
97 kinetic = kinetic * 0.5 / Constants::energyConvert;
99 snap->setTranslationalKineticEnergy(kinetic);
101 return snap->getTranslationalKineticEnergy();
104 RealType Thermo::getRotationalKinetic() {
107 if (!snap->hasRotationalKineticEnergy) {
108 SimInfo::MoleculeIterator miter;
109 vector<StuntDouble*>::iterator iiter;
115 RealType kinetic(0.0);
119 for (sd = mol->beginIntegrableObject(iiter); sd != NULL;
120 sd = mol->nextIntegrableObject(iiter)) {
129 kinetic += angMom[j] * angMom[j] / I(j, j) +
130 angMom[k] * angMom[k] / I(k, k);
132 kinetic += angMom[0] * angMom[0] / I(0, 0) +
133 angMom[1] * angMom[1] / I(1, 1) +
134 angMom[2] * angMom[2] / I(2, 2);
141 MPI_Allreduce(MPI_IN_PLACE, &kinetic, 1, MPI_REALTYPE, MPI_SUM,
145 kinetic = kinetic * 0.5 / Constants::energyConvert;
147 snap->setRotationalKineticEnergy(kinetic);
149 return snap->getRotationalKineticEnergy();
152 RealType Thermo::getKinetic() {
155 if (!snap->hasKineticEnergy) {
156 RealType ke = getTranslationalKinetic() + getRotationalKinetic() +
157 getElectronicKinetic();
159 snap->setKineticEnergy(ke);
161 return snap->getKineticEnergy();
164 RealType Thermo::getPotential() {
169 return snap->getPotentialEnergy();
172 potVec Thermo::getSelectionPotentials() {
177 return snap->getSelectionPotentials();
180 RealType Thermo::getTotalEnergy() {
183 if (!snap->hasTotalEnergy) {
184 snap->setTotalEnergy(this->getKinetic() + this->getPotential());
187 return snap->getTotalEnergy();
193 RealType Thermo::getTemperature() {
196 if (!snap->hasTemperature) {
198 this->getTranslationalKinetic() + this->getRotationalKinetic();
200 RealType temperature =
201 (2.0 * nuclearKE) / (info_->
getNdf() * Constants::kb);
203 snap->setTemperature(temperature);
206 return snap->getTemperature();
209 RealType Thermo::getElectronicKinetic() {
212 if (!snap->hasElectronicKineticEnergy) {
213 SimInfo::MoleculeIterator miter;
214 vector<Atom*>::iterator iiter;
219 RealType kinetic(0.0);
223 for (atom = mol->beginFluctuatingCharge(iiter); atom != NULL;
224 atom = mol->nextFluctuatingCharge(iiter)) {
225 cmass = atom->getChargeMass();
226 cvel = atom->getFlucQVel();
227 kinetic += cmass * cvel * cvel;
232 MPI_Allreduce(MPI_IN_PLACE, &kinetic, 1, MPI_REALTYPE, MPI_SUM,
237 snap->setElectronicKineticEnergy(kinetic);
240 return snap->getElectronicKineticEnergy();
243 RealType Thermo::getElectronicTemperature() {
246 if (!snap->hasElectronicTemperature) {
247 RealType eTemp = (2.0 * this->getElectronicKinetic()) /
250 snap->setElectronicTemperature(eTemp);
253 return snap->getElectronicTemperature();
256 RealType Thermo::getNetCharge() {
259 if (!snap->hasNetCharge) {
260 SimInfo::MoleculeIterator miter;
261 vector<Atom*>::iterator aiter;
265 RealType netCharge(0.0);
269 for (atom = mol->beginAtom(aiter); atom != NULL;
270 atom = mol->nextAtom(aiter)) {
274 if (fca.isFixedCharge()) { charge = fca.getCharge(); }
278 if (fqa.isFluctuatingCharge()) { charge += atom->getFlucQPos(); }
285 MPI_Allreduce(MPI_IN_PLACE, &netCharge, 1, MPI_REALTYPE, MPI_SUM,
289 snap->setNetCharge(netCharge);
292 return snap->getNetCharge();
295 RealType Thermo::getChargeMomentum() {
298 if (!snap->hasChargeMomentum) {
299 SimInfo::MoleculeIterator miter;
300 vector<Atom*>::iterator iiter;
305 RealType momentum(0.0);
309 for (atom = mol->beginFluctuatingCharge(iiter); atom != NULL;
310 atom = mol->nextFluctuatingCharge(iiter)) {
311 cmass = atom->getChargeMass();
312 cvel = atom->getFlucQVel();
314 momentum += cmass * cvel;
319 MPI_Allreduce(MPI_IN_PLACE, &momentum, 1, MPI_REALTYPE, MPI_SUM,
323 snap->setChargeMomentum(momentum);
326 return snap->getChargeMomentum();
329 std::vector<Vector3d> Thermo::getCurrentDensity() {
333 SimInfo::MoleculeIterator miter;
334 std::vector<Atom*>::iterator iiter;
335 std::vector<RigidBody*>::iterator ri;
340 AtomTypeSet::iterator at;
342 std::vector<Vector3d> typeJc(simTypes.size(), V3Zero);
347 for (rb = mol->beginRigidBody(ri); rb != NULL;
348 rb = mol->nextRigidBody(ri)) {
352 for (atom = mol->beginAtom(iiter); atom != NULL;
353 atom = mol->nextAtom(iiter)) {
358 atype = atom->getAtomType();
360 if (fca.isFixedCharge()) q = fca.getCharge();
362 if (fqa.isFluctuatingCharge()) q += atom->getFlucQPos();
365 at = std::find(simTypes.begin(), simTypes.end(), atype);
366 if (at != simTypes.end()) {
367 typeIndex = std::distance(simTypes.begin(), at);
370 if (typeIndex != -1) { typeJc[typeIndex] += q * v; }
376 MPI_Allreduce(MPI_IN_PLACE, &Jc[0], 3, MPI_REALTYPE, MPI_SUM,
378 for (
unsigned int j = 0; j < simTypes.size(); j++) {
379 MPI_Allreduce(MPI_IN_PLACE, &typeJc[j][0], 3, MPI_REALTYPE, MPI_SUM,
384 RealType vol = snap->getVolume();
386 Jc /= (vol * Constants::currentDensityConvert);
387 for (
unsigned int j = 0; j < simTypes.size(); j++) {
388 typeJc[j] /= (vol * Constants::currentDensityConvert);
391 std::vector<Vector3d> result;
394 result.push_back(Jc);
395 for (
unsigned int j = 0; j < simTypes.size(); j++)
396 result.push_back(typeJc[j]);
401 RealType Thermo::getVolume() {
403 return snap->getVolume();
406 RealType Thermo::getPressure() {
409 if (!snap->hasPressure) {
417 pressure = Constants::pressureConvert *
418 (tensor(0, 0) + tensor(1, 1) + tensor(2, 2)) / 3.0;
420 snap->setPressure(pressure);
423 return snap->getPressure();
426 RealType Thermo::getPressure(
Snapshot* snap) {
427 if (!snap->hasPressure) {
434 pressure = Constants::pressureConvert *
435 (tensor(0, 0) + tensor(1, 1) + tensor(2, 2)) / 3.0;
437 snap->setPressure(pressure);
440 return snap->getPressure();
449 if (!snap->hasPressureTensor) {
455 SimInfo::MoleculeIterator i;
456 vector<StuntDouble*>::iterator j;
461 for (sd = mol->beginIntegrableObject(j); sd != NULL;
462 sd = mol->nextIntegrableObject(j)) {
465 p_tens += mass * outProduct(vcom, vcom);
470 MPI_Allreduce(MPI_IN_PLACE, p_tens.getArrayPointer(), 9, MPI_REALTYPE,
471 MPI_SUM, MPI_COMM_WORLD);
474 RealType volume = this->getVolume();
475 Mat3x3d virialTensor = snap->getVirialTensor();
478 (p_tens + Constants::energyConvert * virialTensor) / volume;
480 snap->setPressureTensor(pressureTensor);
482 return snap->getPressureTensor();
489 if (!snap->hasPressureTensor) {
495 SimInfo::MoleculeIterator i;
496 vector<StuntDouble*>::iterator j;
501 for (sd = mol->beginIntegrableObject(j); sd != NULL;
502 sd = mol->nextIntegrableObject(j)) {
505 p_tens += mass * outProduct(vcom, vcom);
510 MPI_Allreduce(MPI_IN_PLACE, p_tens.getArrayPointer(), 9, MPI_REALTYPE,
511 MPI_SUM, MPI_COMM_WORLD);
514 RealType volume = this->getVolume();
515 Mat3x3d virialTensor = snap->getVirialTensor();
518 (p_tens + Constants::energyConvert * virialTensor) / volume;
520 snap->setPressureTensor(pressureTensor);
522 return snap->getPressureTensor();
528 if (!snap->hasSystemDipole) {
529 SimInfo::MoleculeIterator miter;
530 vector<Atom*>::iterator aiter;
543 RealType chargeToC = 1.60217733e-19;
544 RealType angstromToM = 1.0e-10;
545 RealType debyeToCm = 3.33564095198e-30;
549 for (atom = mol->beginAtom(aiter); atom != NULL;
550 atom = mol->nextAtom(aiter)) {
554 if (fca.isFixedCharge()) { charge = fca.getCharge(); }
558 if (fqa.isFluctuatingCharge()) { charge += atom->getFlucQPos(); }
570 }
else if (charge > 0.0) {
576 if (atom->isDipole()) {
577 dipoleVector += atom->getDipole() * debyeToCm;
583 MPI_Allreduce(MPI_IN_PLACE, &pChg, 1, MPI_REALTYPE, MPI_SUM,
585 MPI_Allreduce(MPI_IN_PLACE, &nChg, 1, MPI_REALTYPE, MPI_SUM,
588 MPI_Allreduce(MPI_IN_PLACE, &pCount, 1, MPI_INTEGER, MPI_SUM,
590 MPI_Allreduce(MPI_IN_PLACE, &nCount, 1, MPI_INTEGER, MPI_SUM,
593 MPI_Allreduce(MPI_IN_PLACE, pPos.getArrayPointer(), 3, MPI_REALTYPE,
594 MPI_SUM, MPI_COMM_WORLD);
595 MPI_Allreduce(MPI_IN_PLACE, nPos.getArrayPointer(), 3, MPI_REALTYPE,
596 MPI_SUM, MPI_COMM_WORLD);
598 MPI_Allreduce(MPI_IN_PLACE, dipoleVector.getArrayPointer(), 3,
599 MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
606 RealType chg_value = nChg <= pChg ? nChg : pChg;
609 if (pCount > 0 && nCount > 0) {
615 boxDipole += (pPos - nPos) * chg_value;
616 snap->setSystemDipole(boxDipole);
619 return snap->getSystemDipole();
625 if (!snap->hasSystemQuadrupole) {
626 SimInfo::MoleculeIterator miter;
627 vector<Atom*>::iterator aiter;
635 RealType chargeToC = 1.60217733e-19;
636 RealType angstromToM = 1.0e-10;
637 RealType debyeToCm = 3.33564095198e-30;
641 for (atom = mol->beginAtom(aiter); atom != NULL;
642 atom = mol->nextAtom(aiter)) {
650 if (fca.isFixedCharge()) { charge = fca.getCharge(); }
654 if (fqa.isFluctuatingCharge()) { charge += atom->getFlucQPos(); }
658 qpole += 0.5 * charge * outProduct(ri, ri);
663 dipole = atom->getDipole() * debyeToCm;
664 qpole += 0.5 * outProduct(dipole, ri);
665 qpole += 0.5 * outProduct(ri, dipole);
668 if (ma.isQuadrupole()) {
669 qpole += atom->getQuadrupole() * debyeToCm * angstromToM;
675 MPI_Allreduce(MPI_IN_PLACE, qpole.getArrayPointer(), 9, MPI_REALTYPE,
676 MPI_SUM, MPI_COMM_WORLD);
679 snap->setSystemQuadrupole(qpole);
682 return snap->getSystemQuadrupole();
688 SimInfo::MoleculeIterator miter;
689 vector<StuntDouble*>::iterator iiter;
692 RigidBody::AtomIterator ai;
712 for (sd = mol->beginIntegrableObject(iiter); sd != NULL;
713 sd = mol->nextIntegrableObject(iiter)) {
717 kinetic = mass * (vel[0] * vel[0] + vel[1] * vel[1] + vel[2] * vel[2]);
727 kinetic += angMom[j] * angMom[j] / I(j, j) +
728 angMom[k] * angMom[k] / I(k, k);
730 kinetic += angMom[0] * angMom[0] / I(0, 0) +
731 angMom[1] * angMom[1] / I(1, 1) +
732 angMom[2] * angMom[2] / I(2, 2);
740 for (atom = rb->beginAtom(ai); atom != NULL;
741 atom = rb->nextAtom(ai)) {
742 potential += atom->getParticlePot();
748 potential *= Constants::energyConvert;
750 eatom = (kinetic + potential) / 2.0;
751 heatFluxJc[0] += eatom * vel[0];
752 heatFluxJc[1] += eatom * vel[1];
753 heatFluxJc[2] += eatom * vel[2];
762 MPI_Allreduce(MPI_IN_PLACE, &heatFluxJc[0], 3, MPI_REALTYPE, MPI_SUM,
769 currSnapshot->getConductiveHeatFlux() * Constants::energyConvert;
772 return (heatFluxJv + heatFluxJc) / this->getVolume();
778 if (!snap->hasCOMvel) {
779 SimInfo::MoleculeIterator i;
783 RealType totalMass(0.0);
787 RealType mass = mol->
getMass();
793 MPI_Allreduce(MPI_IN_PLACE, &totalMass, 1, MPI_REALTYPE, MPI_SUM,
795 MPI_Allreduce(MPI_IN_PLACE, comVel.getArrayPointer(), 3, MPI_REALTYPE,
796 MPI_SUM, MPI_COMM_WORLD);
800 snap->setCOMvel(comVel);
802 return snap->getCOMvel();
809 SimInfo::MoleculeIterator i;
813 RealType totalMass(0.0);
817 RealType mass = mol->
getMass();
819 com += mass * mol->
getCom();
823 MPI_Allreduce(MPI_IN_PLACE, &totalMass, 1, MPI_REALTYPE, MPI_SUM,
825 MPI_Allreduce(MPI_IN_PLACE, com.getArrayPointer(), 3, MPI_REALTYPE,
826 MPI_SUM, MPI_COMM_WORLD);
832 return snap->getCOM();
842 if (!(snap->hasCOM && snap->hasCOMvel)) {
843 SimInfo::MoleculeIterator i;
846 RealType totalMass(0.0);
853 RealType mass = mol->
getMass();
855 com += mass * mol->
getCom();
860 MPI_Allreduce(MPI_IN_PLACE, &totalMass, 1, MPI_REALTYPE, MPI_SUM,
862 MPI_Allreduce(MPI_IN_PLACE, com.getArrayPointer(), 3, MPI_REALTYPE,
863 MPI_SUM, MPI_COMM_WORLD);
864 MPI_Allreduce(MPI_IN_PLACE, comVel.getArrayPointer(), 3, MPI_REALTYPE,
865 MPI_SUM, MPI_COMM_WORLD);
871 snap->setCOMvel(comVel);
873 com = snap->getCOM();
874 comVel = snap->getCOMvel();
891 Molecule::IntegrableObjectIterator j;
894 if (!(snap->hasInertiaTensor && snap->hasCOMw)) {
906 SimInfo::MoleculeIterator i;
916 for (sd = mol->beginIntegrableObject(j); sd != NULL;
917 sd = mol->nextIntegrableObject(j)) {
919 v = sd->
getVel() - comVel;
924 xx += r[0] * r[0] * m;
925 yy += r[1] * r[1] * m;
926 zz += r[2] * r[2] * m;
929 xy += r[0] * r[1] * m;
930 xz += r[0] * r[2] * m;
931 yz += r[1] * r[2] * m;
933 angularMomentum +=
cross(r, v) * m;
937 inertiaTensor(0, 0) = yy + zz;
938 inertiaTensor(0, 1) = -xy;
939 inertiaTensor(0, 2) = -xz;
940 inertiaTensor(1, 0) = -xy;
941 inertiaTensor(1, 1) = xx + zz;
942 inertiaTensor(1, 2) = -yz;
943 inertiaTensor(2, 0) = -xz;
944 inertiaTensor(2, 1) = -yz;
945 inertiaTensor(2, 2) = xx + yy;
949 MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
951 MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
954 snap->setCOMw(angularMomentum);
955 snap->setInertiaTensor(inertiaTensor);
958 angularMomentum = snap->getCOMw();
959 inertiaTensor = snap->getInertiaTensor();
967 if (!(snap->hasBoundingBox)) {
968 SimInfo::MoleculeIterator i;
969 Molecule::RigidBodyIterator ri;
970 Molecule::AtomIterator ai;
980 for (rb = mol->beginRigidBody(ri); rb != NULL;
981 rb = mol->nextRigidBody(ri)) {
985 for (atom = mol->beginAtom(ai); atom != NULL;
986 atom = mol->nextAtom(ai)) {
987 pos = atom->getPos();
993 for (
int i = 0; i < 3; i++) {
994 bMax[i] = max(bMax[i], pos[i]);
995 bMin[i] = min(bMin[i], pos[i]);
1003 MPI_Allreduce(MPI_IN_PLACE, &bMax[0], 3, MPI_REALTYPE, MPI_MAX,
1006 MPI_Allreduce(MPI_IN_PLACE, &bMin[0], 3, MPI_REALTYPE, MPI_MIN,
1010 for (
int i = 0; i < 3; i++) {
1011 bBox(i, i) = bMax[i] - bMin[i];
1023 if (!snap->hasCOMw) {
1030 SimInfo::MoleculeIterator i;
1041 thisr = mol->
getCom() - com;
1042 thisp = (mol->
getComVel() - comVel) * thisMass;
1044 angularMomentum +=
cross(thisr, thisp);
1049 MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
1052 snap->setCOMw(angularMomentum);
1055 return snap->getCOMw();
1069 if (!snap->hasGyrationalVolume) {
1073 RealType sysconstants;
1077 geomCnst = 3.0 / 2.0;
1081 det = intTensor.determinant();
1085 4.0 / 3.0 * Constants::PI * pow(sysconstants, geomCnst) * sqrt(det);
1087 snap->setGyrationalVolume(volume);
1089 return snap->getGyrationalVolume();
1095 if (!(snap->hasInertiaTensor && snap->hasGyrationalVolume)) {
1098 RealType sysconstants;
1101 geomCnst = 3.0 / 2.0;
1105 detI = intTensor.determinant();
1109 4.0 / 3.0 * Constants::PI * pow(sysconstants, geomCnst) * sqrt(detI);
1110 snap->setGyrationalVolume(volume);
1112 volume = snap->getGyrationalVolume();
1118 RealType Thermo::getTaggedAtomPairDistance() {
1120 Globals* simParams = info_->getSimParams();
1122 if (simParams->haveTaggedAtomPair() &&
1123 simParams->havePrintTaggedPairDistance()) {
1124 if (simParams->getPrintTaggedPairDistance()) {
1125 pair<int, int> tap = simParams->getTaggedAtomPair();
1129 int mol1 = info_->getGlobalMolMembership(tap.first);
1130 int mol2 = info_->getGlobalMolMembership(tap.second);
1136 if (proc1 == worldRank) {
1142 MPI_Bcast(data, 3, MPI_REALTYPE, proc1, MPI_COMM_WORLD);
1144 MPI_Bcast(data, 3, MPI_REALTYPE, proc1, MPI_COMM_WORLD);
1148 if (proc2 == worldRank) {
1154 MPI_Bcast(data, 3, MPI_REALTYPE, proc2, MPI_COMM_WORLD);
1156 MPI_Bcast(data, 3, MPI_REALTYPE, proc2, MPI_COMM_WORLD);
1163 pos2 = at2->getPos();
1166 currSnapshot->wrapVector(rab);
1174 RealType Thermo::getHullVolume() {
1177 if (!snap->hasHullVolume) {
1180 Globals* simParams = info_->getSimParams();
1181 const std::string ht = simParams->getHULL_Method();
1183 if (ht ==
"Convex") {
1185 }
else if (ht ==
"AlphaShape") {
1186 surfaceMesh_ =
new AlphaHull(simParams->getAlpha());
1193 std::vector<StuntDouble*> localSites_;
1196 SimInfo::MoleculeIterator i;
1197 Molecule::IntegrableObjectIterator j;
1201 for (sd = mol->beginIntegrableObject(j); sd != NULL;
1202 sd = mol->nextIntegrableObject(j)) {
1203 localSites_.push_back(sd);
1208 surfaceMesh_->computeHull(localSites_);
1209 snap->setHullVolume(surfaceMesh_->getVolume());
1211 delete surfaceMesh_;
1214 return snap->getHullVolume();
Compute alpha complex or alpha shape.
AtomType is what OpenMD looks to for unchanging data about an atom.
RealType getMass()
get total mass of this molecule
Vector3d getComVel()
Returns the velocity of center of mass of this molecule.
Vector3d getCom()
Returns the current center of mass position of this molecule.
Real * getArrayPointer()
Returns the pointer of internal array.
int getNdf()
Returns the number of degrees of freedom.
int getNGlobalIntegrableObjects()
Returns the total number of integrable objects (total number of rigid bodies plus the total number of...
int getMolToProc(int globalIndex)
Finds the processor where a molecule resides.
Molecule * beginMolecule(MoleculeIterator &i)
Returns the first molecule in this SimInfo and intialize the iterator.
int getNFluctuatingCharges()
Returns the total number of fluctuating charges that are present.
Molecule * nextMolecule(MoleculeIterator &i)
Returns the next avaliable Molecule based on the iterator.
SnapshotManager * getSnapshotManager()
Returns the snapshot manager.
StuntDouble * getIOIndexToIntegrableObject(int index)
return an integral objects by its global index.
AtomTypeSet getSimulatedAtomTypes()
Returns the set of atom types present in this simulation.
The Snapshot class is a repository storing dynamic data during a Simulation.
Mat3x3d getBoundingBox()
Returns the Bounding Box.
void setBoundingBox(const Mat3x3d &m)
Sets the Bounding Box.
void wrapVector(Vector3d &v)
Wrapping the vector according to periodic boundary condition.
Snapshot * getCurrentSnapshot()
Returns the pointer of current snapshot.
Real determinant() const
Returns the determinant of this matrix.
"Don't move, or you're dead! Stand up! Captain, we've got them!"
Vector3d getVel()
Returns the current velocity of this stuntDouble.
int linearAxis()
Returns the linear axis of the rigidbody, atom and directional atom will always return -1.
RealType getMass()
Returns the mass of this stuntDouble.
virtual Mat3x3d getI()=0
Returns the inertia tensor of this stuntDouble.
bool isLinear()
Tests the if this stuntDouble is a linear rigidbody.
Vector3d getPos()
Returns the current position of this stuntDouble.
RealType getParticlePot()
Returns the current particlePot of this stuntDouble.
bool isRigidBody()
Tests if this stuntDouble is a rigid body.
Vector3d getJ()
Returns the current angular momentum of this stuntDouble (body -fixed).
bool isDirectional()
Tests if this stuntDouble is a directional one.
Mat3x3d getPressureTensor()
gives the pressure tensor in amu*fs^-2*Ang^-1
Vector3d getSystemDipole()
accumulate and return the simulation box dipole moment in C*m
void getComAll(Vector3d &com, Vector3d &comVel)
Returns the center of the mass and Center of Mass velocity of the whole system.
void getInertiaTensor(Mat3x3d &inertiaTensor, Vector3d &angularMomentum)
Returns the inertia tensor and the total angular momentum for for the entire system.
RealType getGyrationalVolume()
Returns volume of system as estimated by an ellipsoid defined by the radii of gyration.
Vector3d getCom()
Returns the center of the mass of the whole system.
Mat3x3d getSystemQuadrupole()
accumulate and return the simulation box dipole moment in debye Angstroms
Vector3d getAngularMomentum()
Returns system angular momentum.
Vector3d getComVel()
Returns the velocity of center of mass of the whole system.
Mat3x3d getBoundingBox()
Returns the Axis-aligned bounding box for the current system.
Real length()
Returns the length of this vector.
Real * getArrayPointer()
Returns the pointer of internal array.
This basic Periodic Table class was originally taken from the data.cpp file in OpenBabel.
Vector3< Real > cross(const Vector3< Real > &v1, const Vector3< Real > &v2)
Returns the cross product of two Vectors.