--- trunk/src/selection/NameFinder.cpp 2005/04/05 23:09:48 452 +++ trunk/src/selection/NameFinder.cpp 2014/01/15 22:26:18 1962 @@ -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,231 +28,401 @@ * 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, 234107 (2008). + * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). + * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ #include "selection/NameFinder.hpp" #include "utils/wildcards.hpp" #include "utils/StringTokenizer.hpp" #include "primitives/Molecule.hpp" #include "utils/StringUtils.hpp" -namespace oopse { +namespace OpenMD { -TreeNode::~TreeNode(){ + TreeNode::~TreeNode(){ std::map::iterator i; for ( i = children.begin(); i != children.end(); ++i) { - i->second->~TreeNode(); + i->second->~TreeNode(); } children.clear(); -} + } - -NameFinder::NameFinder(SimInfo* info) : info_(info), root_(NULL){ - nStuntDouble_ = info_->getNGlobalAtoms() + info_->getNGlobalRigidBodies(); + NameFinder::NameFinder(SimInfo* info) : info_(info), root_(NULL){ + nObjects_.push_back(info_->getNGlobalAtoms()+info_->getNGlobalRigidBodies()); + nObjects_.push_back(info_->getNGlobalBonds()); + nObjects_.push_back(info_->getNGlobalBends()); + nObjects_.push_back(info_->getNGlobalTorsions()); + nObjects_.push_back(info_->getNGlobalInversions()); loadNames(); -} + } - -NameFinder::~NameFinder(){ + NameFinder::~NameFinder(){ delete root_; -} + } -void NameFinder::loadNames() { - - std::map::iterator foundIter; + void NameFinder::loadNames() { SimInfo::MoleculeIterator mi; - Molecule* mol; Molecule::AtomIterator ai; - Atom* atom; Molecule::RigidBodyIterator rbIter; + Molecule::BondIterator bondIter; + Molecule::BendIterator bendIter; + Molecule::TorsionIterator torsionIter; + Molecule::InversionIterator inversionIter; + + Molecule* mol; + Atom* atom; RigidBody* rb; + Bond* bond; + Bend* bend; + Torsion* torsion; + Inversion* inversion; root_ = new TreeNode; - root_->bs.resize(nStuntDouble_); + root_->bs.resize(nObjects_); root_->bs.setAll(); // - for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) { + for (mol = info_->beginMolecule(mi); mol != NULL; + mol = info_->nextMolecule(mi)) { - std::string molName = mol->getMoleculeName(); - TreeNode* currentMolNode = createNode(root_, molName); + std::string molName = mol->getMoleculeName(); + TreeNode* molNode = createNode(root_, molName); - for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { - std::string atomName = atom->getType(); - TreeNode* currentAtomNode = createNode(currentMolNode, atomName); + for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { + std::string atomName = atom->getType(); + TreeNode* atomNode = createNode(molNode, atomName); - currentMolNode->bs.setBitOn(atom->getGlobalIndex()); - currentAtomNode->bs.setBitOn(atom->getGlobalIndex()); - } + molNode->bs.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); + atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); + } - for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { - std::string rbName = rb->getType(); - TreeNode* currentRbNode = createNode(currentMolNode, rbName); + for (rb = mol->beginRigidBody(rbIter); rb != NULL; + rb = mol->nextRigidBody(rbIter)) { + std::string rbName = rb->getType(); + TreeNode* rbNode = createNode(molNode, rbName); - currentMolNode->bs.setBitOn(rb->getGlobalIndex()); - currentRbNode->bs.setBitOn(rb->getGlobalIndex()); + molNode->bs.bitsets_[STUNTDOUBLE].setBitOn(rb->getGlobalIndex()); + rbNode->bs.bitsets_[STUNTDOUBLE].setBitOn(rb->getGlobalIndex()); - //create nodes for atoms belong to this rigidbody - for(atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) { - std::string rbAtomName = atom->getType(); - TreeNode* currentRbAtomNode = createNode(currentRbNode, rbName);; + // COMMENTED OUT because rigid bodies are IntegrableObjects + // (e.g. they are independently mobile, so selecting their + // member atoms will give some odd results if we are computing + // degrees of freedom elsewhere. - currentRbAtomNode->bs.setBitOn(atom->getGlobalIndex()); - } + // //create nodes for atoms belong to this rigidbody + // for(atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) { + // std::string rbAtomName = atom->getType(); + // TreeNode* rbAtomNode = createNode(rbNode, rbAtomName); + // rbAtomNode->bs.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); + // } + } + + for (bond = mol->beginBond(bondIter); bond != NULL; + bond = mol->nextBond(bondIter)) { + + std::string bondName = bond->getName(); + TreeNode* bondNode = createNode(molNode, bondName); + + molNode->bs.bitsets_[BOND].setBitOn(bond->getGlobalIndex()); + bondNode->bs.bitsets_[BOND].setBitOn(bond->getGlobalIndex()); + + std::vector atoms = bond->getAtoms(); + std::vector::iterator ai; + + for (ai = atoms.begin(); ai != atoms.end(); ++ai) { + std::string atomName = (*ai)->getType(); + TreeNode* atomNode = createNode(bondNode, atomName); + atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); } + } + for (bend = mol->beginBend(bendIter); bend != NULL; + bend = mol->nextBend(bendIter)) { + + std::string bendName = bend->getName(); + TreeNode* bendNode = createNode(molNode, bendName); + + molNode->bs.bitsets_[BEND].setBitOn(bend->getGlobalIndex()); + bendNode->bs.bitsets_[BEND].setBitOn(bend->getGlobalIndex()); + + std::vector atoms = bend->getAtoms(); + std::vector::iterator ai; - } + for (ai = atoms.begin(); ai != atoms.end(); ++ai) { + std::string atomName = (*ai)->getType(); + TreeNode* atomNode = createNode(bendNode, atomName); + atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); + } -} + } + for (torsion = mol->beginTorsion(torsionIter); torsion != NULL; + torsion = mol->nextTorsion(torsionIter)) { -TreeNode* NameFinder::createNode(TreeNode* parent, const std::string& name) { + std::string torsionName = torsion->getName(); + TreeNode* torsionNode = createNode(molNode, torsionName); + + molNode->bs.bitsets_[TORSION].setBitOn(torsion->getGlobalIndex()); + torsionNode->bs.bitsets_[TORSION].setBitOn(torsion->getGlobalIndex()); + + std::vector atoms = torsion->getAtoms(); + std::vector::iterator ai; + + for (ai = atoms.begin(); ai != atoms.end(); ++ai) { + std::string atomName = (*ai)->getType(); + TreeNode* atomNode = createNode(torsionNode, atomName); + atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); + } + + } + for (inversion = mol->beginInversion(inversionIter); inversion != NULL; + inversion = mol->nextInversion(inversionIter)) { + + std::string inversionName = inversion->getName(); + TreeNode* inversionNode = createNode(molNode, inversionName); + + molNode->bs.bitsets_[INVERSION].setBitOn(inversion->getGlobalIndex()); + inversionNode->bs.bitsets_[INVERSION].setBitOn(inversion->getGlobalIndex()); + std::vector atoms = inversion->getAtoms(); + std::vector::iterator ai; + + for (ai = atoms.begin(); ai != atoms.end(); ++ai) { + std::string atomName = (*ai)->getType(); + TreeNode* atomNode = createNode(inversionNode, atomName); + atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); + } + } + } + } + + TreeNode* NameFinder::createNode(TreeNode* parent, const std::string& name) { TreeNode* node; std::map::iterator foundIter; foundIter = parent->children.find(name); if ( foundIter == parent->children.end()) { - node = new TreeNode; - node->name = name; - node->bs.resize(nStuntDouble_); - parent->children.insert(std::make_pair(name, node)); + node = new TreeNode; + node->name = name; + node->bs.resize(nObjects_); + parent->children.insert(std::make_pair(name, node)); }else { - node = foundIter->second; + node = foundIter->second; } return node; -} + } -BitSet NameFinder::match(const std::string& name){ - BitSet bs(nStuntDouble_); + SelectionSet NameFinder::match(const std::string& name){ + SelectionSet bs(nObjects_); StringTokenizer tokenizer(name, "."); std::vector names; while(tokenizer.hasMoreTokens()) { - names.push_back(tokenizer.nextToken()); + names.push_back(tokenizer.nextToken()); } int size = names.size(); - switch(size) { - case 1 : - //could be molecule name, atom name and rigidbody name - matchMolecule(names[0], bs); - matchStuntDouble("*", names[0], bs); - - break; - case 2: - //could be molecule.*(include atoms and rigidbodies) or rigidbody.*(atoms belong to rigidbody) - if (!isInteger(names[1])){ - matchRigidAtoms("*", names[0], names[1], bs); - matchStuntDouble(names[0], names[1], bs); - } else { - int internalIndex = lexi_cast(names[1]); - if (internalIndex < 0) { - std::cerr << names[0] << ". " << names[1] << " is an invalid name" << std::endl; - } else { - matchInternalIndex(names[0], internalIndex, bs); - } - } - - break; - case 3: - //must be molecule.rigidbody.* - matchRigidAtoms(names[0], names[1], names[2], bs); - break; - default: - std::cerr << "invalid name: " << name << std::endl; - break; + switch(size) { + case 1 : + //could be molecule name, atom name and rigidbody name + matchMolecule(names[0], bs); + matchStuntDouble("*", names[0], bs); + matchBond("*", names[0], bs); + matchBend("*", names[0], bs); + matchTorsion("*", names[0], bs); + matchInversion("*", names[0], bs); + + break; + case 2: + //could be molecule.*(include atoms and rigidbodies) or rigidbody.*(atoms belong to rigidbody) + + if (!isInteger(names[1])){ + matchRigidAtoms("*", names[0], names[1], bs); + matchStuntDouble(names[0], names[1], bs); + } else { + int internalIndex = lexi_cast(names[1]); + if (internalIndex < 0) { + std::cerr << names[0] << ". " << names[1] << " is an invalid name" << std::endl; + } else { + matchInternalIndex(names[0], internalIndex, bs); + } + } + + break; + case 3: + //must be molecule.rigidbody.* + matchRigidAtoms(names[0], names[1], names[2], bs); + break; + default: + std::cerr << "invalid name: " << name << std::endl; + break; } return bs; -} + } -void NameFinder::matchMolecule(const std::string& molName, BitSet& bs) { + void NameFinder::matchMolecule(const std::string& molName, SelectionSet& bs) { std::vector molNodes = getMatchedChildren(root_, molName); std::vector::iterator i; for( i = molNodes.begin(); i != molNodes.end(); ++i ) { - bs |= (*i)->bs; + bs |= (*i)->bs; } -} + } -void NameFinder::matchStuntDouble(const std::string& molName, const std::string& sdName, BitSet& bs){ + void NameFinder::matchStuntDouble(const std::string& molName, const std::string& sdName, SelectionSet& bs){ std::vector molNodes = getMatchedChildren(root_, molName); std::vector::iterator i; for( i = molNodes.begin(); i != molNodes.end(); ++i ) { - std::vector sdNodes = getMatchedChildren(*i, sdName); - std::vector::iterator j; - for (j = sdNodes.begin(); j != sdNodes.end(); ++j) { - bs |= (*j)->bs; - } + std::vector sdNodes = getMatchedChildren(*i, sdName); + std::vector::iterator j; + for (j = sdNodes.begin(); j != sdNodes.end(); ++j) { + bs |= (*j)->bs; + } } -} + } -void NameFinder::matchRigidAtoms(const std::string& molName, const std::string& rbName, const std::string& rbAtomName, BitSet& bs){ + void NameFinder::matchBond(const std::string& molName, + const std::string& bondName, SelectionSet& bs){ std::vector molNodes = getMatchedChildren(root_, molName); std::vector::iterator i; for( i = molNodes.begin(); i != molNodes.end(); ++i ) { - std::vector rbNodes = getMatchedChildren(*i, rbName); - std::vector::iterator j; - for (j = rbNodes.begin(); j != rbNodes.end(); ++j) { - std::vector rbAtomNodes = getMatchedChildren(*j, rbAtomName); - std::vector::iterator k; - for(k = rbAtomNodes.begin(); k != rbAtomNodes.end(); ++k){ - bs |= (*k)->bs; - } - } + std::vector bondNodes = getMatchedChildren(*i, bondName); + std::vector::iterator j; + for (j = bondNodes.begin(); j != bondNodes.end(); ++j) { + bs |= (*j)->bs; + std::vector bondAtomNodes = getAllChildren(*j); + std::vector::iterator k; + for(k = bondAtomNodes.begin(); k != bondAtomNodes.end(); ++k){ + bs |= (*k)->bs; + } + } } + } -} + void NameFinder::matchBend(const std::string& molName, const std::string& bendName, SelectionSet& bs){ + std::vector molNodes = getMatchedChildren(root_, molName); + std::vector::iterator i; + for( i = molNodes.begin(); i != molNodes.end(); ++i ) { + std::vector bendNodes = getMatchedChildren(*i, bendName); + std::vector::iterator j; + for (j = bendNodes.begin(); j != bendNodes.end(); ++j) { + std::vector bendAtomNodes = getAllChildren(*j); + std::vector::iterator k; + for(k = bendAtomNodes.begin(); k != bendAtomNodes.end(); ++k){ + bs |= (*k)->bs; + } + } + } + } + void NameFinder::matchTorsion(const std::string& molName, const std::string& torsionName, SelectionSet& bs){ + std::vector molNodes = getMatchedChildren(root_, molName); + std::vector::iterator i; + for( i = molNodes.begin(); i != molNodes.end(); ++i ) { + std::vector torsionNodes = getMatchedChildren(*i, torsionName); + std::vector::iterator j; + for (j = torsionNodes.begin(); j != torsionNodes.end(); ++j) { + std::vector torsionAtomNodes = getAllChildren(*j); + std::vector::iterator k; + for(k = torsionAtomNodes.begin(); k != torsionAtomNodes.end(); ++k){ + bs |= (*k)->bs; + } + } + } + } + void NameFinder::matchInversion(const std::string& molName, const std::string& inversionName, SelectionSet& bs){ + std::vector molNodes = getMatchedChildren(root_, molName); + std::vector::iterator i; + for( i = molNodes.begin(); i != molNodes.end(); ++i ) { + std::vector inversionNodes = getMatchedChildren(*i, inversionName); + std::vector::iterator j; + for (j = inversionNodes.begin(); j != inversionNodes.end(); ++j) { + std::vector inversionAtomNodes = getAllChildren(*j); + std::vector::iterator k; + for(k = inversionAtomNodes.begin(); k != inversionAtomNodes.end(); ++k){ + bs |= (*k)->bs; + } + } + } + } + void NameFinder::matchRigidAtoms(const std::string& molName, const std::string& rbName, const std::string& rbAtomName, SelectionSet& bs){ + std::vector molNodes = getMatchedChildren(root_, molName); + std::vector::iterator i; + for( i = molNodes.begin(); i != molNodes.end(); ++i ) { + std::vector rbNodes = getMatchedChildren(*i, rbName); + std::vector::iterator j; + for (j = rbNodes.begin(); j != rbNodes.end(); ++j) { + std::vector rbAtomNodes = getMatchedChildren(*j, rbAtomName); + std::vector::iterator k; + for(k = rbAtomNodes.begin(); k != rbAtomNodes.end(); ++k){ + bs |= (*k)->bs; + } + } + } -std::vector NameFinder::getMatchedChildren(TreeNode* node, const std::string& name) { + } + + std::vector NameFinder::getAllChildren(TreeNode* node) { + std::vector childNodes; + std::map::iterator i; + for (i = node->children.begin(); i != node->children.end(); ++i) { + childNodes.push_back(i->second); + } + return childNodes; + } + + std::vector NameFinder::getMatchedChildren(TreeNode* node, const std::string& name) { std::vector matchedNodes; std::map::iterator i; for (i = node->children.begin(); i != node->children.end(); ++i) { - if (isMatched( i->first, name)) { - matchedNodes.push_back(i->second); - } + if (isMatched( i->first, name)) { + matchedNodes.push_back(i->second); + } } return matchedNodes; -} + } -bool NameFinder::isMatched(const std::string& str, const std::string& wildcard) { - return Wildcard::wildcardfit (wildcard.c_str(), str.c_str()); -} + bool NameFinder::isMatched(const std::string& str, const std::string& wildcard) { + return Wildcard::wildcardfit(wildcard.c_str(), str.c_str()) > 0 ? true : false; + } -void NameFinder::matchInternalIndex(const std::string& name, int internalIndex, BitSet& bs){ + void NameFinder::matchInternalIndex(const std::string& name, int internalIndex, SelectionSet& bs){ - std::map::iterator foundIter; SimInfo::MoleculeIterator mi; Molecule* mol; - for (mol = info_->beginMolecule(mi); mol != NULL; mol = info_->nextMolecule(mi)) { + for (mol = info_->beginMolecule(mi); mol != NULL; + mol = info_->nextMolecule(mi)) { - if (isMatched(mol->getMoleculeName(), name) ) { - int natoms = mol->getNAtoms(); - int nrigidbodies = mol->getNRigidBodies(); - if (internalIndex >= natoms + nrigidbodies) { - continue; - } else if (internalIndex < natoms) { - bs.setBitOn(mol->getAtomAt(internalIndex)->getGlobalIndex()); - continue; - } else if ( internalIndex < natoms + nrigidbodies) { - bs.setBitOn(mol->getRigidBodyAt(internalIndex - natoms)->getGlobalIndex()); - } - } - - } - -} - - bool NameFinder::isInteger(const std::string str) { - for(int i =0; i < str.size(); ++i){ - if (!std::isdigit(str[i])) { - return false; - } + if (isMatched(mol->getMoleculeName(), name) ) { + int natoms = mol->getNAtoms(); + int nrigidbodies = mol->getNRigidBodies(); + if (internalIndex >= natoms + nrigidbodies) { + continue; + } else if (internalIndex < natoms) { + bs.bitsets_[STUNTDOUBLE].setBitOn(mol->getAtomAt(internalIndex)->getGlobalIndex()); + continue; + } else if ( internalIndex < natoms + nrigidbodies) { + bs.bitsets_[STUNTDOUBLE].setBitOn(mol->getRigidBodyAt(internalIndex - natoms)->getGlobalIndex()); + } + } } + } + bool NameFinder::isInteger(const std::string &str) { + for(unsigned int i = 0; i < str.size(); ++i){ + if (!std::isdigit(str[i])) { + return false; + } + } return true; - } - + } }