--- trunk/src/primitives/Molecule.cpp 2005/02/21 15:28:56 372 +++ trunk/src/primitives/Molecule.cpp 2008/07/14 12:35:58 1277 @@ -1,4 +1,4 @@ - /* +/* * Copyright (c) 2005 The University of Notre Dame. All Rights Reserved. * * The University of Notre Dame grants you ("Licensee") a @@ -54,78 +54,87 @@ namespace oopse { #include "utils/simError.h" namespace oopse { -Molecule::Molecule(int stampId, int globalIndex, const std::string& molName) + Molecule::Molecule(int stampId, int globalIndex, const std::string& molName) : stampId_(stampId), globalIndex_(globalIndex), moleculeName_(molName) { - -} - -Molecule::~Molecule() { - - MemoryUtils::deleteVectorOfPointer(atoms_); - MemoryUtils::deleteVectorOfPointer(bonds_); - MemoryUtils::deleteVectorOfPointer(bends_); - MemoryUtils::deleteVectorOfPointer(torsions_); - MemoryUtils::deleteVectorOfPointer(rigidBodies_); - MemoryUtils::deleteVectorOfPointer(cutoffGroups_); - MemoryUtils::deleteVectorOfPointer(constraintPairs_); - MemoryUtils::deleteVectorOfPointer(constraintElems_); - //integrableObjects_ don't own the objects + } + + Molecule::~Molecule() { + + MemoryUtils::deletePointers(atoms_); + MemoryUtils::deletePointers(bonds_); + MemoryUtils::deletePointers(bends_); + MemoryUtils::deletePointers(torsions_); + MemoryUtils::deletePointers(inversions_); + MemoryUtils::deletePointers(rigidBodies_); + MemoryUtils::deletePointers(cutoffGroups_); + MemoryUtils::deletePointers(constraintPairs_); + MemoryUtils::deletePointers(constraintElems_); + // integrableObjects_ don't own the objects integrableObjects_.clear(); -} - -void Molecule::addAtom(Atom* atom) { + } + + void Molecule::addAtom(Atom* atom) { if (std::find(atoms_.begin(), atoms_.end(), atom) == atoms_.end()) { - atoms_.push_back(atom); + atoms_.push_back(atom); } -} - -void Molecule::addBond(Bond* bond) { + } + + void Molecule::addBond(Bond* bond) { if (std::find(bonds_.begin(), bonds_.end(), bond) == bonds_.end()) { - bonds_.push_back(bond); + bonds_.push_back(bond); } -} - -void Molecule::addBend(Bend* bend) { + } + + void Molecule::addBend(Bend* bend) { if (std::find(bends_.begin(), bends_.end(), bend) == bends_.end()) { - bends_.push_back(bend); + bends_.push_back(bend); } -} - -void Molecule::addTorsion(Torsion* torsion) { - if (std::find(torsions_.begin(), torsions_.end(), torsion) == torsions_.end()) { - torsions_.push_back(torsion); + } + + void Molecule::addTorsion(Torsion* torsion) { + if (std::find(torsions_.begin(), torsions_.end(), torsion) == + torsions_.end()) { + torsions_.push_back(torsion); } -} + } -void Molecule::addRigidBody(RigidBody *rb) { - if (std::find(rigidBodies_.begin(), rigidBodies_.end(), rb) == rigidBodies_.end()) { - rigidBodies_.push_back(rb); + void Molecule::addInversion(Inversion* inversion) { + if (std::find(inversions_.begin(), inversions_.end(), inversion) == + inversions_.end()) { + inversions_.push_back(inversion); } -} - -void Molecule::addCutoffGroup(CutoffGroup* cp) { - if (std::find(cutoffGroups_.begin(), cutoffGroups_.end(), cp) == cutoffGroups_.end()) { - cutoffGroups_.push_back(cp); + } + + void Molecule::addRigidBody(RigidBody *rb) { + if (std::find(rigidBodies_.begin(), rigidBodies_.end(), rb) == + rigidBodies_.end()) { + rigidBodies_.push_back(rb); } - -} - -void Molecule::addConstraintPair(ConstraintPair* cp) { - if (std::find(constraintPairs_.begin(), constraintPairs_.end(), cp) == constraintPairs_.end()) { - constraintPairs_.push_back(cp); + } + + void Molecule::addCutoffGroup(CutoffGroup* cp) { + if (std::find(cutoffGroups_.begin(), cutoffGroups_.end(), cp) == + cutoffGroups_.end()) { + cutoffGroups_.push_back(cp); + } + } + + void Molecule::addConstraintPair(ConstraintPair* cp) { + if (std::find(constraintPairs_.begin(), constraintPairs_.end(), cp) == + constraintPairs_.end()) { + constraintPairs_.push_back(cp); + } + } + + void Molecule::addConstraintElem(ConstraintElem* cp) { + if (std::find(constraintElems_.begin(), constraintElems_.end(), cp) == + constraintElems_.end()) { + constraintElems_.push_back(cp); } - -} - -void Molecule::addConstraintElem(ConstraintElem* cp) { - if (std::find(constraintElems_.begin(), constraintElems_.end(), cp) == constraintElems_.end()) { - constraintElems_.push_back(cp); - } - -} - -void Molecule::complete() { + } + + void Molecule::complete() { std::set rigidAtoms; RigidBody* rb; @@ -133,23 +142,27 @@ void Molecule::complete() { for (rb = beginRigidBody(rbIter); rb != NULL; rb = nextRigidBody(rbIter)) { - rigidAtoms.insert(rb->getBeginAtomIter(), rb->getEndAtomIter()); + rigidAtoms.insert(rb->getBeginAtomIter(), rb->getEndAtomIter()); } - + Atom* atom; AtomIterator ai; for (atom = beginAtom(ai); atom != NULL; atom = nextAtom(ai)) { - - if (rigidAtoms.find(*ai) == rigidAtoms.end()) { - //if an atom does not belong to a rigid body, it is an integrable object - integrableObjects_.push_back(*ai); - } - } + + if (rigidAtoms.find(*ai) == rigidAtoms.end()) { + + // If an atom does not belong to a rigid body, it is an + // integrable object + integrableObjects_.push_back(*ai); + } + } + //find all free atoms (which do not belong to rigid bodies) - //performs the "difference" operation from set theory, the output range contains a copy of every - //element that is contained in [allAtoms.begin(), allAtoms.end()) and not contained in - //[rigidAtoms.begin(), rigidAtoms.end()). + // performs the "difference" operation from set theory, the output + // range contains a copy of every element that is contained in + // [allAtoms.begin(), allAtoms.end()) and not contained in + // [rigidAtoms.begin(), rigidAtoms.end()). //std::set_difference(allAtoms.begin(), allAtoms.end(), rigidAtoms.begin(), rigidAtoms.end(), // std::back_inserter(integrableObjects_)); @@ -164,106 +177,117 @@ void Molecule::complete() { integrableObjects_.push_back(rb); } //integrableObjects_.insert(integrableObjects_.end(), rigidBodies_.begin(), rigidBodies_.end()); -} + } -double Molecule::getMass() { + RealType Molecule::getMass() { StuntDouble* sd; std::vector::iterator i; - double mass = 0.0; - - for (sd = beginIntegrableObject(i); sd != NULL; sd = nextIntegrableObject(i)){ - mass += sd->getMass(); + RealType mass = 0.0; + + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ + mass += sd->getMass(); } + + return mass; + } - return mass; - -} - -Vector3d Molecule::getCom() { + Vector3d Molecule::getCom() { StuntDouble* sd; std::vector::iterator i; Vector3d com; - double totalMass = 0; - double mass; + RealType totalMass = 0; + RealType mass; - for (sd = beginIntegrableObject(i); sd != NULL; sd = nextIntegrableObject(i)){ - mass = sd->getMass(); - totalMass += mass; - com += sd->getPos() * mass; + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ + mass = sd->getMass(); + totalMass += mass; + com += sd->getPos() * mass; } - + com /= totalMass; return com; -} + } -void Molecule::moveCom(const Vector3d& delta) { + void Molecule::moveCom(const Vector3d& delta) { StuntDouble* sd; std::vector::iterator i; - for (sd = beginIntegrableObject(i); sd != NULL; sd = nextIntegrableObject(i)){ - sd->setPos(sd->getPos() + delta); - } + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ + sd->setPos(sd->getPos() + delta); + } + } -} - -Vector3d Molecule::getComVel() { + Vector3d Molecule::getComVel() { StuntDouble* sd; std::vector::iterator i; Vector3d velCom; - double totalMass = 0; - double mass; + RealType totalMass = 0; + RealType mass; - for (sd = beginIntegrableObject(i); sd != NULL; sd = nextIntegrableObject(i)){ - mass = sd->getMass(); - totalMass += mass; - velCom += sd->getVel() * mass; + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ + mass = sd->getMass(); + totalMass += mass; + velCom += sd->getVel() * mass; } - + velCom /= totalMass; - + return velCom; -} + } -double Molecule::getPotential() { + RealType Molecule::getPotential() { Bond* bond; Bend* bend; Torsion* torsion; + Inversion* inversion; Molecule::BondIterator bondIter;; Molecule::BendIterator bendIter; Molecule::TorsionIterator torsionIter; + Molecule::InversionIterator inversionIter; - double potential = 0.0; + RealType potential = 0.0; for (bond = beginBond(bondIter); bond != NULL; bond = nextBond(bondIter)) { - potential += bond->getPotential(); + potential += bond->getPotential(); } for (bend = beginBend(bendIter); bend != NULL; bend = nextBend(bendIter)) { - potential += bend->getPotential(); + potential += bend->getPotential(); } - for (torsion = beginTorsion(torsionIter); torsion != NULL; torsion = nextTorsion(torsionIter)) { - potential += torsion->getPotential(); + for (torsion = beginTorsion(torsionIter); torsion != NULL; torsion = + nextTorsion(torsionIter)) { + potential += torsion->getPotential(); } - + + for (inversion = beginInversion(inversionIter); torsion != NULL; + inversion = nextInversion(inversionIter)) { + potential += inversion->getPotential(); + } + return potential; - -} - -std::ostream& operator <<(std::ostream& o, Molecule& mol) { + + } + + std::ostream& operator <<(std::ostream& o, Molecule& mol) { o << std::endl; o << "Molecule " << mol.getGlobalIndex() << "has: " << std::endl; o << mol.getNAtoms() << " atoms" << std::endl; o << mol.getNBonds() << " bonds" << std::endl; o << mol.getNBends() << " bends" << std::endl; o << mol.getNTorsions() << " torsions" << std::endl; + o << mol.getNInversions() << " inversions" << std::endl; o << mol.getNRigidBodies() << " rigid bodies" << std::endl; o << mol.getNIntegrableObjects() << "integrable objects" << std::endl; o << mol.getNCutoffGroups() << "cutoff groups" << std::endl; o << mol.getNConstraintPairs() << "constraint pairs" << std::endl; return o; -} - + } + }//end namespace oopse