--- trunk/src/primitives/Molecule.cpp 2005/04/15 22:04:00 507 +++ branches/development/src/primitives/Molecule.cpp 2011/11/22 20:38:56 1665 @@ -6,19 +6,10 @@ * redistribute this software in source and binary code form, provided * that the following conditions are met: * - * 1. Acknowledgement of the program authors must be made in any - * publication of scientific results based in part on use of the - * program. An acceptable form of acknowledgement is citation of - * the article in which the program was described (Matthew - * A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher - * J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented - * Parallel Simulation Engine for Molecular Dynamics," - * J. Comput. Chem. 26, pp. 252-271 (2005)) - * - * 2. Redistributions of source code must retain the above copyright + * 1. Redistributions of source code must retain the above copyright * notice, this list of conditions and the following disclaimer. * - * 3. Redistributions in binary form must reproduce the above copyright + * 2. Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in the * documentation and/or other materials provided with the * distribution. @@ -37,6 +28,16 @@ * arising out of the use of or inability to use software, even if the * University of Notre Dame has been advised of the possibility of * such damages. + * + * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your + * research, please cite the appropriate papers when you publish your + * work. Good starting points are: + * + * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). + * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). + * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). + * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). + * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ /** @@ -53,78 +54,87 @@ #include "utils/MemoryUtils.hpp" #include "utils/simError.h" -namespace oopse { +namespace OpenMD { Molecule::Molecule(int stampId, int globalIndex, const std::string& molName) : stampId_(stampId), globalIndex_(globalIndex), moleculeName_(molName) { - - } - + } + 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_ don't own the objects integrableObjects_.clear(); } - + void Molecule::addAtom(Atom* atom) { if (std::find(atoms_.begin(), atoms_.end(), atom) == atoms_.end()) { atoms_.push_back(atom); } } - + void Molecule::addBond(Bond* bond) { if (std::find(bonds_.begin(), bonds_.end(), bond) == bonds_.end()) { bonds_.push_back(bond); } } - + void Molecule::addBend(Bend* bend) { if (std::find(bends_.begin(), bends_.end(), bend) == bends_.end()) { bends_.push_back(bend); } } - + void Molecule::addTorsion(Torsion* torsion) { - if (std::find(torsions_.begin(), torsions_.end(), torsion) == torsions_.end()) { + if (std::find(torsions_.begin(), torsions_.end(), torsion) == + torsions_.end()) { torsions_.push_back(torsion); } } + void Molecule::addInversion(Inversion* inversion) { + if (std::find(inversions_.begin(), inversions_.end(), inversion) == + inversions_.end()) { + inversions_.push_back(inversion); + } + } + void Molecule::addRigidBody(RigidBody *rb) { - if (std::find(rigidBodies_.begin(), rigidBodies_.end(), rb) == rigidBodies_.end()) { + if (std::find(rigidBodies_.begin(), rigidBodies_.end(), rb) == + rigidBodies_.end()) { rigidBodies_.push_back(rb); } } - + void Molecule::addCutoffGroup(CutoffGroup* cp) { - if (std::find(cutoffGroups_.begin(), cutoffGroups_.end(), cp) == cutoffGroups_.end()) { + 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()) { + 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()) { + if (std::find(constraintElems_.begin(), constraintElems_.end(), cp) == + constraintElems_.end()) { constraintElems_.push_back(cp); } - } - + void Molecule::complete() { std::set rigidAtoms; @@ -135,21 +145,25 @@ namespace oopse { for (rb = beginRigidBody(rbIter); rb != NULL; rb = nextRigidBody(rbIter)) { 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 + + // 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_)); @@ -166,32 +180,33 @@ namespace oopse { //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)){ + RealType mass = 0.0; + + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ mass += sd->getMass(); } - - return mass; - + + return mass; } 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)){ + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ mass = sd->getMass(); totalMass += mass; com += sd->getPos() * mass; } - + com /= totalMass; return com; @@ -201,40 +216,43 @@ namespace oopse { StuntDouble* sd; std::vector::iterator i; - for (sd = beginIntegrableObject(i); sd != NULL; sd = nextIntegrableObject(i)){ + for (sd = beginIntegrableObject(i); sd != NULL; sd = + nextIntegrableObject(i)){ sd->setPos(sd->getPos() + delta); - } - + } } 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)){ + 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(); @@ -244,14 +262,47 @@ namespace oopse { potential += bend->getPotential(); } - for (torsion = beginTorsion(torsionIter); torsion != NULL; torsion = nextTorsion(torsionIter)) { + 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; + + } + + void Molecule::addProperty(GenericData* genData) { + properties_.addProperty(genData); + } + void Molecule::removeProperty(const std::string& propName) { + properties_.removeProperty(propName); } + void Molecule::clearProperties() { + properties_.clearProperties(); + } + + std::vector Molecule::getPropertyNames() { + return properties_.getPropertyNames(); + } + + std::vector Molecule::getProperties() { + return properties_.getProperties(); + } + + GenericData* Molecule::getPropertyByName(const std::string& propName) { + return properties_.getPropertyByName(propName); + } + + + + std::ostream& operator <<(std::ostream& o, Molecule& mol) { o << std::endl; o << "Molecule " << mol.getGlobalIndex() << "has: " << std::endl; @@ -259,11 +310,12 @@ namespace oopse { 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 + +}//end namespace OpenMD