1 |
/* |
2 |
* Copyright (c) 2005 The University of Notre Dame. All Rights Reserved. |
3 |
* |
4 |
* The University of Notre Dame grants you ("Licensee") a |
5 |
* non-exclusive, royalty free, license to use, modify and |
6 |
* redistribute this software in source and binary code form, provided |
7 |
* that the following conditions are met: |
8 |
* |
9 |
* 1. Redistributions of source code must retain the above copyright |
10 |
* notice, this list of conditions and the following disclaimer. |
11 |
* |
12 |
* 2. Redistributions in binary form must reproduce the above copyright |
13 |
* notice, this list of conditions and the following disclaimer in the |
14 |
* documentation and/or other materials provided with the |
15 |
* distribution. |
16 |
* |
17 |
* This software is provided "AS IS," without a warranty of any |
18 |
* kind. All express or implied conditions, representations and |
19 |
* warranties, including any implied warranty of merchantability, |
20 |
* fitness for a particular purpose or non-infringement, are hereby |
21 |
* excluded. The University of Notre Dame and its licensors shall not |
22 |
* be liable for any damages suffered by licensee as a result of |
23 |
* using, modifying or distributing the software or its |
24 |
* derivatives. In no event will the University of Notre Dame or its |
25 |
* licensors be liable for any lost revenue, profit or data, or for |
26 |
* direct, indirect, special, consequential, incidental or punitive |
27 |
* damages, however caused and regardless of the theory of liability, |
28 |
* arising out of the use of or inability to use software, even if the |
29 |
* University of Notre Dame has been advised of the possibility of |
30 |
* such damages. |
31 |
* |
32 |
* SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your |
33 |
* research, please cite the appropriate papers when you publish your |
34 |
* work. Good starting points are: |
35 |
* |
36 |
* [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). |
37 |
* [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). |
38 |
* [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). |
39 |
* [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). |
40 |
* [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). |
41 |
*/ |
42 |
#include "selection/NameFinder.hpp" |
43 |
#include "utils/wildcards.hpp" |
44 |
#include "utils/StringTokenizer.hpp" |
45 |
#include "primitives/Molecule.hpp" |
46 |
#include "utils/StringUtils.hpp" |
47 |
namespace OpenMD { |
48 |
|
49 |
TreeNode::~TreeNode(){ |
50 |
std::map<std::string, TreeNode*>::iterator i; |
51 |
for ( i = children.begin(); i != children.end(); ++i) { |
52 |
i->second->~TreeNode(); |
53 |
} |
54 |
children.clear(); |
55 |
} |
56 |
|
57 |
NameFinder::NameFinder(SimInfo* info) : info_(info), root_(NULL){ |
58 |
nObjects_.push_back(info_->getNGlobalAtoms()+info_->getNGlobalRigidBodies()); |
59 |
nObjects_.push_back(info_->getNGlobalBonds()); |
60 |
nObjects_.push_back(info_->getNGlobalBends()); |
61 |
nObjects_.push_back(info_->getNGlobalTorsions()); |
62 |
nObjects_.push_back(info_->getNGlobalInversions()); |
63 |
nObjects_.push_back(info_->getNGlobalMolecules()); |
64 |
|
65 |
loadNames(); |
66 |
} |
67 |
|
68 |
NameFinder::~NameFinder(){ |
69 |
delete root_; |
70 |
} |
71 |
|
72 |
void NameFinder::loadNames() { |
73 |
SimInfo::MoleculeIterator mi; |
74 |
Molecule::AtomIterator ai; |
75 |
Molecule::RigidBodyIterator rbIter; |
76 |
Molecule::BondIterator bondIter; |
77 |
Molecule::BendIterator bendIter; |
78 |
Molecule::TorsionIterator torsionIter; |
79 |
Molecule::InversionIterator inversionIter; |
80 |
|
81 |
Molecule* mol; |
82 |
Atom* atom; |
83 |
RigidBody* rb; |
84 |
Bond* bond; |
85 |
Bend* bend; |
86 |
Torsion* torsion; |
87 |
Inversion* inversion; |
88 |
|
89 |
root_ = new TreeNode; |
90 |
root_->bs.resize(nObjects_); |
91 |
root_->bs.setAll(); // |
92 |
|
93 |
for (mol = info_->beginMolecule(mi); mol != NULL; |
94 |
mol = info_->nextMolecule(mi)) { |
95 |
|
96 |
std::string molName = mol->getMoleculeName(); |
97 |
TreeNode* molNode = createNode(root_, molName); |
98 |
molNode->bs.bitsets_[MOLECULE].setBitOn(mol->getGlobalIndex()); |
99 |
|
100 |
for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { |
101 |
std::string atomName = atom->getType(); |
102 |
TreeNode* atomNode = createNode(molNode, atomName); |
103 |
|
104 |
molNode->bs.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); |
105 |
atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); |
106 |
} |
107 |
|
108 |
for (rb = mol->beginRigidBody(rbIter); rb != NULL; |
109 |
rb = mol->nextRigidBody(rbIter)) { |
110 |
std::string rbName = rb->getType(); |
111 |
TreeNode* rbNode = createNode(molNode, rbName); |
112 |
|
113 |
molNode->bs.bitsets_[STUNTDOUBLE].setBitOn(rb->getGlobalIndex()); |
114 |
rbNode->bs.bitsets_[STUNTDOUBLE].setBitOn(rb->getGlobalIndex()); |
115 |
|
116 |
// COMMENTED OUT because rigid bodies are IntegrableObjects |
117 |
// (e.g. they are independently mobile, so selecting their |
118 |
// member atoms will give some odd results if we are computing |
119 |
// degrees of freedom elsewhere. |
120 |
|
121 |
// //create nodes for atoms belong to this rigidbody |
122 |
// for(atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) { |
123 |
// std::string rbAtomName = atom->getType(); |
124 |
// TreeNode* rbAtomNode = createNode(rbNode, rbAtomName); |
125 |
|
126 |
// rbAtomNode->bs.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); |
127 |
// } |
128 |
} |
129 |
|
130 |
for (bond = mol->beginBond(bondIter); bond != NULL; |
131 |
bond = mol->nextBond(bondIter)) { |
132 |
|
133 |
std::string bondName = bond->getName(); |
134 |
TreeNode* bondNode = createNode(molNode, bondName); |
135 |
|
136 |
molNode->bs.bitsets_[BOND].setBitOn(bond->getGlobalIndex()); |
137 |
bondNode->bs.bitsets_[BOND].setBitOn(bond->getGlobalIndex()); |
138 |
|
139 |
std::vector<Atom*> atoms = bond->getAtoms(); |
140 |
std::vector<Atom*>::iterator ai; |
141 |
|
142 |
for (ai = atoms.begin(); ai != atoms.end(); ++ai) { |
143 |
std::string atomName = (*ai)->getType(); |
144 |
TreeNode* atomNode = createNode(bondNode, atomName); |
145 |
atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); |
146 |
} |
147 |
} |
148 |
for (bend = mol->beginBend(bendIter); bend != NULL; |
149 |
bend = mol->nextBend(bendIter)) { |
150 |
|
151 |
std::string bendName = bend->getName(); |
152 |
TreeNode* bendNode = createNode(molNode, bendName); |
153 |
|
154 |
molNode->bs.bitsets_[BEND].setBitOn(bend->getGlobalIndex()); |
155 |
bendNode->bs.bitsets_[BEND].setBitOn(bend->getGlobalIndex()); |
156 |
|
157 |
std::vector<Atom*> atoms = bend->getAtoms(); |
158 |
std::vector<Atom*>::iterator ai; |
159 |
|
160 |
for (ai = atoms.begin(); ai != atoms.end(); ++ai) { |
161 |
std::string atomName = (*ai)->getType(); |
162 |
TreeNode* atomNode = createNode(bendNode, atomName); |
163 |
atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); |
164 |
} |
165 |
|
166 |
} |
167 |
for (torsion = mol->beginTorsion(torsionIter); torsion != NULL; |
168 |
torsion = mol->nextTorsion(torsionIter)) { |
169 |
|
170 |
std::string torsionName = torsion->getName(); |
171 |
TreeNode* torsionNode = createNode(molNode, torsionName); |
172 |
|
173 |
molNode->bs.bitsets_[TORSION].setBitOn(torsion->getGlobalIndex()); |
174 |
torsionNode->bs.bitsets_[TORSION].setBitOn(torsion->getGlobalIndex()); |
175 |
|
176 |
std::vector<Atom*> atoms = torsion->getAtoms(); |
177 |
std::vector<Atom*>::iterator ai; |
178 |
|
179 |
for (ai = atoms.begin(); ai != atoms.end(); ++ai) { |
180 |
std::string atomName = (*ai)->getType(); |
181 |
TreeNode* atomNode = createNode(torsionNode, atomName); |
182 |
atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); |
183 |
} |
184 |
|
185 |
} |
186 |
for (inversion = mol->beginInversion(inversionIter); inversion != NULL; |
187 |
inversion = mol->nextInversion(inversionIter)) { |
188 |
|
189 |
std::string inversionName = inversion->getName(); |
190 |
TreeNode* inversionNode = createNode(molNode, inversionName); |
191 |
|
192 |
molNode->bs.bitsets_[INVERSION].setBitOn(inversion->getGlobalIndex()); |
193 |
inversionNode->bs.bitsets_[INVERSION].setBitOn(inversion->getGlobalIndex()); |
194 |
std::vector<Atom*> atoms = inversion->getAtoms(); |
195 |
std::vector<Atom*>::iterator ai; |
196 |
|
197 |
for (ai = atoms.begin(); ai != atoms.end(); ++ai) { |
198 |
std::string atomName = (*ai)->getType(); |
199 |
TreeNode* atomNode = createNode(inversionNode, atomName); |
200 |
atomNode->bs.bitsets_[STUNTDOUBLE].setBitOn((*ai)->getGlobalIndex()); |
201 |
} |
202 |
} |
203 |
} |
204 |
} |
205 |
|
206 |
TreeNode* NameFinder::createNode(TreeNode* parent, const std::string& name) { |
207 |
TreeNode* node; |
208 |
std::map<std::string, TreeNode*>::iterator foundIter; |
209 |
foundIter = parent->children.find(name); |
210 |
if ( foundIter == parent->children.end()) { |
211 |
node = new TreeNode; |
212 |
node->name = name; |
213 |
node->bs.resize(nObjects_); |
214 |
parent->children.insert(std::make_pair(name, node)); |
215 |
}else { |
216 |
node = foundIter->second; |
217 |
} |
218 |
return node; |
219 |
} |
220 |
|
221 |
SelectionSet NameFinder::match(const std::string& name){ |
222 |
SelectionSet bs(nObjects_); |
223 |
|
224 |
StringTokenizer tokenizer(name, "."); |
225 |
|
226 |
std::vector<std::string> names; |
227 |
while(tokenizer.hasMoreTokens()) { |
228 |
names.push_back(tokenizer.nextToken()); |
229 |
} |
230 |
|
231 |
int size = names.size(); |
232 |
|
233 |
switch(size) { |
234 |
case 1 : |
235 |
//could be molecule name, atom name and rigidbody name |
236 |
matchMolecule(names[0], bs); |
237 |
matchStuntDouble("*", names[0], bs); |
238 |
matchBond("*", names[0], bs); |
239 |
matchBend("*", names[0], bs); |
240 |
matchTorsion("*", names[0], bs); |
241 |
matchInversion("*", names[0], bs); |
242 |
|
243 |
break; |
244 |
case 2: |
245 |
//could be molecule.*(include atoms and rigidbodies) or rigidbody.*(atoms belong to rigidbody) |
246 |
|
247 |
if (!isInteger(names[1])){ |
248 |
matchRigidAtoms("*", names[0], names[1], bs); |
249 |
matchStuntDouble(names[0], names[1], bs); |
250 |
} else { |
251 |
int internalIndex = lexi_cast<int>(names[1]); |
252 |
if (internalIndex < 0) { |
253 |
std::cerr << names[0] << ". " << names[1] << " is an invalid name" << std::endl; |
254 |
} else { |
255 |
matchInternalIndex(names[0], internalIndex, bs); |
256 |
} |
257 |
} |
258 |
|
259 |
break; |
260 |
case 3: |
261 |
//must be molecule.rigidbody.* |
262 |
matchRigidAtoms(names[0], names[1], names[2], bs); |
263 |
break; |
264 |
default: |
265 |
std::cerr << "invalid name: " << name << std::endl; |
266 |
break; |
267 |
} |
268 |
|
269 |
return bs; |
270 |
} |
271 |
|
272 |
void NameFinder::matchMolecule(const std::string& molName, SelectionSet& bs) { |
273 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
274 |
std::vector<TreeNode*>::iterator i; |
275 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
276 |
bs |= (*i)->bs; |
277 |
} |
278 |
} |
279 |
|
280 |
void NameFinder::matchStuntDouble(const std::string& molName, const std::string& sdName, SelectionSet& bs){ |
281 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
282 |
std::vector<TreeNode*>::iterator i; |
283 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
284 |
std::vector<TreeNode*> sdNodes = getMatchedChildren(*i, sdName); |
285 |
std::vector<TreeNode*>::iterator j; |
286 |
for (j = sdNodes.begin(); j != sdNodes.end(); ++j) { |
287 |
bs |= (*j)->bs; |
288 |
} |
289 |
} |
290 |
|
291 |
} |
292 |
|
293 |
void NameFinder::matchBond(const std::string& molName, |
294 |
const std::string& bondName, SelectionSet& bs){ |
295 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
296 |
std::vector<TreeNode*>::iterator i; |
297 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
298 |
std::vector<TreeNode*> bondNodes = getMatchedChildren(*i, bondName); |
299 |
std::vector<TreeNode*>::iterator j; |
300 |
for (j = bondNodes.begin(); j != bondNodes.end(); ++j) { |
301 |
bs |= (*j)->bs; |
302 |
std::vector<TreeNode*> bondAtomNodes = getAllChildren(*j); |
303 |
std::vector<TreeNode*>::iterator k; |
304 |
for(k = bondAtomNodes.begin(); k != bondAtomNodes.end(); ++k){ |
305 |
bs |= (*k)->bs; |
306 |
} |
307 |
} |
308 |
} |
309 |
} |
310 |
|
311 |
void NameFinder::matchBend(const std::string& molName, const std::string& bendName, SelectionSet& bs){ |
312 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
313 |
std::vector<TreeNode*>::iterator i; |
314 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
315 |
std::vector<TreeNode*> bendNodes = getMatchedChildren(*i, bendName); |
316 |
std::vector<TreeNode*>::iterator j; |
317 |
for (j = bendNodes.begin(); j != bendNodes.end(); ++j) { |
318 |
std::vector<TreeNode*> bendAtomNodes = getAllChildren(*j); |
319 |
std::vector<TreeNode*>::iterator k; |
320 |
for(k = bendAtomNodes.begin(); k != bendAtomNodes.end(); ++k){ |
321 |
bs |= (*k)->bs; |
322 |
} |
323 |
} |
324 |
} |
325 |
} |
326 |
void NameFinder::matchTorsion(const std::string& molName, const std::string& torsionName, SelectionSet& bs){ |
327 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
328 |
std::vector<TreeNode*>::iterator i; |
329 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
330 |
std::vector<TreeNode*> torsionNodes = getMatchedChildren(*i, torsionName); |
331 |
std::vector<TreeNode*>::iterator j; |
332 |
for (j = torsionNodes.begin(); j != torsionNodes.end(); ++j) { |
333 |
std::vector<TreeNode*> torsionAtomNodes = getAllChildren(*j); |
334 |
std::vector<TreeNode*>::iterator k; |
335 |
for(k = torsionAtomNodes.begin(); k != torsionAtomNodes.end(); ++k){ |
336 |
bs |= (*k)->bs; |
337 |
} |
338 |
} |
339 |
} |
340 |
} |
341 |
void NameFinder::matchInversion(const std::string& molName, const std::string& inversionName, SelectionSet& bs){ |
342 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
343 |
std::vector<TreeNode*>::iterator i; |
344 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
345 |
std::vector<TreeNode*> inversionNodes = getMatchedChildren(*i, inversionName); |
346 |
std::vector<TreeNode*>::iterator j; |
347 |
for (j = inversionNodes.begin(); j != inversionNodes.end(); ++j) { |
348 |
std::vector<TreeNode*> inversionAtomNodes = getAllChildren(*j); |
349 |
std::vector<TreeNode*>::iterator k; |
350 |
for(k = inversionAtomNodes.begin(); k != inversionAtomNodes.end(); ++k){ |
351 |
bs |= (*k)->bs; |
352 |
} |
353 |
} |
354 |
} |
355 |
} |
356 |
|
357 |
void NameFinder::matchRigidAtoms(const std::string& molName, const std::string& rbName, const std::string& rbAtomName, SelectionSet& bs){ |
358 |
std::vector<TreeNode*> molNodes = getMatchedChildren(root_, molName); |
359 |
std::vector<TreeNode*>::iterator i; |
360 |
for( i = molNodes.begin(); i != molNodes.end(); ++i ) { |
361 |
std::vector<TreeNode*> rbNodes = getMatchedChildren(*i, rbName); |
362 |
std::vector<TreeNode*>::iterator j; |
363 |
for (j = rbNodes.begin(); j != rbNodes.end(); ++j) { |
364 |
std::vector<TreeNode*> rbAtomNodes = getMatchedChildren(*j, rbAtomName); |
365 |
std::vector<TreeNode*>::iterator k; |
366 |
for(k = rbAtomNodes.begin(); k != rbAtomNodes.end(); ++k){ |
367 |
bs |= (*k)->bs; |
368 |
} |
369 |
} |
370 |
} |
371 |
|
372 |
} |
373 |
|
374 |
std::vector<TreeNode*> NameFinder::getAllChildren(TreeNode* node) { |
375 |
std::vector<TreeNode*> childNodes; |
376 |
std::map<std::string, TreeNode*>::iterator i; |
377 |
for (i = node->children.begin(); i != node->children.end(); ++i) { |
378 |
childNodes.push_back(i->second); |
379 |
} |
380 |
return childNodes; |
381 |
} |
382 |
|
383 |
std::vector<TreeNode*> NameFinder::getMatchedChildren(TreeNode* node, const std::string& name) { |
384 |
std::vector<TreeNode*> matchedNodes; |
385 |
std::map<std::string, TreeNode*>::iterator i; |
386 |
for (i = node->children.begin(); i != node->children.end(); ++i) { |
387 |
if (isMatched( i->first, name)) { |
388 |
matchedNodes.push_back(i->second); |
389 |
} |
390 |
} |
391 |
|
392 |
return matchedNodes; |
393 |
} |
394 |
|
395 |
bool NameFinder::isMatched(const std::string& str, const std::string& wildcard) { |
396 |
return Wildcard::wildcardfit(wildcard.c_str(), str.c_str()) > 0 ? true : false; |
397 |
} |
398 |
|
399 |
|
400 |
void NameFinder::matchInternalIndex(const std::string& name, int internalIndex, SelectionSet& bs){ |
401 |
|
402 |
SimInfo::MoleculeIterator mi; |
403 |
Molecule* mol; |
404 |
|
405 |
for (mol = info_->beginMolecule(mi); mol != NULL; |
406 |
mol = info_->nextMolecule(mi)) { |
407 |
|
408 |
if (isMatched(mol->getMoleculeName(), name) ) { |
409 |
int natoms = mol->getNAtoms(); |
410 |
int nrigidbodies = mol->getNRigidBodies(); |
411 |
if (internalIndex >= natoms + nrigidbodies) { |
412 |
continue; |
413 |
} else if (internalIndex < natoms) { |
414 |
bs.bitsets_[STUNTDOUBLE].setBitOn(mol->getAtomAt(internalIndex)->getGlobalIndex()); |
415 |
continue; |
416 |
} else if ( internalIndex < natoms + nrigidbodies) { |
417 |
bs.bitsets_[STUNTDOUBLE].setBitOn(mol->getRigidBodyAt(internalIndex - natoms)->getGlobalIndex()); |
418 |
} |
419 |
} |
420 |
} |
421 |
} |
422 |
|
423 |
bool NameFinder::isInteger(const std::string &str) { |
424 |
for(unsigned int i = 0; i < str.size(); ++i){ |
425 |
if (!std::isdigit(str[i])) { |
426 |
return false; |
427 |
} |
428 |
} |
429 |
return true; |
430 |
} |
431 |
} |