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 1195 by cpuglis, Thu Dec 6 20:04:02 2007 UTC vs.
Revision 1303 by cli2, Mon Oct 13 21:35:22 2008 UTC

# Line 47 | Line 47
47   * @version 1.0
48   */
49    
50 + #include <algorithm>
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"
# Line 89 | Line 91 | namespace oopse {
91      if (bondType) {
92        return bondType;
93      } else {
94 <      //if no exact match found, try wild card match
95 <      return bondTypeCont_.find(keys, wildCardAtomTypeName_);
94 >      AtomType* atype1;
95 >      AtomType* atype2;
96 >      std::vector<std::string> at1key;
97 >      at1key.push_back(at1);
98 >      atype1 = atomTypeCont_.find(at1key);
99 >  
100 >      std::vector<std::string> at2key;
101 >      at2key.push_back(at2);
102 >      atype2 = atomTypeCont_.find(at2key);
103 >
104 >      // query atom types for their chains of responsibility
105 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
106 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
107 >
108 >      std::vector<AtomType*>::iterator i;
109 >      std::vector<AtomType*>::iterator j;
110 >
111 >      int ii = 0;
112 >      int jj = 0;
113 >      int bondTypeScore;
114 >
115 >      std::vector<std::pair<int, std::vector<std::string> > > foundBonds;
116 >
117 >      for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
118 >        jj = 0;
119 >        for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
120 >
121 >          bondTypeScore = ii + jj;
122 >
123 >          std::vector<std::string> myKeys;
124 >          myKeys.push_back((*i)->getName());
125 >          myKeys.push_back((*j)->getName());
126 >
127 >          BondType* bondType = bondTypeCont_.find(myKeys);
128 >          if (bondType) {
129 >            foundBonds.push_back(std::make_pair(bondTypeScore, myKeys));
130 >          }
131 >          jj++;
132 >        }
133 >        ii++;
134 >      }
135 >
136 >
137 >      if (foundBonds.size() > 0) {
138 >        // sort the foundBonds by the score:
139 >        std::sort(foundBonds.begin(), foundBonds.end());
140 >    
141 >        int bestScore = foundBonds[0].first;
142 >        std::vector<std::string> theKeys = foundBonds[0].second;
143 >        
144 >        BondType* bestType = bondTypeCont_.find(theKeys);
145 >        
146 >        return bestType;
147 >      } else {
148 >        //if no exact match found, try wild card match
149 >        return bondTypeCont_.find(keys, wildCardAtomTypeName_);      
150 >      }
151      }
152    }
153 <
153 >  
154    BendType* ForceField::getBendType(const std::string &at1,
155                                      const std::string &at2,
156                                      const std::string &at3) {
# Line 107 | Line 164 | namespace oopse {
164      if (bendType) {
165        return bendType;
166      } else {
167 <      //if no exact match found, try wild card match
168 <      return bendTypeCont_.find(keys, wildCardAtomTypeName_);
167 >
168 >      AtomType* atype1;
169 >      AtomType* atype2;
170 >      AtomType* atype3;
171 >      std::vector<std::string> at1key;
172 >      at1key.push_back(at1);
173 >      atype1 = atomTypeCont_.find(at1key);
174 >  
175 >      std::vector<std::string> at2key;
176 >      at2key.push_back(at2);
177 >      atype2 = atomTypeCont_.find(at2key);
178 >
179 >      std::vector<std::string> at3key;
180 >      at3key.push_back(at3);
181 >      atype3 = atomTypeCont_.find(at3key);
182 >
183 >      // query atom types for their chains of responsibility
184 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
185 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
186 >      std::vector<AtomType*> at3Chain = atype3->allYourBase();
187 >
188 >      std::vector<AtomType*>::iterator i;
189 >      std::vector<AtomType*>::iterator j;
190 >      std::vector<AtomType*>::iterator k;
191 >
192 >      int ii = 0;
193 >      int jj = 0;
194 >      int kk = 0;
195 >      int IKscore;
196 >
197 >      std::vector<tuple3<int, int, std::vector<std::string> > > foundBends;
198 >
199 >      for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
200 >        ii = 0;
201 >        for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
202 >          kk = 0;
203 >          for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
204 >          
205 >            IKscore = ii + kk;
206 >
207 >            std::vector<std::string> myKeys;
208 >            myKeys.push_back((*i)->getName());
209 >            myKeys.push_back((*j)->getName());
210 >            myKeys.push_back((*k)->getName());
211 >
212 >            BendType* bendType = bendTypeCont_.find(myKeys);
213 >            if (bendType) {
214 >              foundBends.push_back( make_tuple3(jj, IKscore, myKeys) );
215 >            }
216 >            kk++;
217 >          }
218 >          ii++;
219 >        }
220 >        jj++;
221 >      }
222 >      
223 >      if (foundBends.size() > 0) {
224 >        std::sort(foundBends.begin(), foundBends.end());
225 >        int jscore = foundBends[0].first;
226 >        int ikscore = foundBends[0].second;
227 >        std::vector<std::string> theKeys = foundBends[0].third;      
228 >        
229 >        BendType* bestType = bendTypeCont_.find(theKeys);  
230 >        return bestType;
231 >      } else {        
232 >        //if no exact match found, try wild card match
233 >        return bendTypeCont_.find(keys, wildCardAtomTypeName_);      
234 >      }
235      }
236    }
237  
# Line 122 | Line 245 | namespace oopse {
245      keys.push_back(at3);    
246      keys.push_back(at4);    
247  
248 +
249 +    //try exact match first
250      TorsionType* torsionType = torsionTypeCont_.find(keys);
251      if (torsionType) {
252        return torsionType;
253      } else {
254 <      //if no exact match found, try wild card match
255 <      return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
256 <    }
257 <    
258 <    return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
254 >
255 >      AtomType* atype1;
256 >      AtomType* atype2;
257 >      AtomType* atype3;
258 >      AtomType* atype4;
259 >      std::vector<std::string> at1key;
260 >      at1key.push_back(at1);
261 >      atype1 = atomTypeCont_.find(at1key);
262 >  
263 >      std::vector<std::string> at2key;
264 >      at2key.push_back(at2);
265 >      atype2 = atomTypeCont_.find(at2key);
266 >
267 >      std::vector<std::string> at3key;
268 >      at3key.push_back(at3);
269 >      atype3 = atomTypeCont_.find(at3key);
270 >
271 >      std::vector<std::string> at4key;
272 >      at4key.push_back(at4);
273 >      atype4 = atomTypeCont_.find(at4key);
274 >
275 >      // query atom types for their chains of responsibility
276 >      std::vector<AtomType*> at1Chain = atype1->allYourBase();
277 >      std::vector<AtomType*> at2Chain = atype2->allYourBase();
278 >      std::vector<AtomType*> at3Chain = atype3->allYourBase();
279 >      std::vector<AtomType*> at4Chain = atype4->allYourBase();
280 >
281 >      std::vector<AtomType*>::iterator i;
282 >      std::vector<AtomType*>::iterator j;
283 >      std::vector<AtomType*>::iterator k;
284 >      std::vector<AtomType*>::iterator l;
285 >
286 >      int ii = 0;
287 >      int jj = 0;
288 >      int kk = 0;
289 >      int ll = 0;
290 >      int ILscore;
291 >      int JKscore;
292 >
293 >      std::vector<tuple3<int, int, std::vector<std::string> > > foundTorsions;
294 >
295 >      for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
296 >        kk = 0;
297 >        for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
298 >          ii = 0;      
299 >          for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
300 >            ll = 0;
301 >            for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
302 >          
303 >              ILscore = ii + ll;
304 >              JKscore = jj + kk;
305 >
306 >              std::vector<std::string> myKeys;
307 >              myKeys.push_back((*i)->getName());
308 >              myKeys.push_back((*j)->getName());
309 >              myKeys.push_back((*k)->getName());
310 >              myKeys.push_back((*l)->getName());
311 >
312 >              TorsionType* torsionType = torsionTypeCont_.find(myKeys);
313 >              if (torsionType) {
314 >                foundTorsions.push_back( make_tuple3(JKscore, ILscore, myKeys) );
315 >              }
316 >              ll++;
317 >            }
318 >            ii++;
319 >          }
320 >          kk++;
321 >        }
322 >        jj++;
323 >      }
324 >      
325 >      if (foundTorsions.size() > 0) {
326 >        std::sort(foundTorsions.begin(), foundTorsions.end());
327 >        int jkscore = foundTorsions[0].first;
328 >        int ilscore = foundTorsions[0].second;
329 >        std::vector<std::string> theKeys = foundTorsions[0].third;
330 >        
331 >        TorsionType* bestType = torsionTypeCont_.find(theKeys);
332 >        return bestType;
333 >      } else {
334 >        //if no exact match found, try wild card match
335 >        return torsionTypeCont_.find(keys, wildCardAtomTypeName_);
336 >      }
337 >    }
338    }
339  
340 +  InversionType* ForceField::getInversionType(const std::string &at1,
341 +                                              const std::string &at2,
342 +                                              const std::string &at3,
343 +                                              const std::string &at4) {
344 +    std::vector<std::string> keys;
345 +    keys.push_back(at1);
346 +    keys.push_back(at2);    
347 +    keys.push_back(at3);    
348 +    keys.push_back(at4);    
349 +
350 +    //try exact match first
351 +    InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(keys);
352 +    if (inversionType) {
353 +      return inversionType;
354 +    } else {
355 +      
356 +      AtomType* atype1;
357 +      AtomType* atype2;
358 +      AtomType* atype3;
359 +      AtomType* atype4;
360 +      std::vector<std::string> at1key;
361 +      at1key.push_back(at1);
362 +      atype1 = atomTypeCont_.find(at1key);
363 +      
364 +      std::vector<std::string> at2key;
365 +      at2key.push_back(at2);
366 +      atype2 = atomTypeCont_.find(at2key);
367 +      
368 +      std::vector<std::string> at3key;
369 +      at3key.push_back(at3);
370 +      atype3 = atomTypeCont_.find(at3key);
371 +      
372 +      std::vector<std::string> at4key;
373 +      at4key.push_back(at4);
374 +      atype4 = atomTypeCont_.find(at4key);
375 +
376 +      // query atom types for their chains of responsibility
377 +      std::vector<AtomType*> at1Chain = atype1->allYourBase();
378 +      std::vector<AtomType*> at2Chain = atype2->allYourBase();
379 +      std::vector<AtomType*> at3Chain = atype3->allYourBase();
380 +      std::vector<AtomType*> at4Chain = atype4->allYourBase();
381 +
382 +      std::vector<AtomType*>::iterator i;
383 +      std::vector<AtomType*>::iterator j;
384 +      std::vector<AtomType*>::iterator k;
385 +      std::vector<AtomType*>::iterator l;
386 +
387 +      int ii = 0;
388 +      int jj = 0;
389 +      int kk = 0;
390 +      int ll = 0;
391 +      int Iscore;
392 +      int JKLscore;
393 +      
394 +      std::vector<tuple3<int, int, std::vector<std::string> > > foundInversions;
395 +      
396 +      for (j = at2Chain.begin(); j != at2Chain.end(); j++) {
397 +        kk = 0;
398 +        for (k = at3Chain.begin(); k != at3Chain.end(); k++) {
399 +          ii = 0;      
400 +          for (i = at1Chain.begin(); i != at1Chain.end(); i++) {
401 +            ll = 0;
402 +            for (l = at4Chain.begin(); l != at4Chain.end(); l++) {
403 +              
404 +              Iscore = ii;
405 +              JKLscore = jj + kk + ll;
406 +              
407 +              std::vector<std::string> myKeys;
408 +              myKeys.push_back((*i)->getName());
409 +              myKeys.push_back((*j)->getName());
410 +              myKeys.push_back((*k)->getName());
411 +              myKeys.push_back((*l)->getName());
412 +              
413 +              InversionType* inversionType = inversionTypeCont_.permutedFindSkippingFirstElement(myKeys);
414 +              if (inversionType) {
415 +                foundInversions.push_back( make_tuple3(Iscore, JKLscore, myKeys) );
416 +              }
417 +              ll++;
418 +            }
419 +            ii++;
420 +          }
421 +          kk++;
422 +        }
423 +        jj++;
424 +      }
425 +        
426 +      if (foundInversions.size() > 0) {
427 +        std::sort(foundInversions.begin(), foundInversions.end());
428 +        int iscore = foundInversions[0].first;
429 +        int jklscore = foundInversions[0].second;
430 +        std::vector<std::string> theKeys = foundInversions[0].third;
431 +        
432 +        InversionType* bestType = inversionTypeCont_.permutedFindSkippingFirstElement(theKeys);
433 +        return bestType;
434 +      } else {
435 +        //if no exact match found, try wild card match
436 +        return inversionTypeCont_.find(keys, wildCardAtomTypeName_);
437 +      }
438 +    }
439 +  }
440 +  
441    NonBondedInteractionType* ForceField::getNonBondedInteractionType(const std::string &at1, const std::string &at2) {
442      std::vector<std::string> keys;
443      keys.push_back(at1);
# Line 177 | Line 482 | namespace oopse {
482      keys.push_back(at4);  
483      return torsionTypeCont_.find(keys);
484    }
485 <
485 >  
486 >  InversionType* ForceField::getExactInversionType(const std::string &at1,
487 >                                                   const std::string &at2,
488 >                                                   const std::string &at3,
489 >                                                   const std::string &at4){
490 >    std::vector<std::string> keys;
491 >    keys.push_back(at1);
492 >    keys.push_back(at2);    
493 >    keys.push_back(at3);    
494 >    keys.push_back(at4);  
495 >    return inversionTypeCont_.find(keys);
496 >  }
497 >  
498    NonBondedInteractionType* ForceField::getExactNonBondedInteractionType(const std::string &at1, const std::string &at2){
499      std::vector<std::string> keys;
500      keys.push_back(at1);
501      keys.push_back(at2);    
502      return nonBondedInteractionTypeCont_.find(keys);
503    }
504 +  
505  
188 –
506    bool ForceField::addAtomType(const std::string &at, AtomType* atomType) {
507      std::vector<std::string> keys;
508      keys.push_back(at);
509      return atomTypeCont_.add(keys, atomType);
510    }
511  
512 +  bool ForceField::replaceAtomType(const std::string &at, AtomType* atomType) {
513 +    std::vector<std::string> keys;
514 +    keys.push_back(at);
515 +    return atomTypeCont_.replace(keys, atomType);
516 +  }
517 +
518    bool ForceField::addBondType(const std::string &at1, const std::string &at2,
519                                 BondType* bondType) {
520      std::vector<std::string> keys;
# Line 222 | Line 545 | namespace oopse {
545      return torsionTypeCont_.add(keys, torsionType);
546    }
547  
548 +  bool ForceField::addInversionType(const std::string &at1,
549 +                                    const std::string &at2,
550 +                                    const std::string &at3,
551 +                                    const std::string &at4,
552 +                                    InversionType* inversionType) {
553 +    std::vector<std::string> keys;
554 +    keys.push_back(at1);
555 +    keys.push_back(at2);    
556 +    keys.push_back(at3);    
557 +    keys.push_back(at4);    
558 +    return inversionTypeCont_.add(keys, inversionType);
559 +  }
560 +  
561    bool ForceField::addNonBondedInteractionType(const std::string &at1,
562                                                 const std::string &at2,
563                                                 NonBondedInteractionType* nbiType) {

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines