ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/selection/SelectionEvaluator.cpp
Revision: 1629
Committed: Wed Sep 14 21:15:17 2011 UTC (13 years, 7 months ago) by gezelter
File size: 14212 byte(s)
Log Message:
Merging changes from old branch into development branch

File Contents

# User Rev Content
1 tim 277 /*
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 gezelter 1390 * 1. Redistributions of source code must retain the above copyright
10 tim 277 * notice, this list of conditions and the following disclaimer.
11     *
12 gezelter 1390 * 2. Redistributions in binary form must reproduce the above copyright
13 tim 277 * 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 gezelter 1390 *
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, 24107 (2008).
39     * [4] Vardeman & Gezelter, in progress (2009).
40 tim 277 */
41    
42 tim 283 #include <stack>
43 tim 277 #include "selection/SelectionEvaluator.hpp"
44 tim 283 #include "primitives/Atom.hpp"
45     #include "primitives/DirectionalAtom.hpp"
46     #include "primitives/RigidBody.hpp"
47     #include "primitives/Molecule.hpp"
48 gezelter 1627 #include "io/ifstrstream.hpp"
49 tim 283
50 gezelter 1390 namespace OpenMD {
51 tim 277
52    
53 gezelter 507 SelectionEvaluator::SelectionEvaluator(SimInfo* si)
54 gezelter 1412 : info(si), nameFinder(info), distanceFinder(info), hullFinder(info),
55     indexFinder(info),
56 cli2 1364 isLoaded_(false){
57     nStuntDouble = info->getNGlobalAtoms() + info->getNGlobalRigidBodies();
58 gezelter 507 }
59 tim 283
60 cli2 1364 bool SelectionEvaluator::loadScript(const std::string& filename,
61     const std::string& script) {
62 tim 288 clearDefinitionsAndLoadPredefined();
63 tim 278 this->filename = filename;
64     this->script = script;
65     if (! compiler.compile(filename, script)) {
66 gezelter 507 error = true;
67     errorMessage = compiler.getErrorMessage();
68 cli2 1364
69     sprintf( painCave.errMsg,
70     "SelectionCompiler Error: %s\n", errorMessage.c_str());
71 gezelter 1390 painCave.severity = OPENMD_ERROR;
72 cli2 1364 painCave.isFatal = 1;
73     simError();
74 gezelter 507 return false;
75 tim 278 }
76    
77     pc = 0;
78     aatoken = compiler.getAatokenCompiled();
79     linenumbers = compiler.getLineNumbers();
80     lineIndices = compiler.getLineIndices();
81 tim 288
82     std::vector<std::vector<Token> >::const_iterator i;
83    
84     isDynamic_ = false;
85     for (i = aatoken.begin(); i != aatoken.end(); ++i) {
86 gezelter 507 if (containDynamicToken(*i)) {
87     isDynamic_ = true;
88     break;
89     }
90 tim 288 }
91    
92     isLoaded_ = true;
93 tim 278 return true;
94 gezelter 507 }
95 tim 278
96 gezelter 507 void SelectionEvaluator::clearState() {
97 tim 278 error = false;
98 tim 281 errorMessage = "";
99 gezelter 507 }
100 tim 278
101 gezelter 507 bool SelectionEvaluator::loadScriptString(const std::string& script) {
102 tim 278 clearState();
103 tim 281 return loadScript("", script);
104 gezelter 507 }
105 tim 278
106 gezelter 507 bool SelectionEvaluator::loadScriptFile(const std::string& filename) {
107 tim 278 clearState();
108     return loadScriptFileInternal(filename);
109 gezelter 507 }
110 tim 278
111 cli2 1364 bool SelectionEvaluator::loadScriptFileInternal(const std::string & filename) {
112     ifstrstream ifs(filename.c_str());
113 tim 288 if (!ifs.is_open()) {
114 gezelter 507 return false;
115 tim 288 }
116 cli2 1364
117 tim 288 const int bufferSize = 65535;
118     char buffer[bufferSize];
119     std::string script;
120     while(ifs.getline(buffer, bufferSize)) {
121 gezelter 507 script += buffer;
122 tim 288 }
123     return loadScript(filename, script);
124 gezelter 507 }
125 cli2 1364
126 gezelter 1390 void SelectionEvaluator::instructionDispatchLoop(OpenMDBitSet& bs){
127 tim 288
128 tim 281 while ( pc < aatoken.size()) {
129 gezelter 507 statement = aatoken[pc++];
130     statementLength = statement.size();
131     Token token = statement[0];
132     switch (token.tok) {
133     case Token::define:
134     define();
135     break;
136     case Token::select:
137     select(bs);
138     break;
139     default:
140     unrecognizedCommand(token);
141     return;
142     }
143 tim 278 }
144 tim 288
145 gezelter 507 }
146 tim 278
147 gezelter 1390 OpenMDBitSet SelectionEvaluator::expression(const std::vector<Token>& code,
148 cli2 1364 int pcStart) {
149 gezelter 1390 OpenMDBitSet bs;
150     std::stack<OpenMDBitSet> stack;
151 cli2 1364
152 tim 288 for (int pc = pcStart; pc < code.size(); ++pc) {
153 tim 278 Token instruction = code[pc];
154 tim 282
155 tim 278 switch (instruction.tok) {
156 tim 281 case Token::expressionBegin:
157 tim 278 break;
158 tim 281 case Token::expressionEnd:
159 tim 282 break;
160 tim 281 case Token::all:
161 gezelter 1390 bs = OpenMDBitSet(nStuntDouble);
162 tim 283 bs.setAll();
163     stack.push(bs);
164 tim 278 break;
165 tim 281 case Token::none:
166 gezelter 1390 bs = OpenMDBitSet(nStuntDouble);
167 tim 283 stack.push(bs);
168 tim 278 break;
169 tim 281 case Token::opOr:
170 tim 283 bs = stack.top();
171     stack.pop();
172     stack.top() |= bs;
173 tim 278 break;
174 tim 281 case Token::opAnd:
175 tim 283 bs = stack.top();
176     stack.pop();
177     stack.top() &= bs;
178 tim 278 break;
179 tim 281 case Token::opNot:
180 tim 283 stack.top().flip();
181 tim 278 break;
182 tim 281 case Token::within:
183 tim 283 withinInstruction(instruction, stack.top());
184 tim 278 break;
185 gezelter 1412 case Token::hull:
186     stack.push(hull());
187     break;
188 gezelter 507 //case Token::selected:
189     // stack.push(getSelectionSet());
190     // break;
191 tim 281 case Token::name:
192 tim 283 stack.push(nameInstruction(boost::any_cast<std::string>(instruction.value)));
193 tim 278 break;
194 tim 295 case Token::index:
195     stack.push(indexInstruction(instruction.value));
196 tim 281 break;
197     case Token::identifier:
198 tim 283 stack.push(lookupValue(boost::any_cast<std::string>(instruction.value)));
199 tim 281 break;
200     case Token::opLT:
201     case Token::opLE:
202     case Token::opGE:
203     case Token::opGT:
204     case Token::opEQ:
205     case Token::opNE:
206 tim 283 stack.push(comparatorInstruction(instruction));
207 tim 278 break;
208     default:
209     unrecognizedExpression();
210     }
211     }
212 tim 283 if (stack.size() != 1)
213 tim 278 evalError("atom expression compiler error - stack over/underflow");
214 tim 283
215     return stack.top();
216 tim 278 }
217    
218    
219    
220 gezelter 1390 OpenMDBitSet SelectionEvaluator::comparatorInstruction(const Token& instruction) {
221 tim 278 int comparator = instruction.tok;
222     int property = instruction.intValue;
223 tim 283 float comparisonValue = boost::any_cast<float>(instruction.value);
224     float propertyValue;
225 gezelter 1390 OpenMDBitSet bs(nStuntDouble);
226 tim 295 bs.clearAll();
227 tim 283
228     SimInfo::MoleculeIterator mi;
229     Molecule* mol;
230     Molecule::AtomIterator ai;
231     Atom* atom;
232     Molecule::RigidBodyIterator rbIter;
233     RigidBody* rb;
234 tim 281
235 cli2 1364 for (mol = info->beginMolecule(mi); mol != NULL;
236     mol = info->nextMolecule(mi)) {
237    
238 gezelter 507 for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
239     compareProperty(atom, bs, property, comparator, comparisonValue);
240     }
241 cli2 1364
242     for (rb = mol->beginRigidBody(rbIter); rb != NULL;
243     rb = mol->nextRigidBody(rbIter)) {
244     compareProperty(rb, bs, property, comparator, comparisonValue);
245     }
246 tim 278 }
247    
248 tim 283 return bs;
249 gezelter 507 }
250 tim 278
251 gezelter 1390 void SelectionEvaluator::compareProperty(StuntDouble* sd, OpenMDBitSet& bs,
252 cli2 1364 int property, int comparator,
253     float comparisonValue) {
254 tim 963 RealType propertyValue = 0.0;
255 gezelter 507 switch (property) {
256     case Token::mass:
257     propertyValue = sd->getMass();
258     break;
259     case Token::charge:
260     if (sd->isAtom()){
261     Atom* atom = static_cast<Atom*>(sd);
262     propertyValue = getCharge(atom);
263     } else if (sd->isRigidBody()) {
264     RigidBody* rb = static_cast<RigidBody*>(sd);
265     RigidBody::AtomIterator ai;
266     Atom* atom;
267     for (atom = rb->beginAtom(ai); atom != NULL; atom = rb->nextAtom(ai)) {
268     propertyValue+= getCharge(atom);
269     }
270     }
271     break;
272 cli2 1360 case Token::x:
273     propertyValue = sd->getPos().x();
274     break;
275     case Token::y:
276     propertyValue = sd->getPos().y();
277     break;
278     case Token::z:
279     propertyValue = sd->getPos().z();
280     break;
281 gezelter 1629 case Token::r:
282     propertyValue = sd->getPos().length();
283     break;
284 gezelter 507 default:
285     unrecognizedAtomProperty(property);
286     }
287 tim 283
288 gezelter 507 bool match = false;
289     switch (comparator) {
290     case Token::opLT:
291     match = propertyValue < comparisonValue;
292     break;
293     case Token::opLE:
294     match = propertyValue <= comparisonValue;
295     break;
296     case Token::opGE:
297     match = propertyValue >= comparisonValue;
298     break;
299     case Token::opGT:
300     match = propertyValue > comparisonValue;
301     break;
302     case Token::opEQ:
303     match = propertyValue == comparisonValue;
304     break;
305     case Token::opNE:
306     match = propertyValue != comparisonValue;
307     break;
308     }
309     if (match)
310     bs.setBitOn(sd->getGlobalIndex());
311 cli2 1364
312 gezelter 507 }
313 tim 283
314 cli2 1364 void SelectionEvaluator::withinInstruction(const Token& instruction,
315 gezelter 1390 OpenMDBitSet& bs){
316 tim 295
317 tim 281 boost::any withinSpec = instruction.value;
318 tim 295 float distance;
319 tim 282 if (withinSpec.type() == typeid(float)){
320 gezelter 507 distance = boost::any_cast<float>(withinSpec);
321 tim 295 } else if (withinSpec.type() == typeid(int)) {
322 gezelter 507 distance = boost::any_cast<int>(withinSpec);
323 tim 295 } else {
324 gezelter 507 evalError("casting error in withinInstruction");
325     bs.clearAll();
326 tim 278 }
327 tim 282
328 tim 295 bs = distanceFinder.find(bs, distance);
329 gezelter 507 }
330 cli2 1364
331 gezelter 507 void SelectionEvaluator::define() {
332 tim 282 assert(statement.size() >= 3);
333 cli2 1364
334 tim 282 std::string variable = boost::any_cast<std::string>(statement[1].value);
335 cli2 1364
336     variables.insert(VariablesType::value_type(variable,
337     expression(statement, 2)));
338 gezelter 507 }
339 cli2 1364
340 tim 278
341 gezelter 507 /** @todo */
342     void SelectionEvaluator::predefine(const std::string& script) {
343 cli2 1364
344 tim 282 if (compiler.compile("#predefine", script)) {
345 gezelter 507 std::vector<std::vector<Token> > aatoken = compiler.getAatokenCompiled();
346     if (aatoken.size() != 1) {
347     evalError("predefinition does not have exactly 1 command:"
348     + script);
349     return;
350     }
351     std::vector<Token> statement = aatoken[0];
352     if (statement.size() > 2) {
353     int tok = statement[1].tok;
354 cli2 1364 if (tok == Token::identifier ||
355     (tok & Token::predefinedset) == Token::predefinedset) {
356 gezelter 507 std::string variable = boost::any_cast<std::string>(statement[1].value);
357     variables.insert(VariablesType::value_type(variable, statement));
358 cli2 1364
359 gezelter 507 } else {
360     evalError("invalid variable name:" + script);
361     }
362     }else {
363     evalError("bad predefinition length:" + script);
364 cli2 1364 }
365 tim 282
366     } else {
367 gezelter 507 evalError("predefined set compile error:" + script +
368     "\ncompile error:" + compiler.getErrorMessage());
369 tim 282 }
370 gezelter 507 }
371 tim 282
372 gezelter 1390 void SelectionEvaluator::select(OpenMDBitSet& bs){
373 tim 288 bs = expression(statement, 1);
374 gezelter 507 }
375 cli2 1364
376 gezelter 1390 OpenMDBitSet SelectionEvaluator::lookupValue(const std::string& variable){
377 cli2 1364
378 gezelter 1390 OpenMDBitSet bs(nStuntDouble);
379 tim 282 std::map<std::string, boost::any>::iterator i = variables.find(variable);
380 tim 295
381 tim 282 if (i != variables.end()) {
382 gezelter 1390 if (i->second.type() == typeid(OpenMDBitSet)) {
383     return boost::any_cast<OpenMDBitSet>(i->second);
384 gezelter 507 } else if (i->second.type() == typeid(std::vector<Token>)){
385     bs = expression(boost::any_cast<std::vector<Token> >(i->second), 2);
386     i->second = bs; /**@todo fixme */
387     return bs;
388     }
389 tim 283 } else {
390 gezelter 507 unrecognizedIdentifier(variable);
391 tim 282 }
392 cli2 1364
393 tim 295 return bs;
394 gezelter 507 }
395 cli2 1364
396 gezelter 1390 OpenMDBitSet SelectionEvaluator::nameInstruction(const std::string& name){
397 cli2 1364 return nameFinder.match(name);
398 gezelter 507 }
399 tim 283
400 gezelter 507 bool SelectionEvaluator::containDynamicToken(const std::vector<Token>& tokens){
401 tim 288 std::vector<Token>::const_iterator i;
402     for (i = tokens.begin(); i != tokens.end(); ++i) {
403 gezelter 507 if (i->tok & Token::dynamic) {
404     return true;
405     }
406 tim 288 }
407 cli2 1364
408 tim 288 return false;
409 gezelter 507 }
410 tim 288
411 gezelter 507 void SelectionEvaluator::clearDefinitionsAndLoadPredefined() {
412 tim 288 variables.clear();
413     //load predefine
414     //predefine();
415 gezelter 507 }
416 tim 288
417 gezelter 1390 OpenMDBitSet SelectionEvaluator::evaluate() {
418     OpenMDBitSet bs(nStuntDouble);
419 tim 288 if (isLoaded_) {
420 gezelter 507 pc = 0;
421     instructionDispatchLoop(bs);
422 tim 288 }
423    
424     return bs;
425 gezelter 507 }
426 tim 295
427 gezelter 1390 OpenMDBitSet SelectionEvaluator::indexInstruction(const boost::any& value) {
428     OpenMDBitSet bs(nStuntDouble);
429 tim 295
430     if (value.type() == typeid(int)) {
431 gezelter 507 int index = boost::any_cast<int>(value);
432     if (index < 0 || index >= bs.size()) {
433     invalidIndex(index);
434     } else {
435     bs = indexFinder.find(index);
436     }
437 tim 295 } else if (value.type() == typeid(std::pair<int, int>)) {
438 gezelter 507 std::pair<int, int> indexRange= boost::any_cast<std::pair<int, int> >(value);
439     assert(indexRange.first <= indexRange.second);
440     if (indexRange.first < 0 || indexRange.second >= bs.size()) {
441     invalidIndexRange(indexRange);
442     }else {
443     bs = indexFinder.find(indexRange.first, indexRange.second);
444     }
445 tim 295 }
446    
447     return bs;
448 gezelter 507 }
449 tim 295
450 tim 432
451 gezelter 1412 OpenMDBitSet SelectionEvaluator::hull() {
452     OpenMDBitSet bs(nStuntDouble);
453    
454     bs = hullFinder.findHull();
455    
456     return bs;
457     }
458    
459    
460    
461 tim 963 RealType SelectionEvaluator::getCharge(Atom* atom) {
462     RealType charge =0.0;
463 tim 432 AtomType* atomType = atom->getAtomType();
464     if (atomType->isCharge()) {
465 gezelter 507 GenericData* data = atomType->getPropertyByName("Charge");
466     if (data != NULL) {
467     DoubleGenericData* doubleData= dynamic_cast<DoubleGenericData*>(data);
468 tim 432
469 gezelter 507 if (doubleData != NULL) {
470     charge = doubleData->getData();
471 tim 432
472 gezelter 507 } else {
473     sprintf( painCave.errMsg,
474     "Can not cast GenericData to DoubleGenericData\n");
475 gezelter 1390 painCave.severity = OPENMD_ERROR;
476 gezelter 507 painCave.isFatal = 1;
477     simError();
478     }
479     }
480 tim 432 }
481    
482     return charge;
483 gezelter 507 }
484 tim 432
485     }

Properties

Name Value
svn:keywords Author Id Revision Date