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