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