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 1151 by chuckv, Fri Jul 6 18:15:03 2007 UTC vs.
branches/development/src/UseTheForce/ForceField.cpp (file contents), Revision 1710 by gezelter, Fri May 18 21:44:02 2012 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]  Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 + * [5]  Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
41   */
42  
43   /**
# Line 47 | Line 48
48   * @version 1.0
49   */
50    
51 + #include <algorithm>
52   #include "UseTheForce/ForceField.hpp"
53   #include "utils/simError.h"
54 < #include "UseTheForce/DarkSide/atype_interface.h"
55 < #include "UseTheForce/DarkSide/fForceOptions_interface.h"
54 < #include "UseTheForce/DarkSide/switcheroo_interface.h"
55 < namespace oopse {
54 > #include "utils/Tuple.hpp"
55 > #include "types/LennardJonesAdapter.hpp"
56  
57 + namespace OpenMD {
58 +
59    ForceField::ForceField() {
60 +
61      char* tempPath;
62      tempPath = getenv("FORCE_PARAM_PATH");
63 <
63 >    
64      if (tempPath == NULL) {
65        //convert a macro from compiler to a string in c++
66        STR_DEFINE(ffPath_, FRC_PATH );
# Line 66 | Line 69 | namespace oopse {
69      }
70    }
71  
72 <
73 <  ForceField::~ForceField() {
74 <    deleteAtypes();
75 <    deleteSwitch();
76 <  }
77 <
72 >  /**
73 >   * getAtomType by string
74 >   *
75 >   * finds the requested atom type in this force field using the string
76 >   * name of the atom type.
77 >   */
78    AtomType* ForceField::getAtomType(const std::string &at) {
79      std::vector<std::string> keys;
80      keys.push_back(at);
81      return atomTypeCont_.find(keys);
82    }
83  
84 <  BondType* ForceField::getBondType(const std::string &at1, const std::string &at2) {
84 >  /**
85 >   * getAtomType by ident
86 >   *
87 >   * finds the requested atom type in this force field using the
88 >   * integer ident instead of the string name of the atom type.
89 >   */
90 >  AtomType* ForceField::getAtomType(int ident) {  
91 >    std::string at = atypeIdentToName.find(ident)->second;
92 >    return getAtomType(at);
93 >  }
94 >
95 >  BondType* ForceField::getBondType(const std::string &at1,
96 >                                    const std::string &at2) {
97      std::vector<std::string> keys;
98      keys.push_back(at1);
99      keys.push_back(at2);    
# Line 88 | Line 103 | namespace oopse {
103      if (bondType) {
104        return bondType;
105      } else {
106 <      //if no exact match found, try wild card match
107 <      return bondTypeCont_.find(keys, wildCardAtomTypeName_);
108 <    }
106 >      AtomType* atype1;
107 >      AtomType* atype2;
108 >      std::vector<std::string> at1key;
109 >      at1key.push_back(at1);
110 >      atype1 = atomTypeCont_.find(at1key);
111 >  
112 >      std::vector<std::string> at2key;
113 >      at2key.push_back(at2);
114 >      atype2 = atomTypeCont_.find(at2key);
115  
116 <  }
116 >      // query atom types for their chains of responsibility
117 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
118 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
119  
120 <  BendType* ForceField::getBendType(const std::string &at1, const std::string &at2,
120 >      std::vector<AtomType*>::iterator i;
121 >      std::vector<AtomType*>::iterator j;
122 >
123 >      int ii = 0;
124 >      int jj = 0;
125 >      int bondTypeScore;
126 >
127 >      std::vector<std::pair<int, std::vector<std::string> > > foundBonds;
128 >
129 >      for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
130 >        jj = 0;
131 >        for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
132 >
133 >          bondTypeScore = ii + jj;
134 >
135 >          std::vector<std::string> myKeys;
136 >          myKeys.push_back((*i)->getName());
137 >          myKeys.push_back((*j)->getName());
138 >
139 >          BondType* bondType = bondTypeCont_.find(myKeys);
140 >          if (bondType) {
141 >            foundBonds.push_back(std::make_pair(bondTypeScore, myKeys));
142 >          }
143 >          jj++;
144 >        }
145 >        ii++;
146 >      }
147 >
148 >
149 >      if (foundBonds.size() > 0) {
150 >        // sort the foundBonds by the score:
151 >        std::sort(foundBonds.begin(), foundBonds.end());
152 >    
153 >        int bestScore = foundBonds[0].first;
154 >        std::vector<std::string> theKeys = foundBonds[0].second;
155 >        
156 >        BondType* bestType = bondTypeCont_.find(theKeys);
157 >        
158 >        return bestType;
159 >      } else {
160 >        //if no exact match found, try wild card match
161 >        return bondTypeCont_.find(keys, wildCardAtomTypeName_);      
162 >      }
163 >    }
164 >  }
165 >  
166 >  BendType* ForceField::getBendType(const std::string &at1,
167 >                                    const std::string &at2,
168                                      const std::string &at3) {
169      std::vector<std::string> keys;
170      keys.push_back(at1);
# Line 106 | Line 176 | namespace oopse {
176      if (bendType) {
177        return bendType;
178      } else {
179 <      //if no exact match found, try wild card match
180 <      return bendTypeCont_.find(keys, wildCardAtomTypeName_);
179 >
180 >      AtomType* atype1;
181 >      AtomType* atype2;
182 >      AtomType* atype3;
183 >      std::vector<std::string> at1key;
184 >      at1key.push_back(at1);
185 >      atype1 = atomTypeCont_.find(at1key);
186 >  
187 >      std::vector<std::string> at2key;
188 >      at2key.push_back(at2);
189 >      atype2 = atomTypeCont_.find(at2key);
190 >
191 >      std::vector<std::string> at3key;
192 >      at3key.push_back(at3);
193 >      atype3 = atomTypeCont_.find(at3key);
194 >
195 >      // query atom types for their chains of responsibility
196 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
197 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
198 >      std::vector<AtomType*> at3Chain = atype3->allYourBase();
199 >
200 >      std::vector<AtomType*>::iterator i;
201 >      std::vector<AtomType*>::iterator j;
202 >      std::vector<AtomType*>::iterator k;
203 >
204 >      int ii = 0;
205 >      int jj = 0;
206 >      int kk = 0;
207 >      int IKscore;
208 >
209 >      std::vector<tuple3<int, int, std::vector<std::string> > > foundBends;
210 >
211 >      for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
212 >        ii = 0;
213 >        for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
214 >          kk = 0;
215 >          for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
216 >          
217 >            IKscore = ii + kk;
218 >
219 >            std::vector<std::string> myKeys;
220 >            myKeys.push_back((*i)->getName());
221 >            myKeys.push_back((*j)->getName());
222 >            myKeys.push_back((*k)->getName());
223 >
224 >            BendType* bendType = bendTypeCont_.find(myKeys);
225 >            if (bendType) {
226 >              foundBends.push_back( make_tuple3(jj, IKscore, myKeys) );
227 >            }
228 >            kk++;
229 >          }
230 >          ii++;
231 >        }
232 >        jj++;
233 >      }
234 >      
235 >      if (foundBends.size() > 0) {
236 >        std::sort(foundBends.begin(), foundBends.end());
237 >        int jscore = foundBends[0].first;
238 >        int ikscore = foundBends[0].second;
239 >        std::vector<std::string> theKeys = foundBends[0].third;      
240 >        
241 >        BendType* bestType = bendTypeCont_.find(theKeys);  
242 >        return bestType;
243 >      } else {        
244 >        //if no exact match found, try wild card match
245 >        return bendTypeCont_.find(keys, wildCardAtomTypeName_);      
246 >      }
247      }
248    }
249  
250 <  TorsionType* ForceField::getTorsionType(const std::string &at1, const std::string &at2,
251 <                                          const std::string &at3, const std::string &at4) {
250 >  TorsionType* ForceField::getTorsionType(const std::string &at1,
251 >                                          const std::string &at2,
252 >                                          const std::string &at3,
253 >                                          const std::string &at4) {
254      std::vector<std::string> keys;
255      keys.push_back(at1);
256      keys.push_back(at2);    
257      keys.push_back(at3);    
258      keys.push_back(at4);    
259  
260 +
261 +    //try exact match first
262      TorsionType* torsionType = torsionTypeCont_.find(keys);
263      if (torsionType) {
264        return torsionType;
265      } else {
126      //if no exact match found, try wild card match
127      return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
128    }
129    
130    return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
266  
267 +      AtomType* atype1;
268 +      AtomType* atype2;
269 +      AtomType* atype3;
270 +      AtomType* atype4;
271 +      std::vector<std::string> at1key;
272 +      at1key.push_back(at1);
273 +      atype1 = atomTypeCont_.find(at1key);
274 +  
275 +      std::vector<std::string> at2key;
276 +      at2key.push_back(at2);
277 +      atype2 = atomTypeCont_.find(at2key);
278 +
279 +      std::vector<std::string> at3key;
280 +      at3key.push_back(at3);
281 +      atype3 = atomTypeCont_.find(at3key);
282 +
283 +      std::vector<std::string> at4key;
284 +      at4key.push_back(at4);
285 +      atype4 = atomTypeCont_.find(at4key);
286 +
287 +      // query atom types for their chains of responsibility
288 +      std::vector<AtomType*> at1Chain = atype1->allYourBase();
289 +      std::vector<AtomType*> at2Chain = atype2->allYourBase();
290 +      std::vector<AtomType*> at3Chain = atype3->allYourBase();
291 +      std::vector<AtomType*> at4Chain = atype4->allYourBase();
292 +
293 +      std::vector<AtomType*>::iterator i;
294 +      std::vector<AtomType*>::iterator j;
295 +      std::vector<AtomType*>::iterator k;
296 +      std::vector<AtomType*>::iterator l;
297 +
298 +      int ii = 0;
299 +      int jj = 0;
300 +      int kk = 0;
301 +      int ll = 0;
302 +      int ILscore;
303 +      int JKscore;
304 +
305 +      std::vector<tuple3<int, int, std::vector<std::string> > > foundTorsions;
306 +
307 +      for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
308 +        kk = 0;
309 +        for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
310 +          ii = 0;      
311 +          for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
312 +            ll = 0;
313 +            for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
314 +          
315 +              ILscore = ii + ll;
316 +              JKscore = jj + kk;
317 +
318 +              std::vector<std::string> myKeys;
319 +              myKeys.push_back((*i)->getName());
320 +              myKeys.push_back((*j)->getName());
321 +              myKeys.push_back((*k)->getName());
322 +              myKeys.push_back((*l)->getName());
323 +
324 +              TorsionType* torsionType = torsionTypeCont_.find(myKeys);
325 +              if (torsionType) {
326 +                foundTorsions.push_back( make_tuple3(JKscore, ILscore, myKeys) );
327 +              }
328 +              ll++;
329 +            }
330 +            ii++;
331 +          }
332 +          kk++;
333 +        }
334 +        jj++;
335 +      }
336 +      
337 +      if (foundTorsions.size() > 0) {
338 +        std::sort(foundTorsions.begin(), foundTorsions.end());
339 +        int jkscore = foundTorsions[0].first;
340 +        int ilscore = foundTorsions[0].second;
341 +        std::vector<std::string> theKeys = foundTorsions[0].third;
342 +        
343 +        TorsionType* bestType = torsionTypeCont_.find(theKeys);
344 +        return bestType;
345 +      } else {
346 +        //if no exact match found, try wild card match
347 +        return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
348 +      }
349 +    }
350    }
351  
352 <  NonBondedInteractionType* ForceField::getNonBondedInteractionType(const std::string &at1, const std::string &at2) {
352 >  InversionType* ForceField::getInversionType(const std::string &at1,
353 >                                              const std::string &at2,
354 >                                              const std::string &at3,
355 >                                              const std::string &at4) {
356      std::vector<std::string> keys;
357      keys.push_back(at1);
358      keys.push_back(at2);    
359 +    keys.push_back(at3);    
360 +    keys.push_back(at4);    
361  
362      //try exact match first
363 +    InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(keys);
364 +    if (inversionType) {
365 +      return inversionType;
366 +    } else {
367 +      
368 +      AtomType* atype1;
369 +      AtomType* atype2;
370 +      AtomType* atype3;
371 +      AtomType* atype4;
372 +      std::vector<std::string> at1key;
373 +      at1key.push_back(at1);
374 +      atype1 = atomTypeCont_.find(at1key);
375 +      
376 +      std::vector<std::string> at2key;
377 +      at2key.push_back(at2);
378 +      atype2 = atomTypeCont_.find(at2key);
379 +      
380 +      std::vector<std::string> at3key;
381 +      at3key.push_back(at3);
382 +      atype3 = atomTypeCont_.find(at3key);
383 +      
384 +      std::vector<std::string> at4key;
385 +      at4key.push_back(at4);
386 +      atype4 = atomTypeCont_.find(at4key);
387 +
388 +      // query atom types for their chains of responsibility
389 +      std::vector<AtomType*> at1Chain = atype1->allYourBase();
390 +      std::vector<AtomType*> at2Chain = atype2->allYourBase();
391 +      std::vector<AtomType*> at3Chain = atype3->allYourBase();
392 +      std::vector<AtomType*> at4Chain = atype4->allYourBase();
393 +
394 +      std::vector<AtomType*>::iterator i;
395 +      std::vector<AtomType*>::iterator j;
396 +      std::vector<AtomType*>::iterator k;
397 +      std::vector<AtomType*>::iterator l;
398 +
399 +      int ii = 0;
400 +      int jj = 0;
401 +      int kk = 0;
402 +      int ll = 0;
403 +      int Iscore;
404 +      int JKLscore;
405 +      
406 +      std::vector<tuple3<int, int, std::vector<std::string> > > foundInversions;
407 +      
408 +      for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
409 +        kk = 0;
410 +        for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
411 +          ii = 0;      
412 +          for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
413 +            ll = 0;
414 +            for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
415 +              
416 +              Iscore = ii;
417 +              JKLscore = jj + kk + ll;
418 +              
419 +              std::vector<std::string> myKeys;
420 +              myKeys.push_back((*i)->getName());
421 +              myKeys.push_back((*j)->getName());
422 +              myKeys.push_back((*k)->getName());
423 +              myKeys.push_back((*l)->getName());
424 +              
425 +              InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(myKeys);
426 +              if (inversionType) {
427 +                foundInversions.push_back( make_tuple3(Iscore, JKLscore, myKeys) );
428 +              }
429 +              ll++;
430 +            }
431 +            ii++;
432 +          }
433 +          kk++;
434 +        }
435 +        jj++;
436 +      }
437 +        
438 +      if (foundInversions.size() > 0) {
439 +        std::sort(foundInversions.begin(), foundInversions.end());
440 +        int iscore = foundInversions[0].first;
441 +        int jklscore = foundInversions[0].second;
442 +        std::vector<std::string> theKeys = foundInversions[0].third;
443 +        
444 +        InversionType* bestType = inversionTypeCont_.permutedFindSkippingFirstElement(theKeys);
445 +        return bestType;
446 +      } else {
447 +        //if no exact match found, try wild card match
448 +        return inversionTypeCont_.find(keys, wildCardAtomTypeName_);
449 +      }
450 +    }
451 +  }
452 +  
453 +  NonBondedInteractionType* ForceField::getNonBondedInteractionType(const std::string &at1, const std::string &at2) {
454 +    
455 +    std::vector<std::string> keys;
456 +    keys.push_back(at1);
457 +    keys.push_back(at2);    
458 +    
459 +    //try exact match first
460      NonBondedInteractionType* nbiType = nonBondedInteractionTypeCont_.find(keys);
461      if (nbiType) {
462        return nbiType;
463      } else {
464 <      //if no exact match found, try wild card match
465 <      return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
464 >      AtomType* atype1;
465 >      AtomType* atype2;
466 >      std::vector<std::string> at1key;
467 >      at1key.push_back(at1);
468 >      atype1 = atomTypeCont_.find(at1key);
469 >      
470 >      std::vector<std::string> at2key;
471 >      at2key.push_back(at2);
472 >      atype2 = atomTypeCont_.find(at2key);
473 >      
474 >      // query atom types for their chains of responsibility
475 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
476 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
477 >      
478 >      std::vector<AtomType*>::iterator i;
479 >      std::vector<AtomType*>::iterator j;
480 >      
481 >      int ii = 0;
482 >      int jj = 0;
483 >      int nbiTypeScore;
484 >      
485 >      std::vector<std::pair<int, std::vector<std::string> > > foundNBI;
486 >      
487 >      for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
488 >        jj = 0;
489 >        for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
490 >          
491 >          nbiTypeScore = ii + jj;
492 >          
493 >          std::vector<std::string> myKeys;
494 >          myKeys.push_back((*i)->getName());
495 >          myKeys.push_back((*j)->getName());
496 >          
497 >          NonBondedInteractionType* nbiType = nonBondedInteractionTypeCont_.find(myKeys);
498 >          if (nbiType) {
499 >            foundNBI.push_back(std::make_pair(nbiTypeScore, myKeys));
500 >          }
501 >          jj++;
502 >        }
503 >        ii++;
504 >      }
505 >      
506 >      
507 >      if (foundNBI.size() > 0) {
508 >        // sort the foundNBI by the score:
509 >        std::sort(foundNBI.begin(), foundNBI.end());
510 >        
511 >        int bestScore = foundNBI[0].first;
512 >        std::vector<std::string> theKeys = foundNBI[0].second;
513 >        
514 >        NonBondedInteractionType* bestType = nonBondedInteractionTypeCont_.find(theKeys);        
515 >        return bestType;
516 >      } else {
517 >        //if no exact match found, try wild card match
518 >        return nonBondedInteractionTypeCont_.find(keys, wildCardAtomTypeName_);
519 >      }
520      }
147
521    }
522 <
523 <  BondType* ForceField::getExactBondType(const std::string &at1, const std::string &at2){
522 >  
523 >  BondType* ForceField::getExactBondType(const std::string &at1,
524 >                                         const std::string &at2){
525      std::vector<std::string> keys;
526      keys.push_back(at1);
527      keys.push_back(at2);    
528      return bondTypeCont_.find(keys);
529    }
530 <
531 <  BendType* ForceField::getExactBendType(const std::string &at1, const std::string &at2,
530 >  
531 >  BendType* ForceField::getExactBendType(const std::string &at1,
532 >                                         const std::string &at2,
533                                           const std::string &at3){
534      std::vector<std::string> keys;
535      keys.push_back(at1);
# Line 162 | Line 537 | namespace oopse {
537      keys.push_back(at3);    
538      return bendTypeCont_.find(keys);
539    }
540 <
541 <  TorsionType* ForceField::getExactTorsionType(const std::string &at1, const std::string &at2,
542 <                                               const std::string &at3, const std::string &at4){
540 >  
541 >  TorsionType* ForceField::getExactTorsionType(const std::string &at1,
542 >                                               const std::string &at2,
543 >                                               const std::string &at3,
544 >                                               const std::string &at4){
545      std::vector<std::string> keys;
546      keys.push_back(at1);
547      keys.push_back(at2);    
# Line 172 | Line 549 | namespace oopse {
549      keys.push_back(at4);  
550      return torsionTypeCont_.find(keys);
551    }
552 <
552 >  
553 >  InversionType* ForceField::getExactInversionType(const std::string &at1,
554 >                                                   const std::string &at2,
555 >                                                   const std::string &at3,
556 >                                                   const std::string &at4){
557 >    std::vector<std::string> keys;
558 >    keys.push_back(at1);
559 >    keys.push_back(at2);    
560 >    keys.push_back(at3);    
561 >    keys.push_back(at4);  
562 >    return inversionTypeCont_.find(keys);
563 >  }
564 >  
565    NonBondedInteractionType* ForceField::getExactNonBondedInteractionType(const std::string &at1, const std::string &at2){
566      std::vector<std::string> keys;
567      keys.push_back(at1);
568      keys.push_back(at2);    
569      return nonBondedInteractionTypeCont_.find(keys);
570    }
571 +  
572  
183
573    bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
574      std::vector<std::string> keys;
575      keys.push_back(at);
576 +    atypeIdentToName[atomType->getIdent()] = at;
577      return atomTypeCont_.add(keys, atomType);
578    }
579  
580 <  bool ForceField::addBondType(const std::string &at1, const std::string &at2, BondType* bondType) {
580 >  bool ForceField::replaceAtomType(const std::string &at, AtomType* atomType) {
581      std::vector<std::string> keys;
582 +    keys.push_back(at);
583 +    atypeIdentToName[atomType->getIdent()] = at;
584 +    return atomTypeCont_.replace(keys, atomType);
585 +  }
586 +
587 +  bool ForceField::addBondType(const std::string &at1, const std::string &at2,
588 +                               BondType* bondType) {
589 +    std::vector<std::string> keys;
590      keys.push_back(at1);
591      keys.push_back(at2);    
592 <    return bondTypeCont_.add(keys, bondType);
195 <
592 >    return bondTypeCont_.add(keys, bondType);    
593    }
594 <
594 >  
595    bool ForceField::addBendType(const std::string &at1, const std::string &at2,
596                                 const std::string &at3, BendType* bendType) {
597      std::vector<std::string> keys;
# Line 203 | Line 600 | namespace oopse {
600      keys.push_back(at3);    
601      return bendTypeCont_.add(keys, bendType);
602    }
603 <
604 <  bool ForceField::addTorsionType(const std::string &at1, const std::string &at2,
605 <                                  const std::string &at3, const std::string &at4, TorsionType* torsionType) {
603 >  
604 >  bool ForceField::addTorsionType(const std::string &at1,
605 >                                  const std::string &at2,
606 >                                  const std::string &at3,
607 >                                  const std::string &at4,
608 >                                  TorsionType* torsionType) {
609      std::vector<std::string> keys;
610      keys.push_back(at1);
611      keys.push_back(at2);    
# Line 214 | Line 614 | namespace oopse {
614      return torsionTypeCont_.add(keys, torsionType);
615    }
616  
617 <  bool ForceField::addNonBondedInteractionType(const std::string &at1, const std::string &at2, NonBondedInteractionType* nbiType) {
617 >  bool ForceField::addInversionType(const std::string &at1,
618 >                                    const std::string &at2,
619 >                                    const std::string &at3,
620 >                                    const std::string &at4,
621 >                                    InversionType* inversionType) {
622      std::vector<std::string> keys;
623      keys.push_back(at1);
624      keys.push_back(at2);    
625 +    keys.push_back(at3);    
626 +    keys.push_back(at4);    
627 +    return inversionTypeCont_.add(keys, inversionType);
628 +  }
629 +  
630 +  bool ForceField::addNonBondedInteractionType(const std::string &at1,
631 +                                               const std::string &at2,
632 +                                               NonBondedInteractionType* nbiType) {
633 +    std::vector<std::string> keys;
634 +    keys.push_back(at1);
635 +    keys.push_back(at2);    
636      return nonBondedInteractionTypeCont_.add(keys, nbiType);
222
637    }
638 <
638 >  
639    RealType ForceField::getRcutFromAtomType(AtomType* at) {
226    /**@todo */
227    GenericData* data;
640      RealType rcut = 0.0;
641 <
642 <    if (at->isLennardJones()) {
643 <      data = at->getPropertyByName("LennardJones");
644 <      if (data != NULL) {
233 <        LJParamGenericData* ljData = dynamic_cast<LJParamGenericData*>(data);
234 <
235 <        if (ljData != NULL) {
236 <          LJParam ljParam = ljData->getData();
237 <
238 <          //by default use 2.5*sigma as cutoff radius
239 <          rcut = 2.5 * ljParam.sigma;
240 <                
241 <        } else {
242 <          sprintf( painCave.errMsg,
243 <                   "Can not cast GenericData to LJParam\n");
244 <          painCave.severity = OOPSE_ERROR;
245 <          painCave.isFatal = 1;
246 <          simError();          
247 <        }            
248 <      } else {
249 <        sprintf( painCave.errMsg, "Can not find Parameters for LennardJones\n");
250 <        painCave.severity = OOPSE_ERROR;
251 <        painCave.isFatal = 1;
252 <        simError();          
253 <      }
641 >    
642 >    LennardJonesAdapter lja = LennardJonesAdapter(at);
643 >    if (lja.isLennardJones()) {
644 >      rcut = 2.5 * lja.getSigma();
645      }
255
646      return rcut;    
647    }
648 +  
649  
259
650    ifstrstream* ForceField::openForceFieldFile(const std::string& filename) {
651      std::string forceFieldFilename(filename);
652      ifstrstream* ffStream = new ifstrstream();
# Line 278 | Line 668 | namespace oopse {
668                   "\tHave you tried setting the FORCE_PARAM_PATH environment "
669                   "variable?\n",
670                   forceFieldFilename.c_str() );
671 <        painCave.severity = OOPSE_ERROR;
671 >        painCave.severity = OPENMD_ERROR;
672          painCave.isFatal = 1;
673          simError();
674        }
675      }  
286
676      return ffStream;
288
677    }
678  
679 <  void ForceField::setFortranForceOptions(){
292 <    ForceOptions theseFortranOptions;
293 <    forceFieldOptions_.makeFortranOptions(theseFortranOptions);
294 <    setfForceOptions(&theseFortranOptions);
295 <  }
296 < } //end namespace oopse
679 > } //end namespace OpenMD

Comparing:
trunk/src/UseTheForce/ForceField.cpp (property svn:keywords), Revision 1151 by chuckv, Fri Jul 6 18:15:03 2007 UTC vs.
branches/development/src/UseTheForce/ForceField.cpp (property svn:keywords), Revision 1710 by gezelter, Fri May 18 21:44:02 2012 UTC

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

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines