ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/selection/SelectionEvaluator.cpp
(Generate patch)

Comparing trunk/src/selection/SelectionEvaluator.cpp (file contents):
Revision 413 by tim, Wed Mar 9 17:30:29 2005 UTC vs.
Revision 1513 by gezelter, Tue Oct 19 18:40:54 2010 UTC

# Line 6 | Line 6
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.
# Line 37 | Line 28
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, 24107 (2008).          
39 + * [4]  Vardeman & Gezelter, in progress (2009).                        
40   */
41  
42   #include <stack>
# Line 45 | Line 45
45   #include "primitives/DirectionalAtom.hpp"
46   #include "primitives/RigidBody.hpp"
47   #include "primitives/Molecule.hpp"
48 + #include "io/basic_ifstrstream.hpp"
49  
50 < namespace oopse {
50 > namespace OpenMD {
51  
52  
53 < SelectionEvaluator::SelectionEvaluator(SimInfo* si)
54 <    : info(si), nameFinder(info), distanceFinder(info), indexFinder(info), isLoaded_(false){
55 <    
56 <    nStuntDouble = info->getNGlobalAtoms() + info->getNRigidBodies();
57 < }            
53 >  SelectionEvaluator::SelectionEvaluator(SimInfo* si)
54 >    : info(si), nameFinder(info), distanceFinder(info), hullFinder(info),
55 >      indexFinder(info),
56 >      isLoaded_(false){    
57 >      nStuntDouble = info->getNGlobalAtoms() + info->getNGlobalRigidBodies();
58 >    }            
59  
60 < bool SelectionEvaluator::loadScript(const std::string& filename, const std::string& script) {
60 >  bool SelectionEvaluator::loadScript(const std::string& filename,
61 >                                      const std::string& script) {
62      clearDefinitionsAndLoadPredefined();
63      this->filename = filename;
64      this->script = script;
65      if (! compiler.compile(filename, script)) {
66 <        error = true;
67 <        errorMessage = compiler.getErrorMessage();
68 <        std::cerr << "SelectionCompiler Error: " << errorMessage << std::endl;
69 <        return false;
66 >      error = true;
67 >      errorMessage = compiler.getErrorMessage();
68 >
69 >      sprintf( painCave.errMsg,
70 >               "SelectionCompiler Error: %s\n", errorMessage.c_str());
71 >      painCave.severity = OPENMD_ERROR;
72 >      painCave.isFatal = 1;
73 >      simError();
74 >      return false;
75      }
76  
77      pc = 0;
# Line 75 | Line 83 | bool SelectionEvaluator::loadScript(const std::string&
83  
84      isDynamic_ = false;
85      for (i = aatoken.begin(); i != aatoken.end(); ++i) {
86 <        if (containDynamicToken(*i)) {
87 <            isDynamic_ = true;
88 <            break;
89 <        }
86 >      if (containDynamicToken(*i)) {
87 >        isDynamic_ = true;
88 >        break;
89 >      }
90      }
91  
92      isLoaded_ = true;
93      return true;
94 < }
94 >  }
95  
96 < void SelectionEvaluator::clearState() {
89 <    //for (int i = scriptLevelMax; --i >= 0; )
90 <    //    stack[i].clear();
91 <    //scriptLevel = 0;
96 >  void SelectionEvaluator::clearState() {
97      error = false;
98      errorMessage = "";
99 < }
99 >  }
100  
101 < bool SelectionEvaluator::loadScriptString(const std::string& script) {
101 >  bool SelectionEvaluator::loadScriptString(const std::string& script) {
102      clearState();
103      return loadScript("", script);
104 < }
104 >  }
105  
106 < bool SelectionEvaluator::loadScriptFile(const std::string& filename) {
106 >  bool SelectionEvaluator::loadScriptFile(const std::string& filename) {
107      clearState();
108      return loadScriptFileInternal(filename);
109 < }
109 >  }
110  
111 < bool SelectionEvaluator::loadScriptFileInternal(const  std::string & filename) {
112 <  std::ifstream ifs(filename.c_str());
111 >  bool SelectionEvaluator::loadScriptFileInternal(const std::string & filename) {
112 >    ifstrstream ifs(filename.c_str());
113      if (!ifs.is_open()) {
114 <        return false;
114 >      return false;
115      }
116 <
116 >    
117      const int bufferSize = 65535;
118      char buffer[bufferSize];
119      std::string script;
120      while(ifs.getline(buffer, bufferSize)) {
121 <        script += buffer;
121 >      script += buffer;
122      }
123      return loadScript(filename, script);
124 < }
125 <
126 < void SelectionEvaluator::instructionDispatchLoop(BitSet& bs){
124 >  }
125 >  
126 >  void SelectionEvaluator::instructionDispatchLoop(OpenMDBitSet& bs){
127      
128      while ( pc < aatoken.size()) {
129 <        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 <        }
129 >      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      }
144  
145 < }
145 >  }
146  
147 <  BitSet SelectionEvaluator::expression(const std::vector<Token>& code, int pcStart) {
148 <    BitSet bs;
149 <    std::stack<BitSet> stack;
150 <    
147 >  OpenMDBitSet SelectionEvaluator::expression(const std::vector<Token>& code,
148 >                                             int pcStart) {
149 >    OpenMDBitSet bs;
150 >    std::stack<OpenMDBitSet> stack;
151 >  
152      for (int pc = pcStart; pc < code.size(); ++pc) {
153        Token instruction = code[pc];
154  
# Line 152 | Line 158 | void SelectionEvaluator::instructionDispatchLoop(BitSe
158        case Token::expressionEnd:
159          break;
160        case Token::all:
161 <        bs = BitSet(nStuntDouble);
161 >        bs = OpenMDBitSet(nStuntDouble);
162          bs.setAll();
163          stack.push(bs);            
164          break;
165        case Token::none:
166 <        bs = BitSet(nStuntDouble);
166 >        bs = OpenMDBitSet(nStuntDouble);
167          stack.push(bs);            
168          break;
169        case Token::opOr:
# Line 174 | Line 180 | void SelectionEvaluator::instructionDispatchLoop(BitSe
180          stack.top().flip();
181          break;
182        case Token::within:
177
183          withinInstruction(instruction, stack.top());
184          break;
185 <      //case Token::selected:
186 <      //  stack.push(getSelectionSet());
187 <      //  break;
185 >      case Token::hull:
186 >        stack.push(hull());
187 >        break;
188 >        //case Token::selected:
189 >        //  stack.push(getSelectionSet());
190 >        //  break;
191        case Token::name:
192          stack.push(nameInstruction(boost::any_cast<std::string>(instruction.value)));
193          break;
# Line 209 | Line 217 | void SelectionEvaluator::instructionDispatchLoop(BitSe
217  
218  
219  
220 < BitSet SelectionEvaluator::comparatorInstruction(const Token& instruction) {
220 >  OpenMDBitSet SelectionEvaluator::comparatorInstruction(const Token& instruction) {
221      int comparator = instruction.tok;
222      int property = instruction.intValue;
223      float comparisonValue = boost::any_cast<float>(instruction.value);
224      float propertyValue;
225 <    BitSet bs(nStuntDouble);
225 >    OpenMDBitSet bs(nStuntDouble);
226      bs.clearAll();
227      
228      SimInfo::MoleculeIterator mi;
# Line 223 | Line 231 | BitSet SelectionEvaluator::comparatorInstruction(const
231      Atom* atom;
232      Molecule::RigidBodyIterator rbIter;
233      RigidBody* rb;
226    
227    for (mol = info->beginMolecule(mi); mol != NULL; mol = info->nextMolecule(mi)) {
234  
235 <        for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
236 <            compareProperty(atom, bs, property, comparator, comparisonValue);
237 <        }
238 <        
239 <        //change the positions of atoms which belong to the rigidbodies
240 <        for (rb = mol->beginRigidBody(rbIter); rb != NULL; rb = mol->nextRigidBody(rbIter)) {
241 <            compareProperty(rb, bs, property, comparator, comparisonValue);
242 <        }        
235 >    for (mol = info->beginMolecule(mi); mol != NULL;
236 >         mol = info->nextMolecule(mi)) {
237 >
238 >      for(atom = mol->beginAtom(ai); atom != NULL; atom = mol->nextAtom(ai)) {
239 >        compareProperty(atom, bs, property, comparator, comparisonValue);
240 >      }
241 >    
242 >      for (rb = mol->beginRigidBody(rbIter); rb != NULL;
243 >           rb = mol->nextRigidBody(rbIter)) {
244 >        compareProperty(rb, bs, property, comparator, comparisonValue);
245 >      }
246      }
247  
248      return bs;
249 < }
249 >  }
250  
251 < void SelectionEvaluator::compareProperty(StuntDouble* sd, BitSet& bs, int property, int comparator, float comparisonValue) {
252 <        double propertyValue;
253 <        switch (property) {
254 <        case Token::mass:
255 <            propertyValue = sd->getMass();
256 <            break;
257 <        case Token::charge:
258 <            return;
259 <            //break;
260 <        case Token::dipole:
261 <            return;
262 <            //break;
263 <        default:
264 <            unrecognizedAtomProperty(property);
265 <        }
251 >  void SelectionEvaluator::compareProperty(StuntDouble* sd, OpenMDBitSet& bs,
252 >                                           int property, int comparator,
253 >                                           float comparisonValue) {
254 >    RealType propertyValue = 0.0;
255 >    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 >    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 >    case Token::r:
282 >      propertyValue = sd->getPos().length();
283 >      break;
284 >    default:
285 >      unrecognizedAtomProperty(property);
286 >    }
287          
288 <        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 <
312 < }
288 >    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 >    
312 >  }
313  
314 < void SelectionEvaluator::withinInstruction(const Token& instruction, BitSet& bs){
314 >  void SelectionEvaluator::withinInstruction(const Token& instruction,
315 >                                             OpenMDBitSet& bs){
316      
317      boost::any withinSpec = instruction.value;
318      float distance;
319      if (withinSpec.type() == typeid(float)){
320 <        distance = boost::any_cast<float>(withinSpec);
320 >      distance = boost::any_cast<float>(withinSpec);
321      } else if (withinSpec.type() == typeid(int)) {
322 <        distance = boost::any_cast<int>(withinSpec);    
322 >      distance = boost::any_cast<int>(withinSpec);    
323      } else {
324 <        evalError("casting error in withinInstruction");
325 <        bs.clearAll();
324 >      evalError("casting error in withinInstruction");
325 >      bs.clearAll();
326      }
327      
328      bs = distanceFinder.find(bs, distance);            
329 < }
330 <
331 < void SelectionEvaluator::define() {
329 >  }
330 >  
331 >  void SelectionEvaluator::define() {
332      assert(statement.size() >= 3);
333 <
333 >    
334      std::string variable = boost::any_cast<std::string>(statement[1].value);
335 +    
336 +    variables.insert(VariablesType::value_type(variable,
337 +                                               expression(statement, 2)));
338 +  }
339 +  
340  
341 <    variables.insert(VariablesType::value_type(variable, expression(statement, 2)));
342 < }
343 <
308 <
309 < /** @todo */
310 < void SelectionEvaluator::predefine(const std::string& script) {
311 <
341 >  /** @todo */
342 >  void SelectionEvaluator::predefine(const std::string& script) {
343 >    
344      if (compiler.compile("#predefine", script)) {
345 <        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 <            if (tok == Token::identifier || (tok & Token::predefinedset) == Token::predefinedset) {
355 <                std::string variable = boost::any_cast<std::string>(statement[1].value);
356 <                variables.insert(VariablesType::value_type(variable, statement));
357 <
358 <            } else {
359 <                evalError("invalid variable name:" + script);
360 <            }
361 <        }else {
362 <            evalError("bad predefinition length:" + script);
363 <        }
364 <
345 >      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 >        if (tok == Token::identifier ||
355 >            (tok & Token::predefinedset) == Token::predefinedset) {
356 >          std::string variable = boost::any_cast<std::string>(statement[1].value);
357 >          variables.insert(VariablesType::value_type(variable, statement));
358 >          
359 >        } else {
360 >          evalError("invalid variable name:" + script);
361 >        }
362 >      }else {
363 >        evalError("bad predefinition length:" + script);
364 >      }      
365          
366      } else {
367 <        evalError("predefined set compile error:" + script +
368 <          "\ncompile error:" + compiler.getErrorMessage());
367 >      evalError("predefined set compile error:" + script +
368 >                "\ncompile error:" + compiler.getErrorMessage());
369      }
370 +  }
371  
372 < }
340 <
341 < void SelectionEvaluator::select(BitSet& bs){
372 >  void SelectionEvaluator::select(OpenMDBitSet& bs){
373      bs = expression(statement, 1);
374 < }
375 <
376 < BitSet SelectionEvaluator::lookupValue(const std::string& variable){
377 <
378 <    BitSet bs(nStuntDouble);
374 >  }
375 >  
376 >  OpenMDBitSet SelectionEvaluator::lookupValue(const std::string& variable){
377 >    
378 >    OpenMDBitSet bs(nStuntDouble);
379      std::map<std::string, boost::any>::iterator i = variables.find(variable);
380      
381      if (i != variables.end()) {
382 <        if (i->second.type() == typeid(BitSet)) {
383 <            return boost::any_cast<BitSet>(i->second);
384 <        } 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 <        }
382 >      if (i->second.type() == typeid(OpenMDBitSet)) {
383 >        return boost::any_cast<OpenMDBitSet>(i->second);
384 >      } 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      } else {
390 <        unrecognizedIdentifier(variable);
390 >      unrecognizedIdentifier(variable);
391      }
392 <
392 >    
393      return bs;
394 < }
394 >  }
395 >  
396 >  OpenMDBitSet SelectionEvaluator::nameInstruction(const std::string& name){    
397 >    return nameFinder.match(name);    
398 >  }    
399  
400 < 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){
400 >  bool SelectionEvaluator::containDynamicToken(const std::vector<Token>& tokens){
401      std::vector<Token>::const_iterator i;
402      for (i = tokens.begin(); i != tokens.end(); ++i) {
403 <        if (i->tok & Token::dynamic) {
404 <            return true;
405 <        }
403 >      if (i->tok & Token::dynamic) {
404 >        return true;
405 >      }
406      }
407 <
407 >    
408      return false;
409 < }    
409 >  }    
410  
411 < void SelectionEvaluator::clearDefinitionsAndLoadPredefined() {
411 >  void SelectionEvaluator::clearDefinitionsAndLoadPredefined() {
412      variables.clear();
413      //load predefine
414      //predefine();
415 < }
415 >  }
416  
417 < BitSet SelectionEvaluator::evaluate() {
418 <    BitSet bs(nStuntDouble);
417 >  OpenMDBitSet SelectionEvaluator::evaluate() {
418 >    OpenMDBitSet bs(nStuntDouble);
419      if (isLoaded_) {
420 <        pc = 0;
421 <        instructionDispatchLoop(bs);
420 >      pc = 0;
421 >      instructionDispatchLoop(bs);
422      }
423  
424      return bs;
425 < }
425 >  }
426  
427 < BitSet SelectionEvaluator::indexInstruction(const boost::any& value) {
428 <    BitSet bs(nStuntDouble);
427 >  OpenMDBitSet SelectionEvaluator::indexInstruction(const boost::any& value) {
428 >    OpenMDBitSet bs(nStuntDouble);
429  
430      if (value.type() == typeid(int)) {
431 <        int index = boost::any_cast<int>(value);
432 <        if (index < 0 || index >= bs.size()) {
433 <            invalidIndex(index);
434 <        } else {
435 <            indexFinder.find(index);
436 <        }
431 >      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      } else if (value.type() == typeid(std::pair<int, int>)) {
438 <        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 <            indexFinder.find(indexRange.first, indexRange.second);
444 <        }
438 >      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      }
446  
447      return bs;
448 < }
448 >  }
449  
450 +
451 +  OpenMDBitSet SelectionEvaluator::hull() {
452 +    OpenMDBitSet bs(nStuntDouble);
453 +    
454 +    bs = hullFinder.findHull();
455 +    
456 +    return bs;
457 +  }
458 +
459 +
460 +
461 +  RealType SelectionEvaluator::getCharge(Atom* atom) {
462 +    RealType charge =0.0;
463 +    AtomType* atomType = atom->getAtomType();
464 +    if (atomType->isCharge()) {
465 +      GenericData* data = atomType->getPropertyByName("Charge");
466 +      if (data != NULL) {
467 +        DoubleGenericData* doubleData= dynamic_cast<DoubleGenericData*>(data);
468 +
469 +        if (doubleData != NULL) {
470 +          charge = doubleData->getData();
471 +
472 +        } else {
473 +          sprintf( painCave.errMsg,
474 +                   "Can not cast GenericData to DoubleGenericData\n");
475 +          painCave.severity = OPENMD_ERROR;
476 +          painCave.isFatal = 1;
477 +          simError();          
478 +        }
479 +      }
480 +    }
481 +
482 +    return charge;
483 +  }
484 +
485   }

Comparing trunk/src/selection/SelectionEvaluator.cpp (property svn:keywords):
Revision 413 by tim, Wed Mar 9 17:30:29 2005 UTC vs.
Revision 1513 by gezelter, Tue Oct 19 18:40:54 2010 UTC

# Line 0 | Line 1
1 + Author Id Revision Date

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines