6 |
|
* redistribute this software in source and binary code form, provided |
7 |
|
* that the following conditions are met: |
8 |
|
* |
9 |
< |
* 1. Acknowledgement of the program authors must be made in any |
10 |
< |
* publication of scientific results based in part on use of the |
11 |
< |
* program. An acceptable form of acknowledgement is citation of |
12 |
< |
* the article in which the program was described (Matthew |
13 |
< |
* A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher |
14 |
< |
* J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented |
15 |
< |
* Parallel Simulation Engine for Molecular Dynamics," |
16 |
< |
* J. Comput. Chem. 26, pp. 252-271 (2005)) |
17 |
< |
* |
18 |
< |
* 2. Redistributions of source code must retain the above copyright |
9 |
> |
* 1. Redistributions of source code must retain the above copyright |
10 |
|
* notice, this list of conditions and the following disclaimer. |
11 |
|
* |
12 |
< |
* 3. Redistributions in binary form must reproduce the above copyright |
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. |
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 |
|
|
43 |
|
#include <stack> |
46 |
|
#include "primitives/DirectionalAtom.hpp" |
47 |
|
#include "primitives/RigidBody.hpp" |
48 |
|
#include "primitives/Molecule.hpp" |
49 |
+ |
#include "io/ifstrstream.hpp" |
50 |
+ |
#include "types/FixedChargeAdapter.hpp" |
51 |
+ |
#include "types/FluctuatingChargeAdapter.hpp" |
52 |
|
|
49 |
– |
namespace oopse { |
53 |
|
|
54 |
+ |
namespace OpenMD { |
55 |
|
|
52 |
– |
SelectionEvaluator::SelectionEvaluator(SimInfo* si) |
53 |
– |
: info(si), nameFinder(info), distanceFinder(info), indexFinder(info), isLoaded_(false){ |
54 |
– |
|
55 |
– |
nStuntDouble = info->getNGlobalAtoms() + info->getNRigidBodies(); |
56 |
– |
} |
56 |
|
|
57 |
< |
bool SelectionEvaluator::loadScript(const std::string& filename, const std::string& script) { |
57 |
> |
SelectionEvaluator::SelectionEvaluator(SimInfo* si) |
58 |
> |
: info(si), nameFinder(info), distanceFinder(info), hullFinder(info), |
59 |
> |
indexFinder(info), hasSurfaceArea_(false), |
60 |
> |
isLoaded_(false){ |
61 |
> |
nObjects.push_back(info->getNGlobalAtoms() + info->getNGlobalRigidBodies()); |
62 |
> |
nObjects.push_back(info->getNGlobalBonds()); |
63 |
> |
nObjects.push_back(info->getNGlobalBends()); |
64 |
> |
nObjects.push_back(info->getNGlobalTorsions()); |
65 |
> |
nObjects.push_back(info->getNGlobalInversions()); |
66 |
> |
nObjects.push_back(info->getNGlobalMolecules()); |
67 |
> |
} |
68 |
> |
|
69 |
> |
bool SelectionEvaluator::loadScript(const std::string& filename, |
70 |
> |
const std::string& script) { |
71 |
|
clearDefinitionsAndLoadPredefined(); |
72 |
|
this->filename = filename; |
73 |
|
this->script = script; |
74 |
|
if (! compiler.compile(filename, script)) { |
75 |
< |
error = true; |
76 |
< |
errorMessage = compiler.getErrorMessage(); |
77 |
< |
std::cerr << "SelectionCompiler Error: " << errorMessage << std::endl; |
78 |
< |
return false; |
75 |
> |
error = true; |
76 |
> |
errorMessage = compiler.getErrorMessage(); |
77 |
> |
|
78 |
> |
sprintf( painCave.errMsg, |
79 |
> |
"SelectionCompiler Error: %s\n", errorMessage.c_str()); |
80 |
> |
painCave.severity = OPENMD_ERROR; |
81 |
> |
painCave.isFatal = 1; |
82 |
> |
simError(); |
83 |
> |
return false; |
84 |
|
} |
85 |
|
|
86 |
|
pc = 0; |
92 |
|
|
93 |
|
isDynamic_ = false; |
94 |
|
for (i = aatoken.begin(); i != aatoken.end(); ++i) { |
95 |
< |
if (containDynamicToken(*i)) { |
96 |
< |
isDynamic_ = true; |
97 |
< |
break; |
98 |
< |
} |
95 |
> |
if (containDynamicToken(*i)) { |
96 |
> |
isDynamic_ = true; |
97 |
> |
break; |
98 |
> |
} |
99 |
|
} |
100 |
|
|
101 |
|
isLoaded_ = true; |
102 |
|
return true; |
103 |
< |
} |
103 |
> |
} |
104 |
|
|
105 |
< |
void SelectionEvaluator::clearState() { |
89 |
< |
//for (int i = scriptLevelMax; --i >= 0; ) |
90 |
< |
// stack[i].clear(); |
91 |
< |
//scriptLevel = 0; |
105 |
> |
void SelectionEvaluator::clearState() { |
106 |
|
error = false; |
107 |
|
errorMessage = ""; |
108 |
< |
} |
108 |
> |
} |
109 |
|
|
110 |
< |
bool SelectionEvaluator::loadScriptString(const std::string& script) { |
110 |
> |
bool SelectionEvaluator::loadScriptString(const std::string& script) { |
111 |
|
clearState(); |
112 |
|
return loadScript("", script); |
113 |
< |
} |
113 |
> |
} |
114 |
|
|
115 |
< |
bool SelectionEvaluator::loadScriptFile(const std::string& filename) { |
115 |
> |
bool SelectionEvaluator::loadScriptFile(const std::string& filename) { |
116 |
|
clearState(); |
117 |
|
return loadScriptFileInternal(filename); |
118 |
< |
} |
118 |
> |
} |
119 |
|
|
120 |
< |
bool SelectionEvaluator::loadScriptFileInternal(const std::string & filename) { |
121 |
< |
std::ifstream ifs(filename.c_str()); |
120 |
> |
bool SelectionEvaluator::loadScriptFileInternal(const std::string & filename) { |
121 |
> |
ifstrstream ifs(filename.c_str()); |
122 |
|
if (!ifs.is_open()) { |
123 |
< |
return false; |
123 |
> |
return false; |
124 |
|
} |
125 |
< |
|
125 |
> |
|
126 |
|
const int bufferSize = 65535; |
127 |
|
char buffer[bufferSize]; |
128 |
|
std::string script; |
129 |
|
while(ifs.getline(buffer, bufferSize)) { |
130 |
< |
script += buffer; |
130 |
> |
script += buffer; |
131 |
|
} |
132 |
|
return loadScript(filename, script); |
133 |
< |
} |
134 |
< |
|
135 |
< |
void SelectionEvaluator::instructionDispatchLoop(BitSet& bs){ |
133 |
> |
} |
134 |
> |
|
135 |
> |
void SelectionEvaluator::instructionDispatchLoop(SelectionSet& bs){ |
136 |
|
|
137 |
|
while ( pc < aatoken.size()) { |
138 |
< |
statement = aatoken[pc++]; |
139 |
< |
statementLength = statement.size(); |
140 |
< |
Token token = statement[0]; |
141 |
< |
switch (token.tok) { |
142 |
< |
case Token::define: |
143 |
< |
define(); |
144 |
< |
break; |
145 |
< |
case Token::select: |
146 |
< |
select(bs); |
147 |
< |
break; |
148 |
< |
default: |
149 |
< |
unrecognizedCommand(token); |
150 |
< |
return; |
151 |
< |
} |
138 |
> |
statement = aatoken[pc++]; |
139 |
> |
statementLength = statement.size(); |
140 |
> |
Token token = statement[0]; |
141 |
> |
switch (token.tok) { |
142 |
> |
case Token::define: |
143 |
> |
define(); |
144 |
> |
break; |
145 |
> |
case Token::select: |
146 |
> |
select(bs); |
147 |
> |
break; |
148 |
> |
default: |
149 |
> |
unrecognizedCommand(token); |
150 |
> |
return; |
151 |
> |
} |
152 |
|
} |
153 |
|
|
154 |
< |
} |
154 |
> |
} |
155 |
|
|
156 |
< |
BitSet SelectionEvaluator::expression(const std::vector<Token>& code, int pcStart) { |
143 |
< |
BitSet bs; |
144 |
< |
std::stack<BitSet> stack; |
156 |
> |
void SelectionEvaluator::instructionDispatchLoop(SelectionSet& bs, int frame){ |
157 |
|
|
158 |
< |
for (int pc = pcStart; pc < code.size(); ++pc) { |
158 |
> |
while ( pc < aatoken.size()) { |
159 |
> |
statement = aatoken[pc++]; |
160 |
> |
statementLength = statement.size(); |
161 |
> |
Token token = statement[0]; |
162 |
> |
switch (token.tok) { |
163 |
> |
case Token::define: |
164 |
> |
define(); |
165 |
> |
break; |
166 |
> |
case Token::select: |
167 |
> |
select(bs, frame); |
168 |
> |
break; |
169 |
> |
default: |
170 |
> |
unrecognizedCommand(token); |
171 |
> |
return; |
172 |
> |
} |
173 |
> |
} |
174 |
> |
|
175 |
> |
} |
176 |
> |
|
177 |
> |
SelectionSet SelectionEvaluator::expression(const std::vector<Token>& code, |
178 |
> |
int pcStart) { |
179 |
> |
SelectionSet bs = createSelectionSets(); |
180 |
> |
std::stack<SelectionSet> stack; |
181 |
> |
vector<int> bsSize = bs.size(); |
182 |
> |
|
183 |
> |
for (unsigned int pc = pcStart; pc < code.size(); ++pc) { |
184 |
|
Token instruction = code[pc]; |
185 |
|
|
186 |
|
switch (instruction.tok) { |
189 |
|
case Token::expressionEnd: |
190 |
|
break; |
191 |
|
case Token::all: |
192 |
< |
bs = BitSet(nStuntDouble); |
193 |
< |
bs.setAll(); |
192 |
> |
bs = allInstruction(); |
193 |
> |
stack.push(bs); |
194 |
> |
break; |
195 |
> |
case Token::none: |
196 |
> |
bs = createSelectionSets(); |
197 |
> |
stack.push(bs); |
198 |
> |
break; |
199 |
> |
case Token::opOr: |
200 |
> |
bs= stack.top(); |
201 |
> |
stack.pop(); |
202 |
> |
stack.top() |= bs; |
203 |
> |
break; |
204 |
> |
case Token::opAnd: |
205 |
> |
bs = stack.top(); |
206 |
> |
stack.pop(); |
207 |
> |
stack.top() &= bs; |
208 |
> |
break; |
209 |
> |
case Token::opNot: |
210 |
> |
stack.top().flip(); |
211 |
> |
break; |
212 |
> |
case Token::within: |
213 |
> |
withinInstruction(instruction, stack.top()); |
214 |
> |
break; |
215 |
> |
case Token::hull: |
216 |
> |
stack.push(hull()); |
217 |
> |
break; |
218 |
> |
//case Token::selected: |
219 |
> |
// stack.push(getSelectionSet()); |
220 |
> |
// break; |
221 |
> |
case Token::name: |
222 |
> |
stack.push(nameInstruction(boost::any_cast<std::string>(instruction.value))); |
223 |
> |
break; |
224 |
> |
case Token::index: |
225 |
> |
stack.push(indexInstruction(instruction.value)); |
226 |
> |
break; |
227 |
> |
case Token::identifier: |
228 |
> |
stack.push(lookupValue(boost::any_cast<std::string>(instruction.value))); |
229 |
> |
break; |
230 |
> |
case Token::opLT: |
231 |
> |
case Token::opLE: |
232 |
> |
case Token::opGE: |
233 |
> |
case Token::opGT: |
234 |
> |
case Token::opEQ: |
235 |
> |
case Token::opNE: |
236 |
> |
stack.push(comparatorInstruction(instruction)); |
237 |
> |
break; |
238 |
> |
default: |
239 |
> |
unrecognizedExpression(); |
240 |
> |
} |
241 |
> |
} |
242 |
> |
if (stack.size() != 1) |
243 |
> |
evalError("atom expression compiler error - stack over/underflow"); |
244 |
> |
|
245 |
> |
return stack.top(); |
246 |
> |
} |
247 |
> |
|
248 |
> |
|
249 |
> |
SelectionSet SelectionEvaluator::expression(const std::vector<Token>& code, |
250 |
> |
int pcStart, int frame) { |
251 |
> |
SelectionSet bs = createSelectionSets(); |
252 |
> |
std::stack<SelectionSet> stack; |
253 |
> |
|
254 |
> |
for (unsigned int pc = pcStart; pc < code.size(); ++pc) { |
255 |
> |
Token instruction = code[pc]; |
256 |
> |
|
257 |
> |
switch (instruction.tok) { |
258 |
> |
case Token::expressionBegin: |
259 |
> |
break; |
260 |
> |
case Token::expressionEnd: |
261 |
> |
break; |
262 |
> |
case Token::all: |
263 |
> |
bs = allInstruction(); |
264 |
|
stack.push(bs); |
265 |
|
break; |
266 |
|
case Token::none: |
267 |
< |
bs = BitSet(nStuntDouble); |
267 |
> |
bs = SelectionSet(nObjects); |
268 |
|
stack.push(bs); |
269 |
|
break; |
270 |
|
case Token::opOr: |
281 |
|
stack.top().flip(); |
282 |
|
break; |
283 |
|
case Token::within: |
284 |
< |
|
178 |
< |
withinInstruction(instruction, stack.top()); |
284 |
> |
withinInstruction(instruction, stack.top(), frame); |
285 |
|
break; |
286 |
< |
//case Token::selected: |
287 |
< |
// stack.push(getSelectionSet()); |
288 |
< |
// break; |
286 |
> |
case Token::hull: |
287 |
> |
stack.push(hull(frame)); |
288 |
> |
break; |
289 |
> |
//case Token::selected: |
290 |
> |
// stack.push(getSelectionSet()); |
291 |
> |
// break; |
292 |
|
case Token::name: |
293 |
|
stack.push(nameInstruction(boost::any_cast<std::string>(instruction.value))); |
294 |
|
break; |
304 |
|
case Token::opGT: |
305 |
|
case Token::opEQ: |
306 |
|
case Token::opNE: |
307 |
< |
stack.push(comparatorInstruction(instruction)); |
307 |
> |
stack.push(comparatorInstruction(instruction, frame)); |
308 |
|
break; |
309 |
|
default: |
310 |
|
unrecognizedExpression(); |
318 |
|
|
319 |
|
|
320 |
|
|
321 |
< |
BitSet SelectionEvaluator::comparatorInstruction(const Token& instruction) { |
321 |
> |
SelectionSet SelectionEvaluator::comparatorInstruction(const Token& instruction) { |
322 |
|
int comparator = instruction.tok; |
323 |
|
int property = instruction.intValue; |
324 |
|
float comparisonValue = boost::any_cast<float>(instruction.value); |
325 |
< |
float propertyValue; |
217 |
< |
BitSet bs(nStuntDouble); |
325 |
> |
SelectionSet bs = createSelectionSets(); |
326 |
|
bs.clearAll(); |
327 |
|
|
328 |
|
SimInfo::MoleculeIterator mi; |
331 |
|
Atom* atom; |
332 |
|
Molecule::RigidBodyIterator rbIter; |
333 |
|
RigidBody* rb; |
334 |
+ |
|
335 |
+ |
for (mol = info->beginMolecule(mi); mol != NULL; |
336 |
+ |
mol = info->nextMolecule(mi)) { |
337 |
+ |
|
338 |
+ |
for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { |
339 |
+ |
compareProperty(atom, bs, property, comparator, comparisonValue); |
340 |
+ |
} |
341 |
+ |
|
342 |
+ |
for (rb = mol->beginRigidBody(rbIter); rb != NULL; |
343 |
+ |
rb = mol->nextRigidBody(rbIter)) { |
344 |
+ |
compareProperty(rb, bs, property, comparator, comparisonValue); |
345 |
+ |
} |
346 |
+ |
} |
347 |
+ |
|
348 |
+ |
return bs; |
349 |
+ |
} |
350 |
+ |
|
351 |
+ |
SelectionSet SelectionEvaluator::comparatorInstruction(const Token& instruction, int frame) { |
352 |
+ |
int comparator = instruction.tok; |
353 |
+ |
int property = instruction.intValue; |
354 |
+ |
float comparisonValue = boost::any_cast<float>(instruction.value); |
355 |
+ |
SelectionSet bs = createSelectionSets(); |
356 |
+ |
bs.clearAll(); |
357 |
|
|
358 |
< |
for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) { |
358 |
> |
SimInfo::MoleculeIterator mi; |
359 |
> |
Molecule* mol; |
360 |
> |
Molecule::AtomIterator ai; |
361 |
> |
Atom* atom; |
362 |
> |
Molecule::RigidBodyIterator rbIter; |
363 |
> |
RigidBody* rb; |
364 |
|
|
365 |
< |
for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { |
366 |
< |
compareProperty(atom, bs, property, comparator, comparisonValue); |
367 |
< |
} |
368 |
< |
|
369 |
< |
//change the positions of atoms which belong to the rigidbodies |
370 |
< |
for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) { |
371 |
< |
compareProperty(rb, bs, property, comparator, comparisonValue); |
372 |
< |
} |
365 |
> |
for (mol = info->beginMolecule(mi); mol != NULL; |
366 |
> |
mol = info->nextMolecule(mi)) { |
367 |
> |
|
368 |
> |
for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { |
369 |
> |
compareProperty(atom, bs, property, comparator, comparisonValue, frame); |
370 |
> |
} |
371 |
> |
|
372 |
> |
for (rb = mol->beginRigidBody(rbIter); rb != NULL; |
373 |
> |
rb = mol->nextRigidBody(rbIter)) { |
374 |
> |
compareProperty(rb, bs, property, comparator, comparisonValue, frame); |
375 |
> |
} |
376 |
|
} |
377 |
|
|
378 |
|
return bs; |
379 |
< |
} |
379 |
> |
} |
380 |
|
|
381 |
< |
void SelectionEvaluator::compareProperty(StuntDouble* sd, BitSet& bs, int property, int comparator, float comparisonValue) { |
382 |
< |
double propertyValue; |
383 |
< |
switch (property) { |
384 |
< |
case Token::mass: |
385 |
< |
propertyValue = sd->getMass(); |
386 |
< |
break; |
387 |
< |
case Token::charge: |
388 |
< |
return; |
389 |
< |
//break; |
390 |
< |
case Token::dipole: |
391 |
< |
return; |
392 |
< |
//break; |
393 |
< |
default: |
394 |
< |
unrecognizedAtomProperty(property); |
395 |
< |
} |
381 |
> |
void SelectionEvaluator::compareProperty(StuntDouble* sd, SelectionSet& bs, |
382 |
> |
int property, int comparator, |
383 |
> |
float comparisonValue) { |
384 |
> |
RealType propertyValue = 0.0; |
385 |
> |
Vector3d pos; |
386 |
> |
|
387 |
> |
switch (property) { |
388 |
> |
case Token::mass: |
389 |
> |
propertyValue = sd->getMass(); |
390 |
> |
break; |
391 |
> |
case Token::charge: |
392 |
> |
if (sd->isAtom()){ |
393 |
> |
Atom* atom = static_cast<Atom*>(sd); |
394 |
> |
propertyValue = getCharge(atom); |
395 |
> |
} else if (sd->isRigidBody()) { |
396 |
> |
RigidBody* rb = static_cast<RigidBody*>(sd); |
397 |
> |
RigidBody::AtomIterator ai; |
398 |
> |
Atom* atom; |
399 |
> |
for (atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) { |
400 |
> |
propertyValue+= getCharge(atom); |
401 |
> |
} |
402 |
> |
} |
403 |
> |
break; |
404 |
> |
case Token::x: |
405 |
> |
propertyValue = sd->getPos().x(); |
406 |
> |
break; |
407 |
> |
case Token::y: |
408 |
> |
propertyValue = sd->getPos().y(); |
409 |
> |
break; |
410 |
> |
case Token::z: |
411 |
> |
propertyValue = sd->getPos().z(); |
412 |
> |
break; |
413 |
> |
case Token::wrappedX: |
414 |
> |
pos = sd->getPos(); |
415 |
> |
info->getSnapshotManager()->getCurrentSnapshot()->wrapVector(pos); |
416 |
> |
propertyValue = pos.x(); |
417 |
> |
break; |
418 |
> |
case Token::wrappedY: |
419 |
> |
pos = sd->getPos(); |
420 |
> |
info->getSnapshotManager()->getCurrentSnapshot()->wrapVector(pos); |
421 |
> |
propertyValue = pos.y(); |
422 |
> |
break; |
423 |
> |
case Token::wrappedZ: |
424 |
> |
pos = sd->getPos(); |
425 |
> |
info->getSnapshotManager()->getCurrentSnapshot()->wrapVector(pos); |
426 |
> |
propertyValue = pos.z(); |
427 |
> |
break; |
428 |
> |
case Token::r: |
429 |
> |
propertyValue = sd->getPos().length(); |
430 |
> |
break; |
431 |
> |
default: |
432 |
> |
unrecognizedAtomProperty(property); |
433 |
> |
} |
434 |
|
|
435 |
< |
bool match = false; |
436 |
< |
switch (comparator) { |
437 |
< |
case Token::opLT: |
438 |
< |
match = propertyValue < comparisonValue; |
439 |
< |
break; |
440 |
< |
case Token::opLE: |
441 |
< |
match = propertyValue <= comparisonValue; |
442 |
< |
break; |
443 |
< |
case Token::opGE: |
444 |
< |
match = propertyValue >= comparisonValue; |
445 |
< |
break; |
446 |
< |
case Token::opGT: |
447 |
< |
match = propertyValue > comparisonValue; |
448 |
< |
break; |
449 |
< |
case Token::opEQ: |
450 |
< |
match = propertyValue == comparisonValue; |
451 |
< |
break; |
452 |
< |
case Token::opNE: |
453 |
< |
match = propertyValue != comparisonValue; |
454 |
< |
break; |
455 |
< |
} |
279 |
< |
if (match) |
280 |
< |
bs.setBitOn(sd->getGlobalIndex()); |
435 |
> |
bool match = false; |
436 |
> |
switch (comparator) { |
437 |
> |
case Token::opLT: |
438 |
> |
match = propertyValue < comparisonValue; |
439 |
> |
break; |
440 |
> |
case Token::opLE: |
441 |
> |
match = propertyValue <= comparisonValue; |
442 |
> |
break; |
443 |
> |
case Token::opGE: |
444 |
> |
match = propertyValue >= comparisonValue; |
445 |
> |
break; |
446 |
> |
case Token::opGT: |
447 |
> |
match = propertyValue > comparisonValue; |
448 |
> |
break; |
449 |
> |
case Token::opEQ: |
450 |
> |
match = propertyValue == comparisonValue; |
451 |
> |
break; |
452 |
> |
case Token::opNE: |
453 |
> |
match = propertyValue != comparisonValue; |
454 |
> |
break; |
455 |
> |
} |
456 |
|
|
457 |
< |
} |
457 |
> |
if (match) |
458 |
> |
bs.bitsets_[STUNTDOUBLE].setBitOn(sd->getGlobalIndex()); |
459 |
|
|
460 |
< |
void SelectionEvaluator::withinInstruction(const Token& instruction, BitSet& bs){ |
460 |
> |
|
461 |
> |
} |
462 |
> |
|
463 |
> |
void SelectionEvaluator::compareProperty(StuntDouble* sd, SelectionSet& bs, |
464 |
> |
int property, int comparator, |
465 |
> |
float comparisonValue, int frame) { |
466 |
> |
RealType propertyValue = 0.0; |
467 |
> |
Vector3d pos; |
468 |
> |
switch (property) { |
469 |
> |
case Token::mass: |
470 |
> |
propertyValue = sd->getMass(); |
471 |
> |
break; |
472 |
> |
case Token::charge: |
473 |
> |
if (sd->isAtom()){ |
474 |
> |
Atom* atom = static_cast<Atom*>(sd); |
475 |
> |
propertyValue = getCharge(atom,frame); |
476 |
> |
} else if (sd->isRigidBody()) { |
477 |
> |
RigidBody* rb = static_cast<RigidBody*>(sd); |
478 |
> |
RigidBody::AtomIterator ai; |
479 |
> |
Atom* atom; |
480 |
> |
for (atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) { |
481 |
> |
propertyValue+= getCharge(atom,frame); |
482 |
> |
} |
483 |
> |
} |
484 |
> |
break; |
485 |
> |
case Token::x: |
486 |
> |
propertyValue = sd->getPos(frame).x(); |
487 |
> |
break; |
488 |
> |
case Token::y: |
489 |
> |
propertyValue = sd->getPos(frame).y(); |
490 |
> |
break; |
491 |
> |
case Token::z: |
492 |
> |
propertyValue = sd->getPos(frame).z(); |
493 |
> |
break; |
494 |
> |
case Token::wrappedX: |
495 |
> |
pos = sd->getPos(frame); |
496 |
> |
info->getSnapshotManager()->getSnapshot(frame)->wrapVector(pos); |
497 |
> |
propertyValue = pos.x(); |
498 |
> |
break; |
499 |
> |
case Token::wrappedY: |
500 |
> |
pos = sd->getPos(frame); |
501 |
> |
info->getSnapshotManager()->getSnapshot(frame)->wrapVector(pos); |
502 |
> |
propertyValue = pos.y(); |
503 |
> |
break; |
504 |
> |
case Token::wrappedZ: |
505 |
> |
pos = sd->getPos(frame); |
506 |
> |
info->getSnapshotManager()->getSnapshot(frame)->wrapVector(pos); |
507 |
> |
propertyValue = pos.z(); |
508 |
> |
break; |
509 |
> |
|
510 |
> |
case Token::r: |
511 |
> |
propertyValue = sd->getPos(frame).length(); |
512 |
> |
break; |
513 |
> |
default: |
514 |
> |
unrecognizedAtomProperty(property); |
515 |
> |
} |
516 |
> |
|
517 |
> |
bool match = false; |
518 |
> |
switch (comparator) { |
519 |
> |
case Token::opLT: |
520 |
> |
match = propertyValue < comparisonValue; |
521 |
> |
break; |
522 |
> |
case Token::opLE: |
523 |
> |
match = propertyValue <= comparisonValue; |
524 |
> |
break; |
525 |
> |
case Token::opGE: |
526 |
> |
match = propertyValue >= comparisonValue; |
527 |
> |
break; |
528 |
> |
case Token::opGT: |
529 |
> |
match = propertyValue > comparisonValue; |
530 |
> |
break; |
531 |
> |
case Token::opEQ: |
532 |
> |
match = propertyValue == comparisonValue; |
533 |
> |
break; |
534 |
> |
case Token::opNE: |
535 |
> |
match = propertyValue != comparisonValue; |
536 |
> |
break; |
537 |
> |
} |
538 |
> |
if (match) |
539 |
> |
bs.bitsets_[STUNTDOUBLE].setBitOn(sd->getGlobalIndex()); |
540 |
> |
|
541 |
|
|
542 |
+ |
} |
543 |
+ |
|
544 |
+ |
void SelectionEvaluator::withinInstruction(const Token& instruction, |
545 |
+ |
SelectionSet& bs){ |
546 |
+ |
|
547 |
|
boost::any withinSpec = instruction.value; |
548 |
|
float distance; |
549 |
|
if (withinSpec.type() == typeid(float)){ |
550 |
< |
distance = boost::any_cast<float>(withinSpec); |
550 |
> |
distance = boost::any_cast<float>(withinSpec); |
551 |
|
} else if (withinSpec.type() == typeid(int)) { |
552 |
< |
distance = boost::any_cast<int>(withinSpec); |
552 |
> |
distance = boost::any_cast<int>(withinSpec); |
553 |
|
} else { |
554 |
< |
evalError("casting error in withinInstruction"); |
555 |
< |
bs.clearAll(); |
554 |
> |
evalError("casting error in withinInstruction"); |
555 |
> |
bs.clearAll(); |
556 |
|
} |
557 |
|
|
558 |
|
bs = distanceFinder.find(bs, distance); |
559 |
< |
} |
559 |
> |
} |
560 |
|
|
561 |
< |
void SelectionEvaluator::define() { |
561 |
> |
void SelectionEvaluator::withinInstruction(const Token& instruction, |
562 |
> |
SelectionSet& bs, int frame){ |
563 |
> |
|
564 |
> |
boost::any withinSpec = instruction.value; |
565 |
> |
float distance; |
566 |
> |
if (withinSpec.type() == typeid(float)){ |
567 |
> |
distance = boost::any_cast<float>(withinSpec); |
568 |
> |
} else if (withinSpec.type() == typeid(int)) { |
569 |
> |
distance = boost::any_cast<int>(withinSpec); |
570 |
> |
} else { |
571 |
> |
evalError("casting error in withinInstruction"); |
572 |
> |
bs.clearAll(); |
573 |
> |
} |
574 |
> |
|
575 |
> |
bs = distanceFinder.find(bs, distance, frame); |
576 |
> |
} |
577 |
> |
|
578 |
> |
void SelectionEvaluator::define() { |
579 |
|
assert(statement.size() >= 3); |
580 |
< |
|
580 |
> |
|
581 |
|
std::string variable = boost::any_cast<std::string>(statement[1].value); |
582 |
+ |
|
583 |
+ |
variables.insert(VariablesType::value_type(variable, |
584 |
+ |
expression(statement, 2))); |
585 |
+ |
} |
586 |
+ |
|
587 |
|
|
588 |
< |
variables.insert(VariablesType::value_type(variable, expression(statement, 2))); |
589 |
< |
} |
590 |
< |
|
308 |
< |
|
309 |
< |
/** @todo */ |
310 |
< |
void SelectionEvaluator::predefine(const std::string& script) { |
311 |
< |
|
588 |
> |
/** @todo */ |
589 |
> |
void SelectionEvaluator::predefine(const std::string& script) { |
590 |
> |
|
591 |
|
if (compiler.compile("#predefine", script)) { |
592 |
< |
std::vector<std::vector<Token> > aatoken = compiler.getAatokenCompiled(); |
593 |
< |
if (aatoken.size() != 1) { |
594 |
< |
evalError("predefinition does not have exactly 1 command:" |
595 |
< |
+ script); |
596 |
< |
return; |
597 |
< |
} |
598 |
< |
std::vector<Token> statement = aatoken[0]; |
599 |
< |
if (statement.size() > 2) { |
600 |
< |
int tok = statement[1].tok; |
601 |
< |
if (tok == Token::identifier || (tok & Token::predefinedset) == Token::predefinedset) { |
602 |
< |
std::string variable = boost::any_cast<std::string>(statement[1].value); |
603 |
< |
variables.insert(VariablesType::value_type(variable, statement)); |
604 |
< |
|
605 |
< |
} else { |
606 |
< |
evalError("invalid variable name:" + script); |
607 |
< |
} |
608 |
< |
}else { |
609 |
< |
evalError("bad predefinition length:" + script); |
610 |
< |
} |
611 |
< |
|
592 |
> |
std::vector<std::vector<Token> > aatoken = compiler.getAatokenCompiled(); |
593 |
> |
if (aatoken.size() != 1) { |
594 |
> |
evalError("predefinition does not have exactly 1 command:" |
595 |
> |
+ script); |
596 |
> |
return; |
597 |
> |
} |
598 |
> |
std::vector<Token> statement = aatoken[0]; |
599 |
> |
if (statement.size() > 2) { |
600 |
> |
int tok = statement[1].tok; |
601 |
> |
if (tok == Token::identifier || |
602 |
> |
(tok & Token::predefinedset) == Token::predefinedset) { |
603 |
> |
std::string variable = boost::any_cast<std::string>(statement[1].value); |
604 |
> |
variables.insert(VariablesType::value_type(variable, statement)); |
605 |
> |
|
606 |
> |
} else { |
607 |
> |
evalError("invalid variable name:" + script); |
608 |
> |
} |
609 |
> |
}else { |
610 |
> |
evalError("bad predefinition length:" + script); |
611 |
> |
} |
612 |
|
|
613 |
|
} else { |
614 |
< |
evalError("predefined set compile error:" + script + |
615 |
< |
"\ncompile error:" + compiler.getErrorMessage()); |
614 |
> |
evalError("predefined set compile error:" + script + |
615 |
> |
"\ncompile error:" + compiler.getErrorMessage()); |
616 |
|
} |
617 |
+ |
} |
618 |
|
|
619 |
< |
} |
340 |
< |
|
341 |
< |
void SelectionEvaluator::select(BitSet& bs){ |
619 |
> |
void SelectionEvaluator::select(SelectionSet& bs){ |
620 |
|
bs = expression(statement, 1); |
621 |
< |
} |
621 |
> |
} |
622 |
|
|
623 |
< |
BitSet SelectionEvaluator::lookupValue(const std::string& variable){ |
624 |
< |
|
625 |
< |
BitSet bs(nStuntDouble); |
623 |
> |
void SelectionEvaluator::select(SelectionSet& bs, int frame){ |
624 |
> |
bs = expression(statement, 1, frame); |
625 |
> |
} |
626 |
> |
|
627 |
> |
SelectionSet SelectionEvaluator::lookupValue(const std::string& variable){ |
628 |
> |
|
629 |
> |
SelectionSet bs = createSelectionSets(); |
630 |
|
std::map<std::string, boost::any>::iterator i = variables.find(variable); |
631 |
|
|
632 |
|
if (i != variables.end()) { |
633 |
< |
if (i->second.type() == typeid(BitSet)) { |
634 |
< |
return boost::any_cast<BitSet>(i->second); |
635 |
< |
} else if (i->second.type() == typeid(std::vector<Token>)){ |
636 |
< |
bs = expression(boost::any_cast<std::vector<Token> >(i->second), 2); |
637 |
< |
i->second = bs; /**@todo fixme */ |
638 |
< |
return bs; |
639 |
< |
} |
633 |
> |
if (i->second.type() == typeid(SelectionSet)) { |
634 |
> |
return boost::any_cast<SelectionSet>(i->second); |
635 |
> |
} else if (i->second.type() == typeid(std::vector<Token>)){ |
636 |
> |
bs = expression(boost::any_cast<std::vector<Token> >(i->second), 2); |
637 |
> |
i->second = bs; /**@todo fixme */ |
638 |
> |
return bs; |
639 |
> |
} |
640 |
|
} else { |
641 |
< |
unrecognizedIdentifier(variable); |
641 |
> |
unrecognizedIdentifier(variable); |
642 |
|
} |
643 |
< |
|
643 |
> |
|
644 |
|
return bs; |
645 |
< |
} |
645 |
> |
} |
646 |
> |
|
647 |
> |
SelectionSet SelectionEvaluator::nameInstruction(const std::string& name){ |
648 |
> |
return nameFinder.match(name); |
649 |
> |
} |
650 |
|
|
651 |
< |
BitSet SelectionEvaluator::nameInstruction(const std::string& name){ |
366 |
< |
|
367 |
< |
return nameFinder.match(name); |
368 |
< |
|
369 |
< |
} |
370 |
< |
|
371 |
< |
bool SelectionEvaluator::containDynamicToken(const std::vector<Token>& tokens){ |
651 |
> |
bool SelectionEvaluator::containDynamicToken(const std::vector<Token>& tokens){ |
652 |
|
std::vector<Token>::const_iterator i; |
653 |
|
for (i = tokens.begin(); i != tokens.end(); ++i) { |
654 |
< |
if (i->tok & Token::dynamic) { |
655 |
< |
return true; |
656 |
< |
} |
654 |
> |
if (i->tok & Token::dynamic) { |
655 |
> |
return true; |
656 |
> |
} |
657 |
|
} |
658 |
< |
|
658 |
> |
|
659 |
|
return false; |
660 |
< |
} |
660 |
> |
} |
661 |
|
|
662 |
< |
void SelectionEvaluator::clearDefinitionsAndLoadPredefined() { |
662 |
> |
void SelectionEvaluator::clearDefinitionsAndLoadPredefined() { |
663 |
|
variables.clear(); |
664 |
|
//load predefine |
665 |
|
//predefine(); |
666 |
< |
} |
666 |
> |
} |
667 |
|
|
668 |
< |
BitSet SelectionEvaluator::evaluate() { |
669 |
< |
BitSet bs(nStuntDouble); |
668 |
> |
SelectionSet SelectionEvaluator::createSelectionSets() { |
669 |
> |
SelectionSet ss(nObjects); |
670 |
> |
return ss; |
671 |
> |
} |
672 |
> |
|
673 |
> |
SelectionSet SelectionEvaluator::evaluate() { |
674 |
> |
SelectionSet bs = createSelectionSets(); |
675 |
|
if (isLoaded_) { |
676 |
< |
pc = 0; |
677 |
< |
instructionDispatchLoop(bs); |
676 |
> |
pc = 0; |
677 |
> |
instructionDispatchLoop(bs); |
678 |
|
} |
679 |
+ |
return bs; |
680 |
+ |
} |
681 |
|
|
682 |
+ |
SelectionSet SelectionEvaluator::evaluate(int frame) { |
683 |
+ |
SelectionSet bs = createSelectionSets(); |
684 |
+ |
if (isLoaded_) { |
685 |
+ |
pc = 0; |
686 |
+ |
instructionDispatchLoop(bs, frame); |
687 |
+ |
} |
688 |
|
return bs; |
689 |
< |
} |
689 |
> |
} |
690 |
|
|
691 |
< |
BitSet SelectionEvaluator::indexInstruction(const boost::any& value) { |
692 |
< |
BitSet bs(nStuntDouble); |
691 |
> |
SelectionSet SelectionEvaluator::indexInstruction(const boost::any& value) { |
692 |
> |
SelectionSet bs = createSelectionSets(); |
693 |
|
|
694 |
|
if (value.type() == typeid(int)) { |
695 |
< |
int index = boost::any_cast<int>(value); |
696 |
< |
if (index < 0 || index >= bs.size()) { |
697 |
< |
invalidIndex(index); |
698 |
< |
} else { |
699 |
< |
indexFinder.find(index); |
700 |
< |
} |
695 |
> |
int index = boost::any_cast<int>(value); |
696 |
> |
if (index < 0 || index >= bs.bitsets_[STUNTDOUBLE].size()) { |
697 |
> |
invalidIndex(index); |
698 |
> |
} else { |
699 |
> |
bs = indexFinder.find(index); |
700 |
> |
} |
701 |
|
} else if (value.type() == typeid(std::pair<int, int>)) { |
702 |
< |
std::pair<int, int> indexRange= boost::any_cast<std::pair<int, int> >(value); |
703 |
< |
assert(indexRange.first <= indexRange.second); |
704 |
< |
if (indexRange.first < 0 || indexRange.second >= bs.size()) { |
705 |
< |
invalidIndexRange(indexRange); |
706 |
< |
}else { |
707 |
< |
indexFinder.find(indexRange.first, indexRange.second); |
708 |
< |
} |
702 |
> |
std::pair<int, int> indexRange= boost::any_cast<std::pair<int, int> >(value); |
703 |
> |
assert(indexRange.first <= indexRange.second); |
704 |
> |
if (indexRange.first < 0 || |
705 |
> |
indexRange.second >= bs.bitsets_[STUNTDOUBLE].size()) { |
706 |
> |
invalidIndexRange(indexRange); |
707 |
> |
}else { |
708 |
> |
bs = indexFinder.find(indexRange.first, indexRange.second); |
709 |
> |
} |
710 |
|
} |
711 |
|
|
712 |
|
return bs; |
713 |
< |
} |
713 |
> |
} |
714 |
|
|
715 |
+ |
SelectionSet SelectionEvaluator::allInstruction() { |
716 |
+ |
SelectionSet ss = createSelectionSets(); |
717 |
+ |
|
718 |
+ |
SimInfo::MoleculeIterator mi; |
719 |
+ |
Molecule::AtomIterator ai; |
720 |
+ |
Molecule::RigidBodyIterator rbIter; |
721 |
+ |
Molecule::BondIterator bondIter; |
722 |
+ |
Molecule::BendIterator bendIter; |
723 |
+ |
Molecule::TorsionIterator torsionIter; |
724 |
+ |
Molecule::InversionIterator inversionIter; |
725 |
+ |
|
726 |
+ |
Molecule* mol; |
727 |
+ |
Atom* atom; |
728 |
+ |
RigidBody* rb; |
729 |
+ |
Bond* bond; |
730 |
+ |
Bend* bend; |
731 |
+ |
Torsion* torsion; |
732 |
+ |
Inversion* inversion; |
733 |
+ |
|
734 |
+ |
// Doing the loop insures that we're actually on this processor. |
735 |
+ |
|
736 |
+ |
for (mol = info->beginMolecule(mi); mol != NULL; |
737 |
+ |
mol = info->nextMolecule(mi)) { |
738 |
+ |
for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) { |
739 |
+ |
ss.bitsets_[STUNTDOUBLE].setBitOn(atom->getGlobalIndex()); |
740 |
+ |
} |
741 |
+ |
for (rb = mol->beginRigidBody(rbIter); rb != NULL; |
742 |
+ |
rb = mol->nextRigidBody(rbIter)) { |
743 |
+ |
ss.bitsets_[STUNTDOUBLE].setBitOn(rb->getGlobalIndex()); |
744 |
+ |
} |
745 |
+ |
for (bond = mol->beginBond(bondIter); bond != NULL; |
746 |
+ |
bond = mol->nextBond(bondIter)) { |
747 |
+ |
ss.bitsets_[BOND].setBitOn(bond->getGlobalIndex()); |
748 |
+ |
} |
749 |
+ |
for (bend = mol->beginBend(bendIter); bend != NULL; |
750 |
+ |
bend = mol->nextBend(bendIter)) { |
751 |
+ |
ss.bitsets_[BEND].setBitOn(bend->getGlobalIndex()); |
752 |
+ |
} |
753 |
+ |
for (torsion = mol->beginTorsion(torsionIter); torsion != NULL; |
754 |
+ |
torsion = mol->nextTorsion(torsionIter)) { |
755 |
+ |
ss.bitsets_[TORSION].setBitOn(torsion->getGlobalIndex()); |
756 |
+ |
} |
757 |
+ |
for (inversion = mol->beginInversion(inversionIter); inversion != NULL; |
758 |
+ |
inversion = mol->nextInversion(inversionIter)) { |
759 |
+ |
ss.bitsets_[INVERSION].setBitOn(inversion->getGlobalIndex()); |
760 |
+ |
} |
761 |
+ |
ss.bitsets_[MOLECULE].setBitOn(mol->getGlobalIndex()); |
762 |
+ |
} |
763 |
+ |
|
764 |
+ |
return ss; |
765 |
+ |
} |
766 |
+ |
|
767 |
+ |
SelectionSet SelectionEvaluator::hull() { |
768 |
+ |
SelectionSet bs = createSelectionSets(); |
769 |
+ |
|
770 |
+ |
bs = hullFinder.findHull(); |
771 |
+ |
surfaceArea_ = hullFinder.getSurfaceArea(); |
772 |
+ |
hasSurfaceArea_ = true; |
773 |
+ |
return bs; |
774 |
+ |
} |
775 |
+ |
|
776 |
+ |
|
777 |
+ |
SelectionSet SelectionEvaluator::hull(int frame) { |
778 |
+ |
SelectionSet bs = createSelectionSets(); |
779 |
+ |
|
780 |
+ |
bs = hullFinder.findHull(frame); |
781 |
+ |
|
782 |
+ |
return bs; |
783 |
+ |
} |
784 |
+ |
|
785 |
+ |
RealType SelectionEvaluator::getCharge(Atom* atom) { |
786 |
+ |
RealType charge = 0.0; |
787 |
+ |
AtomType* atomType = atom->getAtomType(); |
788 |
+ |
|
789 |
+ |
FixedChargeAdapter fca = FixedChargeAdapter(atomType); |
790 |
+ |
if ( fca.isFixedCharge() ) { |
791 |
+ |
charge = fca.getCharge(); |
792 |
+ |
} |
793 |
+ |
|
794 |
+ |
FluctuatingChargeAdapter fqa = FluctuatingChargeAdapter(atomType); |
795 |
+ |
if ( fqa.isFluctuatingCharge() ) { |
796 |
+ |
charge += atom->getFlucQPos(); |
797 |
+ |
} |
798 |
+ |
return charge; |
799 |
+ |
} |
800 |
+ |
|
801 |
+ |
RealType SelectionEvaluator::getCharge(Atom* atom, int frame) { |
802 |
+ |
RealType charge = 0.0; |
803 |
+ |
AtomType* atomType = atom->getAtomType(); |
804 |
+ |
|
805 |
+ |
FixedChargeAdapter fca = FixedChargeAdapter(atomType); |
806 |
+ |
if ( fca.isFixedCharge() ) { |
807 |
+ |
charge = fca.getCharge(); |
808 |
+ |
} |
809 |
+ |
|
810 |
+ |
FluctuatingChargeAdapter fqa = FluctuatingChargeAdapter(atomType); |
811 |
+ |
if ( fqa.isFluctuatingCharge() ) { |
812 |
+ |
charge += atom->getFlucQPos(frame); |
813 |
+ |
} |
814 |
+ |
return charge; |
815 |
+ |
} |
816 |
+ |
|
817 |
|
} |