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

Comparing:
trunk/src/UseTheForce/ForceField.cpp (file contents), Revision 1275 by cli2, Fri Jul 4 20:54:29 2008 UTC vs.
branches/development/src/UseTheForce/ForceField.cpp (file contents), Revision 1629 by gezelter, Wed Sep 14 21:15:17 2011 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   /**
# Line 51 | Line 51
51   #include "UseTheForce/ForceField.hpp"
52   #include "utils/simError.h"
53   #include "utils/Tuple.hpp"
54 < #include "UseTheForce/DarkSide/atype_interface.h"
55 < #include "UseTheForce/DarkSide/fForceOptions_interface.h"
56 < #include "UseTheForce/DarkSide/switcheroo_interface.h"
57 < namespace oopse {
54 > namespace OpenMD {
55  
56    ForceField::ForceField() {
57 +
58      char* tempPath;
59      tempPath = getenv("FORCE_PARAM_PATH");
60 <
60 >    
61      if (tempPath == NULL) {
62        //convert a macro from compiler to a string in c++
63        STR_DEFINE(ffPath_, FRC_PATH );
# Line 68 | Line 66 | namespace oopse {
66      }
67    }
68  
69 <
70 <  ForceField::~ForceField() {
71 <    deleteAtypes();
72 <    deleteSwitch();
73 <  }
74 <
69 >  /**
70 >   * getAtomType by string
71 >   *
72 >   * finds the requested atom type in this force field using the string
73 >   * name of the atom type.
74 >   */
75    AtomType* ForceField::getAtomType(const std::string &at) {
76      std::vector<std::string> keys;
77      keys.push_back(at);
78      return atomTypeCont_.find(keys);
79    }
80  
81 +  /**
82 +   * getAtomType by ident
83 +   *
84 +   * finds the requested atom type in this force field using the
85 +   * integer ident instead of the string name of the atom type.
86 +   */
87 +  AtomType* ForceField::getAtomType(int ident) {  
88 +    std::string at = atypeIdentToName.find(ident)->second;
89 +    return getAtomType(at);
90 +  }
91 +
92    BondType* ForceField::getBondType(const std::string &at1,
93                                      const std::string &at2) {
94      std::vector<std::string> keys;
# Line 133 | Line 142 | namespace oopse {
142          ii++;
143        }
144  
136      // sort the foundBonds by the score:
145  
146 <      std::sort(foundBonds.begin(), foundBonds.end());
147 <      
148 <      int bestScore = foundBonds[0].first;
149 <      std::vector<std::string> theKeys = foundBonds[0].second;
150 <
151 <      std::cout << "best matching bond = " << theKeys[0] << "\t" << theKeys[1]  << "\t(score = "<< bestScore << ")\n";      
152 <      BondType* bestType = bondTypeCont_.find(theKeys);
153 <      if (bestType)
154 <        return bestType;
155 <      else {
156 <        //if no exact match found, try wild card match
157 <        return bondTypeCont_.find(keys, wildCardAtomTypeName_);      
146 >      if (foundBonds.size() > 0) {
147 >        // sort the foundBonds by the score:
148 >        std::sort(foundBonds.begin(), foundBonds.end());
149 >    
150 >        int bestScore = foundBonds[0].first;
151 >        std::vector<std::string> theKeys = foundBonds[0].second;
152 >        
153 >        BondType* bestType = bondTypeCont_.find(theKeys);
154 >        
155 >        return bestType;
156 >      } else {
157 >        //if no exact match found, try wild card match
158 >        return bondTypeCont_.find(keys, wildCardAtomTypeName_);      
159        }
160      }
161    }
# Line 220 | Line 229 | namespace oopse {
229          jj++;
230        }
231        
232 <      std::sort(foundBends.begin(), foundBends.end());
233 <
234 <      int jscore = foundBends[0].first;
235 <      int ikscore = foundBends[0].second;
236 <      std::vector<std::string> theKeys = foundBends[0].third;
237 <
238 <      std::cout << "best matching bend = " << theKeys[0] << "\t" <<theKeys[1]  << "\t" << theKeys[2] << "\t(scores = "<< jscore << "\t" << ikscore << ")\n";      
239 <
240 <      BendType* bestType = bendTypeCont_.find(theKeys);  
232 <      if (bestType)
233 <        return bestType;
234 <      else {
235 <      
232 >      if (foundBends.size() > 0) {
233 >        std::sort(foundBends.begin(), foundBends.end());
234 >        int jscore = foundBends[0].first;
235 >        int ikscore = foundBends[0].second;
236 >        std::vector<std::string> theKeys = foundBends[0].third;      
237 >        
238 >        BendType* bestType = bendTypeCont_.find(theKeys);  
239 >        return bestType;
240 >      } else {        
241          //if no exact match found, try wild card match
242          return bendTypeCont_.find(keys, wildCardAtomTypeName_);      
243        }
244      }
245    }
246  
242
247    TorsionType* ForceField::getTorsionType(const std::string &at1,
248                                            const std::string &at2,
249                                            const std::string &at3,
# Line 327 | Line 331 | namespace oopse {
331          jj++;
332        }
333        
334 <      std::sort(foundTorsions.begin(), foundTorsions.end());
335 <
336 <      int jkscore = foundTorsions[0].first;
337 <      int ilscore = foundTorsions[0].second;
338 <      std::vector<std::string> theKeys = foundTorsions[0].third;
339 <
340 <      std::cout << "best matching torsion = " << theKeys[0] << "\t" <<theKeys[1]  << "\t" << theKeys[2] << "\t" << theKeys[3] << "\t(scores = "<< jkscore << "\t" << ilscore << ")\n";
341 <
338 <      
339 <      TorsionType* bestType = torsionTypeCont_.find(theKeys);
340 <      if (bestType) {
341 <        return bestType;
334 >      if (foundTorsions.size() > 0) {
335 >        std::sort(foundTorsions.begin(), foundTorsions.end());
336 >        int jkscore = foundTorsions[0].first;
337 >        int ilscore = foundTorsions[0].second;
338 >        std::vector<std::string> theKeys = foundTorsions[0].third;
339 >        
340 >        TorsionType* bestType = torsionTypeCont_.find(theKeys);
341 >        return bestType;
342        } else {
343          //if no exact match found, try wild card match
344          return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
# Line 357 | Line 357 | namespace oopse {
357      keys.push_back(at4);    
358  
359      //try exact match first
360 <    InversionType* inversionType = inversionTypeCont_.find(keys);
360 >    InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(keys);
361      if (inversionType) {
362        return inversionType;
363      } else {
# Line 419 | Line 419 | namespace oopse {
419                myKeys.push_back((*k)->getName());
420                myKeys.push_back((*l)->getName());
421                
422 <              InversionType* inversionType = inversionTypeCont_.find(myKeys);
422 >              InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(myKeys);
423                if (inversionType) {
424                  foundInversions.push_back( make_tuple3(Iscore, JKLscore, myKeys) );
425                }
# Line 431 | Line 431 | namespace oopse {
431          }
432          jj++;
433        }
434 <      
435 <      std::sort(foundInversions.begin(), foundInversions.end());
436 <      
437 <      int iscore = foundInversions[0].first;
438 <      int jklscore = foundInversions[0].second;
439 <      std::vector<std::string> theKeys = foundInversions[0].third;
440 <      
441 <      std::cout << "best matching inversion = " << theKeys[0] << "\t" <<theKeys[1]  << "\t" << theKeys[2] << "\t" << theKeys[3] << "\t(scores = "<< iscore << "\t" << jklscore << ")\n";
442 <      
443 <      
444 <      InversionType* bestType = inversionTypeCont_.find(theKeys);
445 <      if (bestType) {
446 <        return bestType;
434 >        
435 >      if (foundInversions.size() > 0) {
436 >        std::sort(foundInversions.begin(), foundInversions.end());
437 >        int iscore = foundInversions[0].first;
438 >        int jklscore = foundInversions[0].second;
439 >        std::vector<std::string> theKeys = foundInversions[0].third;
440 >        
441 >        InversionType* bestType = inversionTypeCont_.permutedFindSkippingFirstElement(theKeys);
442 >        return bestType;
443        } else {
444          //if no exact match found, try wild card match
445          return inversionTypeCont_.find(keys, wildCardAtomTypeName_);
# Line 452 | Line 448 | namespace oopse {
448    }
449    
450    NonBondedInteractionType* ForceField::getNonBondedInteractionType(const std::string &at1, const std::string &at2) {
451 +    
452      std::vector<std::string> keys;
453      keys.push_back(at1);
454      keys.push_back(at2);    
# Line 461 | Line 458 | namespace oopse {
458      if (nbiType) {
459        return nbiType;
460      } else {
461 <      //if no exact match found, try wild card match
462 <      return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
463 <    }    
461 >      AtomType* atype1;
462 >      AtomType* atype2;
463 >      std::vector<std::string> at1key;
464 >      at1key.push_back(at1);
465 >      atype1 = atomTypeCont_.find(at1key);
466 >      
467 >      std::vector<std::string> at2key;
468 >      at2key.push_back(at2);
469 >      atype2 = atomTypeCont_.find(at2key);
470 >      
471 >      // query atom types for their chains of responsibility
472 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
473 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
474 >      
475 >      std::vector<AtomType*>::iterator i;
476 >      std::vector<AtomType*>::iterator j;
477 >      
478 >      int ii = 0;
479 >      int jj = 0;
480 >      int nbiTypeScore;
481 >      
482 >      std::vector<std::pair<int, std::vector<std::string> > > foundNBI;
483 >      
484 >      for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
485 >        jj = 0;
486 >        for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
487 >          
488 >          nbiTypeScore = ii + jj;
489 >          
490 >          std::vector<std::string> myKeys;
491 >          myKeys.push_back((*i)->getName());
492 >          myKeys.push_back((*j)->getName());
493 >          
494 >          NonBondedInteractionType* nbiType = nonBondedInteractionTypeCont_.find(myKeys);
495 >          if (nbiType) {
496 >            foundNBI.push_back(std::make_pair(nbiTypeScore, myKeys));
497 >          }
498 >          jj++;
499 >        }
500 >        ii++;
501 >      }
502 >      
503 >      
504 >      if (foundNBI.size() > 0) {
505 >        // sort the foundNBI by the score:
506 >        std::sort(foundNBI.begin(), foundNBI.end());
507 >        
508 >        int bestScore = foundNBI[0].first;
509 >        std::vector<std::string> theKeys = foundNBI[0].second;
510 >        
511 >        NonBondedInteractionType* bestType = nonBondedInteractionTypeCont_.find(theKeys);        
512 >        return bestType;
513 >      } else {
514 >        //if no exact match found, try wild card match
515 >        return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
516 >      }
517 >    }
518    }
519    
520    BondType* ForceField::getExactBondType(const std::string &at1,
# Line 519 | Line 570 | namespace oopse {
570    bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
571      std::vector<std::string> keys;
572      keys.push_back(at);
573 +    atypeIdentToName[atomType->getIdent()] = at;
574      return atomTypeCont_.add(keys, atomType);
575    }
576  
577 +  bool ForceField::replaceAtomType(const std::string &at, AtomType* atomType) {
578 +    std::vector<std::string> keys;
579 +    keys.push_back(at);
580 +    atypeIdentToName[atomType->getIdent()] = at;
581 +    return atomTypeCont_.replace(keys, atomType);
582 +  }
583 +
584    bool ForceField::addBondType(const std::string &at1, const std::string &at2,
585                                 BondType* bondType) {
586      std::vector<std::string> keys;
# Line 593 | Line 652 | namespace oopse {
652          } else {
653            sprintf( painCave.errMsg,
654                     "Can not cast GenericData to LJParam\n");
655 <          painCave.severity = OOPSE_ERROR;
655 >          painCave.severity = OPENMD_ERROR;
656            painCave.isFatal = 1;
657            simError();          
658          }            
659        } else {
660          sprintf( painCave.errMsg, "Can not find Parameters for LennardJones\n");
661 <        painCave.severity = OOPSE_ERROR;
661 >        painCave.severity = OPENMD_ERROR;
662          painCave.isFatal = 1;
663          simError();          
664        }
# Line 629 | Line 688 | namespace oopse {
688                   "\tHave you tried setting the FORCE_PARAM_PATH environment "
689                   "variable?\n",
690                   forceFieldFilename.c_str() );
691 <        painCave.severity = OOPSE_ERROR;
691 >        painCave.severity = OPENMD_ERROR;
692          painCave.isFatal = 1;
693          simError();
694        }
# Line 637 | Line 696 | namespace oopse {
696      return ffStream;
697    }
698  
699 <  void ForceField::setFortranForceOptions(){
641 <    ForceOptions theseFortranOptions;
642 <    forceFieldOptions_.makeFortranOptions(theseFortranOptions);
643 <    setfForceOptions(&theseFortranOptions);
644 <  }
645 < } //end namespace oopse
699 > } //end namespace OpenMD

Comparing:
trunk/src/UseTheForce/ForceField.cpp (property svn:keywords), Revision 1275 by cli2, Fri Jul 4 20:54:29 2008 UTC vs.
branches/development/src/UseTheForce/ForceField.cpp (property svn:keywords), Revision 1629 by gezelter, Wed Sep 14 21:15:17 2011 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines