43 |
|
|
44 |
|
namespace OpenMD { |
45 |
|
|
46 |
< |
InteractionManager* InteractionManager::_instance = NULL; |
47 |
< |
SimInfo* InteractionManager::info_ = NULL; |
48 |
< |
bool InteractionManager::initialized_ = false; |
46 |
> |
InteractionManager::InteractionManager() { |
47 |
|
|
48 |
< |
RealType InteractionManager::rCut_ = 0.0; |
49 |
< |
RealType InteractionManager::rSwitch_ = 0.0; |
50 |
< |
RealType InteractionManager::skinThickness_ = 0.0; |
51 |
< |
RealType InteractionManager::listRadius_ = 0.0; |
52 |
< |
CutoffMethod InteractionManager::cutoffMethod_ = SHIFTED_FORCE; |
53 |
< |
SwitchingFunctionType InteractionManager::sft_ = cubic; |
54 |
< |
RealType InteractionManager::vdwScale_[4] = {1.0, 0.0, 0.0, 0.0}; |
55 |
< |
RealType InteractionManager::electrostaticScale_[4] = {1.0, 0.0, 0.0, 0.0}; |
56 |
< |
|
57 |
< |
map<int, AtomType*> InteractionManager::typeMap_; |
60 |
< |
map<pair<AtomType*, AtomType*>, set<NonBondedInteraction*> > InteractionManager::interactions_; |
61 |
< |
|
62 |
< |
LJ* InteractionManager::lj_ = new LJ(); |
63 |
< |
GB* InteractionManager::gb_ = new GB(); |
64 |
< |
Sticky* InteractionManager::sticky_ = new Sticky(); |
65 |
< |
Morse* InteractionManager::morse_ = new Morse(); |
66 |
< |
EAM* InteractionManager::eam_ = new EAM(); |
67 |
< |
SC* InteractionManager::sc_ = new SC(); |
68 |
< |
Electrostatic* InteractionManager::electrostatic_ = new Electrostatic(); |
69 |
< |
MAW* InteractionManager::maw_ = new MAW(); |
70 |
< |
SwitchingFunction* InteractionManager::switcher_ = new SwitchingFunction(); |
71 |
< |
|
72 |
< |
InteractionManager* InteractionManager::Instance() { |
73 |
< |
if (!_instance) { |
74 |
< |
_instance = new InteractionManager(); |
75 |
< |
} |
76 |
< |
return _instance; |
48 |
> |
initialized_ = false; |
49 |
> |
|
50 |
> |
lj_ = new LJ(); |
51 |
> |
gb_ = new GB(); |
52 |
> |
sticky_ = new Sticky(); |
53 |
> |
morse_ = new Morse(); |
54 |
> |
eam_ = new EAM(); |
55 |
> |
sc_ = new SC(); |
56 |
> |
electrostatic_ = new Electrostatic(); |
57 |
> |
maw_ = new MAW(); |
58 |
|
} |
59 |
|
|
60 |
|
void InteractionManager::initialize() { |
70 |
|
electrostatic_->setForceField(forceField_); |
71 |
|
maw_->setForceField(forceField_); |
72 |
|
|
92 |
– |
ForceFieldOptions& fopts = forceField_->getForceFieldOptions(); |
93 |
– |
|
94 |
– |
// Force fields can set options on how to scale van der Waals and electrostatic |
95 |
– |
// interactions for atoms connected via bonds, bends and torsions |
96 |
– |
// in this case the topological distance between atoms is: |
97 |
– |
// 0 = topologically unconnected |
98 |
– |
// 1 = bonded together |
99 |
– |
// 2 = connected via a bend |
100 |
– |
// 3 = connected via a torsion |
101 |
– |
|
102 |
– |
vdwScale_[0] = 1.0; |
103 |
– |
vdwScale_[1] = fopts.getvdw12scale(); |
104 |
– |
vdwScale_[2] = fopts.getvdw13scale(); |
105 |
– |
vdwScale_[3] = fopts.getvdw14scale(); |
106 |
– |
|
107 |
– |
electrostaticScale_[0] = 1.0; |
108 |
– |
electrostaticScale_[1] = fopts.getelectrostatic12scale(); |
109 |
– |
electrostaticScale_[2] = fopts.getelectrostatic13scale(); |
110 |
– |
electrostaticScale_[3] = fopts.getelectrostatic14scale(); |
111 |
– |
|
73 |
|
ForceField::AtomTypeContainer* atomTypes = forceField_->getAtomTypes(); |
74 |
|
ForceField::AtomTypeContainer::MapTypeIterator i1, i2; |
75 |
|
AtomType* atype1; |
263 |
|
} |
264 |
|
} |
265 |
|
} |
305 |
– |
|
306 |
– |
setupCutoffs(); |
307 |
– |
setupSwitching(); |
308 |
– |
setupNeighborlists(); |
266 |
|
|
310 |
– |
//int ljsp = cutoffMethod_ == SHIFTED_POTENTIAL ? 1 : 0; |
311 |
– |
//int ljsf = cutoffMethod_ == SHIFTED_FORCE ? 1 : 0; |
312 |
– |
//notifyFortranCutoffs(&rCut_, &rSwitch_, &ljsp, &ljsf); |
313 |
– |
//notifyFortranSkinThickness(&skinThickness_); |
314 |
– |
|
267 |
|
initialized_ = true; |
316 |
– |
} |
317 |
– |
|
318 |
– |
/** |
319 |
– |
* setupCutoffs |
320 |
– |
* |
321 |
– |
* Sets the values of cutoffRadius and cutoffMethod |
322 |
– |
* |
323 |
– |
* cutoffRadius : realType |
324 |
– |
* If the cutoffRadius was explicitly set, use that value. |
325 |
– |
* If the cutoffRadius was not explicitly set: |
326 |
– |
* Are there electrostatic atoms? Use 12.0 Angstroms. |
327 |
– |
* No electrostatic atoms? Poll the atom types present in the |
328 |
– |
* simulation for suggested cutoff values (e.g. 2.5 * sigma). |
329 |
– |
* Use the maximum suggested value that was found. |
330 |
– |
* |
331 |
– |
* cutoffMethod : (one of HARD, SWITCHED, SHIFTED_FORCE, SHIFTED_POTENTIAL) |
332 |
– |
* If cutoffMethod was explicitly set, use that choice. |
333 |
– |
* If cutoffMethod was not explicitly set, use SHIFTED_FORCE |
334 |
– |
*/ |
335 |
– |
void InteractionManager::setupCutoffs() { |
336 |
– |
|
337 |
– |
Globals* simParams_ = info_->getSimParams(); |
338 |
– |
|
339 |
– |
if (simParams_->haveCutoffRadius()) { |
340 |
– |
rCut_ = simParams_->getCutoffRadius(); |
341 |
– |
} else { |
342 |
– |
if (info_->usesElectrostaticAtoms()) { |
343 |
– |
sprintf(painCave.errMsg, |
344 |
– |
"InteractionManager::setupCutoffs: No value was set for the cutoffRadius.\n" |
345 |
– |
"\tOpenMD will use a default value of 12.0 angstroms" |
346 |
– |
"\tfor the cutoffRadius.\n"); |
347 |
– |
painCave.isFatal = 0; |
348 |
– |
painCave.severity = OPENMD_INFO; |
349 |
– |
simError(); |
350 |
– |
rCut_ = 12.0; |
351 |
– |
} else { |
352 |
– |
RealType thisCut; |
353 |
– |
set<AtomType*>::iterator i; |
354 |
– |
set<AtomType*> atomTypes; |
355 |
– |
atomTypes = info_->getSimulatedAtomTypes(); |
356 |
– |
for (i = atomTypes.begin(); i != atomTypes.end(); ++i) { |
357 |
– |
thisCut = getSuggestedCutoffRadius((*i)); |
358 |
– |
rCut_ = max(thisCut, rCut_); |
359 |
– |
} |
360 |
– |
sprintf(painCave.errMsg, |
361 |
– |
"InteractionManager::setupCutoffs: No value was set for the cutoffRadius.\n" |
362 |
– |
"\tOpenMD will use %lf angstroms.\n", |
363 |
– |
rCut_); |
364 |
– |
painCave.isFatal = 0; |
365 |
– |
painCave.severity = OPENMD_INFO; |
366 |
– |
simError(); |
367 |
– |
} |
368 |
– |
} |
369 |
– |
|
370 |
– |
map<string, CutoffMethod> stringToCutoffMethod; |
371 |
– |
stringToCutoffMethod["HARD"] = HARD; |
372 |
– |
stringToCutoffMethod["SWITCHED"] = SWITCHED; |
373 |
– |
stringToCutoffMethod["SHIFTED_POTENTIAL"] = SHIFTED_POTENTIAL; |
374 |
– |
stringToCutoffMethod["SHIFTED_FORCE"] = SHIFTED_FORCE; |
375 |
– |
|
376 |
– |
if (simParams_->haveCutoffMethod()) { |
377 |
– |
string cutMeth = toUpperCopy(simParams_->getCutoffMethod()); |
378 |
– |
map<string, CutoffMethod>::iterator i; |
379 |
– |
i = stringToCutoffMethod.find(cutMeth); |
380 |
– |
if (i == stringToCutoffMethod.end()) { |
381 |
– |
sprintf(painCave.errMsg, |
382 |
– |
"InteractionManager::setupCutoffs: Could not find chosen cutoffMethod %s\n" |
383 |
– |
"\tShould be one of: " |
384 |
– |
"HARD, SWITCHED, SHIFTED_POTENTIAL, or SHIFTED_FORCE\n", |
385 |
– |
cutMeth.c_str()); |
386 |
– |
painCave.isFatal = 1; |
387 |
– |
painCave.severity = OPENMD_ERROR; |
388 |
– |
simError(); |
389 |
– |
} else { |
390 |
– |
cutoffMethod_ = i->second; |
391 |
– |
} |
392 |
– |
} else { |
393 |
– |
sprintf(painCave.errMsg, |
394 |
– |
"InteractionManager::setupCutoffs: No value was set for the cutoffMethod.\n" |
395 |
– |
"\tOpenMD will use SHIFTED_FORCE.\n"); |
396 |
– |
painCave.isFatal = 0; |
397 |
– |
painCave.severity = OPENMD_INFO; |
398 |
– |
simError(); |
399 |
– |
cutoffMethod_ = SHIFTED_FORCE; |
400 |
– |
} |
268 |
|
} |
402 |
– |
|
403 |
– |
|
404 |
– |
/** |
405 |
– |
* setupSwitching |
406 |
– |
* |
407 |
– |
* Sets the values of switchingRadius and |
408 |
– |
* If the switchingRadius was explicitly set, use that value (but check it) |
409 |
– |
* If the switchingRadius was not explicitly set: use 0.85 * cutoffRadius_ |
410 |
– |
*/ |
411 |
– |
void InteractionManager::setupSwitching() { |
412 |
– |
Globals* simParams_ = info_->getSimParams(); |
413 |
– |
|
414 |
– |
if (simParams_->haveSwitchingRadius()) { |
415 |
– |
rSwitch_ = simParams_->getSwitchingRadius(); |
416 |
– |
if (rSwitch_ > rCut_) { |
417 |
– |
sprintf(painCave.errMsg, |
418 |
– |
"InteractionManager::setupSwitching: switchingRadius (%f) is larger than cutoffRadius(%f)\n", |
419 |
– |
rSwitch_, rCut_); |
420 |
– |
painCave.isFatal = 1; |
421 |
– |
painCave.severity = OPENMD_ERROR; |
422 |
– |
simError(); |
423 |
– |
} |
424 |
– |
} else { |
425 |
– |
rSwitch_ = 0.85 * rCut_; |
426 |
– |
sprintf(painCave.errMsg, |
427 |
– |
"InteractionManager::setupSwitching: No value was set for the switchingRadius.\n" |
428 |
– |
"\tOpenMD will use a default value of 85 percent of the cutoffRadius.\n" |
429 |
– |
"\tswitchingRadius = %f. for this simulation\n", rSwitch_); |
430 |
– |
painCave.isFatal = 0; |
431 |
– |
painCave.severity = OPENMD_WARNING; |
432 |
– |
simError(); |
433 |
– |
} |
434 |
– |
|
435 |
– |
if (simParams_->haveSwitchingFunctionType()) { |
436 |
– |
string funcType = simParams_->getSwitchingFunctionType(); |
437 |
– |
toUpper(funcType); |
438 |
– |
if (funcType == "CUBIC") { |
439 |
– |
sft_ = cubic; |
440 |
– |
} else { |
441 |
– |
if (funcType == "FIFTH_ORDER_POLYNOMIAL") { |
442 |
– |
sft_ = fifth_order_poly; |
443 |
– |
} else { |
444 |
– |
// throw error |
445 |
– |
sprintf( painCave.errMsg, |
446 |
– |
"InteractionManager::setupSwitching : Unknown switchingFunctionType. (Input file specified %s .)\n" |
447 |
– |
"\tswitchingFunctionType must be one of: " |
448 |
– |
"\"cubic\" or \"fifth_order_polynomial\".", |
449 |
– |
funcType.c_str() ); |
450 |
– |
painCave.isFatal = 1; |
451 |
– |
painCave.severity = OPENMD_ERROR; |
452 |
– |
simError(); |
453 |
– |
} |
454 |
– |
} |
455 |
– |
} |
456 |
– |
|
457 |
– |
switcher_->setSwitchType(sft_); |
458 |
– |
switcher_->setSwitch(rSwitch_, rCut_); |
459 |
– |
} |
460 |
– |
|
461 |
– |
/** |
462 |
– |
* setupNeighborlists |
463 |
– |
* |
464 |
– |
* If the skinThickness was explicitly set, use that value (but check it) |
465 |
– |
* If the skinThickness was not explicitly set: use 1.0 angstroms |
466 |
– |
*/ |
467 |
– |
void InteractionManager::setupNeighborlists() { |
468 |
– |
|
469 |
– |
Globals* simParams_ = info_->getSimParams(); |
269 |
|
|
471 |
– |
if (simParams_->haveSkinThickness()) { |
472 |
– |
skinThickness_ = simParams_->getSkinThickness(); |
473 |
– |
} else { |
474 |
– |
skinThickness_ = 1.0; |
475 |
– |
sprintf(painCave.errMsg, |
476 |
– |
"InteractionManager::setupNeighborlists: No value was set for the skinThickness.\n" |
477 |
– |
"\tOpenMD will use a default value of %f Angstroms\n" |
478 |
– |
"\tfor this simulation\n", skinThickness_); |
479 |
– |
painCave.severity = OPENMD_INFO; |
480 |
– |
painCave.isFatal = 0; |
481 |
– |
simError(); |
482 |
– |
} |
483 |
– |
|
484 |
– |
listRadius_ = rCut_ + skinThickness_; |
485 |
– |
} |
486 |
– |
|
487 |
– |
|
270 |
|
void InteractionManager::doPrePair(InteractionData idat){ |
271 |
|
|
272 |
|
if (!initialized_) initialize(); |
273 |
|
|
274 |
|
set<NonBondedInteraction*>::iterator it; |
275 |
|
|
276 |
< |
for (it = interactions_[ *(idat.atypes) ].begin(); |
277 |
< |
it != interactions_[ *(idat.atypes) ].end(); ++it){ |
276 |
> |
for (it = interactions_[ idat.atypes ].begin(); |
277 |
> |
it != interactions_[ idat.atypes ].end(); ++it){ |
278 |
|
if ((*it)->getFamily() == METALLIC_FAMILY) { |
279 |
|
dynamic_cast<MetallicInteraction*>(*it)->calcDensity(idat); |
280 |
|
} |
305 |
|
|
306 |
|
set<NonBondedInteraction*>::iterator it; |
307 |
|
|
308 |
< |
for (it = interactions_[ *(idat.atypes) ].begin(); |
309 |
< |
it != interactions_[ *(idat.atypes) ].end(); ++it) |
308 |
> |
for (it = interactions_[ idat.atypes ].begin(); |
309 |
> |
it != interactions_[ idat.atypes ].end(); ++it) |
310 |
|
(*it)->calcForce(idat); |
311 |
|
|
312 |
|
return; |
318 |
|
|
319 |
|
set<NonBondedInteraction*>::iterator it; |
320 |
|
|
321 |
< |
for (it = interactions_[ *(idat.atypes) ].begin(); |
322 |
< |
it != interactions_[ *(idat.atypes) ].end(); ++it){ |
321 |
> |
for (it = interactions_[ idat.atypes ].begin(); |
322 |
> |
it != interactions_[ idat.atypes ].end(); ++it){ |
323 |
|
if ((*it)->getFamily() == ELECTROSTATIC_FAMILY) { |
324 |
|
dynamic_cast<ElectrostaticInteraction*>(*it)->calcSkipCorrection(idat); |
325 |
|
} |