| 40 |
|
*/ |
| 41 |
|
|
| 42 |
|
#include "nonbonded/InteractionManager.hpp" |
| 43 |
– |
#include "UseTheForce/doForces_interface.h" |
| 43 |
|
|
| 44 |
|
namespace OpenMD { |
| 45 |
|
|
| 49 |
|
|
| 50 |
|
RealType InteractionManager::rCut_ = 0.0; |
| 51 |
|
RealType InteractionManager::rSwitch_ = 0.0; |
| 53 |
– |
RealType InteractionManager::skinThickness_ = 0.0; |
| 54 |
– |
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}; |
| 303 |
|
|
| 304 |
|
setupCutoffs(); |
| 305 |
|
setupSwitching(); |
| 309 |
– |
setupNeighborlists(); |
| 306 |
|
|
| 307 |
< |
int ljsp = cutoffMethod_ == SHIFTED_POTENTIAL ? 1 : 0; |
| 308 |
< |
int ljsf = cutoffMethod_ == SHIFTED_FORCE ? 1 : 0; |
| 309 |
< |
notifyFortranCutoffs(&rCut_, &rSwitch_, &ljsp, &ljsf); |
| 314 |
< |
notifyFortranSkinThickness(&skinThickness_); |
| 307 |
> |
//int ljsp = cutoffMethod_ == SHIFTED_POTENTIAL ? 1 : 0; |
| 308 |
> |
//int ljsf = cutoffMethod_ == SHIFTED_FORCE ? 1 : 0; |
| 309 |
> |
//notifyFortranCutoffs(&rCut_, &rSwitch_, &ljsp, &ljsf); |
| 310 |
|
|
| 311 |
|
initialized_ = true; |
| 312 |
|
} |
| 454 |
|
switcher_->setSwitch(rSwitch_, rCut_); |
| 455 |
|
} |
| 456 |
|
|
| 457 |
< |
/** |
| 463 |
< |
* setupNeighborlists |
| 464 |
< |
* |
| 465 |
< |
* If the skinThickness was explicitly set, use that value (but check it) |
| 466 |
< |
* If the skinThickness was not explicitly set: use 1.0 angstroms |
| 467 |
< |
*/ |
| 468 |
< |
void InteractionManager::setupNeighborlists() { |
| 469 |
< |
|
| 470 |
< |
Globals* simParams_ = info_->getSimParams(); |
| 471 |
< |
|
| 472 |
< |
if (simParams_->haveSkinThickness()) { |
| 473 |
< |
skinThickness_ = simParams_->getSkinThickness(); |
| 474 |
< |
} else { |
| 475 |
< |
skinThickness_ = 1.0; |
| 476 |
< |
sprintf(painCave.errMsg, |
| 477 |
< |
"InteractionManager::setupNeighborlists: No value was set for the skinThickness.\n" |
| 478 |
< |
"\tOpenMD will use a default value of %f Angstroms\n" |
| 479 |
< |
"\tfor this simulation\n", skinThickness_); |
| 480 |
< |
painCave.severity = OPENMD_INFO; |
| 481 |
< |
painCave.isFatal = 0; |
| 482 |
< |
simError(); |
| 483 |
< |
} |
| 484 |
< |
|
| 485 |
< |
listRadius_ = rCut_ + skinThickness_; |
| 486 |
< |
} |
| 487 |
< |
|
| 488 |
< |
|
| 489 |
< |
void InteractionManager::doPrePair(int *atid1, int *atid2, RealType *rij, RealType *rho_i_at_j, RealType *rho_j_at_i){ |
| 457 |
> |
void InteractionManager::doPrePair(InteractionData idat){ |
| 458 |
|
|
| 459 |
|
if (!initialized_) initialize(); |
| 460 |
< |
|
| 493 |
< |
DensityData ddat; |
| 494 |
< |
|
| 495 |
< |
ddat.atype1 = typeMap_[*atid1]; |
| 496 |
< |
ddat.atype2 = typeMap_[*atid2]; |
| 497 |
< |
ddat.rij = *rij; |
| 498 |
< |
ddat.rho_i_at_j = *rho_i_at_j; |
| 499 |
< |
ddat.rho_j_at_i = *rho_j_at_i; |
| 500 |
< |
|
| 501 |
< |
pair<AtomType*, AtomType*> key = make_pair(ddat.atype1, ddat.atype2); |
| 460 |
> |
|
| 461 |
|
set<NonBondedInteraction*>::iterator it; |
| 462 |
|
|
| 463 |
< |
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it){ |
| 463 |
> |
for (it = interactions_[ *(idat.atypes) ].begin(); |
| 464 |
> |
it != interactions_[ *(idat.atypes) ].end(); ++it){ |
| 465 |
|
if ((*it)->getFamily() == METALLIC_FAMILY) { |
| 466 |
< |
dynamic_cast<MetallicInteraction*>(*it)->calcDensity(ddat); |
| 466 |
> |
dynamic_cast<MetallicInteraction*>(*it)->calcDensity(idat); |
| 467 |
|
} |
| 468 |
|
} |
| 469 |
|
|
| 470 |
|
return; |
| 471 |
|
} |
| 472 |
|
|
| 473 |
< |
void InteractionManager::doPreForce(int *atid, RealType *rho, RealType *frho, RealType *dfrhodrho){ |
| 473 |
> |
void InteractionManager::doPreForce(SelfData sdat){ |
| 474 |
|
|
| 475 |
|
if (!initialized_) initialize(); |
| 476 |
< |
|
| 477 |
< |
FunctionalData fdat; |
| 518 |
< |
|
| 519 |
< |
fdat.atype = typeMap_[*atid]; |
| 520 |
< |
fdat.rho = *rho; |
| 521 |
< |
fdat.frho = *frho; |
| 522 |
< |
fdat.dfrhodrho = *dfrhodrho; |
| 523 |
< |
|
| 524 |
< |
pair<AtomType*, AtomType*> key = make_pair(fdat.atype, fdat.atype); |
| 476 |
> |
|
| 477 |
> |
pair<AtomType*, AtomType*> key = make_pair(sdat.atype, sdat.atype); |
| 478 |
|
set<NonBondedInteraction*>::iterator it; |
| 479 |
|
|
| 480 |
|
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it){ |
| 481 |
|
if ((*it)->getFamily() == METALLIC_FAMILY) { |
| 482 |
< |
dynamic_cast<MetallicInteraction*>(*it)->calcFunctional(fdat); |
| 482 |
> |
dynamic_cast<MetallicInteraction*>(*it)->calcFunctional(sdat); |
| 483 |
|
} |
| 484 |
|
} |
| 485 |
|
|
| 486 |
|
return; |
| 487 |
|
} |
| 488 |
|
|
| 489 |
< |
void InteractionManager::doPair(int *atid1, int *atid2, RealType *d, RealType *r, RealType *r2, RealType *sw, int *topoDist, RealType *A1, RealType *A2, RealType *eFrame1, RealType *eFrame2, RealType *vpair, RealType *pot, RealType *f1, RealType *t1, RealType *t2, RealType *rho1, RealType *rho2, RealType *dfrho1, RealType *dfrho2, RealType *fshift1, RealType *fshift2){ |
| 489 |
> |
void InteractionManager::doPair(InteractionData idat){ |
| 490 |
|
|
| 491 |
|
if (!initialized_) initialize(); |
| 492 |
< |
|
| 540 |
< |
InteractionData idat; |
| 541 |
< |
|
| 542 |
< |
idat.atype1 = typeMap_[*atid1]; |
| 543 |
< |
idat.atype2 = typeMap_[*atid2]; |
| 544 |
< |
idat.d = Vector3d(d); |
| 545 |
< |
idat.rij = *r; |
| 546 |
< |
idat.r2 = *r2; |
| 547 |
< |
idat.sw = *sw; |
| 548 |
< |
idat.vdwMult = vdwScale_[*topoDist]; |
| 549 |
< |
idat.electroMult = electrostaticScale_[*topoDist]; |
| 550 |
< |
for (int i = 0; i < 4; i++) { |
| 551 |
< |
idat.vpair[i] = vpair[i]; |
| 552 |
< |
idat.pot[i] = pot[i]; |
| 553 |
< |
} |
| 554 |
< |
idat.f1 = Vector3d(f1[0], f1[1], f1[2]); |
| 555 |
< |
idat.eFrame1 = Mat3x3d(eFrame1); |
| 556 |
< |
idat.eFrame2 = Mat3x3d(eFrame2); |
| 557 |
< |
idat.A1 = RotMat3x3d(A1); |
| 558 |
< |
idat.A2 = RotMat3x3d(A2); |
| 559 |
< |
idat.t1 = Vector3d(t1); |
| 560 |
< |
idat.t2 = Vector3d(t2); |
| 561 |
< |
idat.rho1 = *rho1; |
| 562 |
< |
idat.rho2 = *rho2; |
| 563 |
< |
idat.dfrho1 = *dfrho1; |
| 564 |
< |
idat.dfrho2 = *dfrho2; |
| 565 |
< |
idat.fshift1 = *fshift1; |
| 566 |
< |
idat.fshift2 = *fshift2; |
| 567 |
< |
|
| 568 |
< |
pair<AtomType*, AtomType*> key = make_pair(idat.atype1, idat.atype2); |
| 492 |
> |
|
| 493 |
|
set<NonBondedInteraction*>::iterator it; |
| 494 |
|
|
| 495 |
< |
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it) |
| 495 |
> |
for (it = interactions_[ *(idat.atypes) ].begin(); |
| 496 |
> |
it != interactions_[ *(idat.atypes) ].end(); ++it) |
| 497 |
|
(*it)->calcForce(idat); |
| 498 |
|
|
| 574 |
– |
for (int i = 0; i < 4; i++) { |
| 575 |
– |
vpair[i] = idat.vpair[i]; |
| 576 |
– |
pot[i] = idat.pot[i]; |
| 577 |
– |
} |
| 578 |
– |
|
| 579 |
– |
for (int i = 0; i < 3; i++) { |
| 580 |
– |
f1[i] = idat.f1[i]; |
| 581 |
– |
t1[i] = idat.t1[i]; |
| 582 |
– |
t2[i] = idat.t2[i]; |
| 583 |
– |
} |
| 584 |
– |
|
| 499 |
|
return; |
| 500 |
|
} |
| 501 |
|
|
| 502 |
< |
void InteractionManager::doSkipCorrection(int *atid1, int *atid2, RealType *d, RealType *r, RealType *skippedCharge1, RealType *skippedCharge2, RealType *sw, RealType *electroMult, RealType *pot, RealType *vpair, RealType *f1, RealType *eFrame1, RealType *eFrame2, RealType *t1, RealType *t2){ |
| 502 |
> |
void InteractionManager::doSkipCorrection(InteractionData idat){ |
| 503 |
|
|
| 504 |
< |
if (!initialized_) initialize(); |
| 504 |
> |
if (!initialized_) initialize(); |
| 505 |
|
|
| 592 |
– |
SkipCorrectionData skdat; |
| 593 |
– |
|
| 594 |
– |
skdat.atype1 = typeMap_[*atid1]; |
| 595 |
– |
skdat.atype2 = typeMap_[*atid2]; |
| 596 |
– |
skdat.d = Vector3d(d); |
| 597 |
– |
skdat.rij = *r; |
| 598 |
– |
skdat.skippedCharge1 = *skippedCharge1; |
| 599 |
– |
skdat.skippedCharge2 = *skippedCharge2; |
| 600 |
– |
skdat.sw = *sw; |
| 601 |
– |
skdat.electroMult = *electroMult; |
| 602 |
– |
for (int i = 0; i < 4; i++) { |
| 603 |
– |
skdat.vpair[i] = vpair[i]; |
| 604 |
– |
skdat.pot[i] = pot[i]; |
| 605 |
– |
} |
| 606 |
– |
skdat.f1 = Vector3d(f1); |
| 607 |
– |
skdat.eFrame1 = Mat3x3d(eFrame1); |
| 608 |
– |
skdat.eFrame2 = Mat3x3d(eFrame2); |
| 609 |
– |
skdat.t1 = Vector3d(t1); |
| 610 |
– |
skdat.t2 = Vector3d(t2); |
| 611 |
– |
|
| 612 |
– |
pair<AtomType*, AtomType*> key = make_pair(skdat.atype1, skdat.atype2); |
| 506 |
|
set<NonBondedInteraction*>::iterator it; |
| 507 |
|
|
| 508 |
< |
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it){ |
| 508 |
> |
for (it = interactions_[ *(idat.atypes) ].begin(); |
| 509 |
> |
it != interactions_[ *(idat.atypes) ].end(); ++it){ |
| 510 |
|
if ((*it)->getFamily() == ELECTROSTATIC_FAMILY) { |
| 511 |
< |
dynamic_cast<ElectrostaticInteraction*>(*it)->calcSkipCorrection(skdat); |
| 511 |
> |
dynamic_cast<ElectrostaticInteraction*>(*it)->calcSkipCorrection(idat); |
| 512 |
|
} |
| 513 |
|
} |
| 514 |
|
|
| 621 |
– |
for (int i = 0; i < 4; i++) { |
| 622 |
– |
vpair[i] = skdat.vpair[i]; |
| 623 |
– |
pot[i] = skdat.pot[i]; |
| 624 |
– |
} |
| 625 |
– |
|
| 626 |
– |
for (int i = 0; i < 3; i++) { |
| 627 |
– |
f1[i] = skdat.f1[i]; |
| 628 |
– |
t1[i] = skdat.t1[i]; |
| 629 |
– |
t2[i] = skdat.t2[i]; |
| 630 |
– |
} |
| 631 |
– |
|
| 515 |
|
return; |
| 516 |
|
} |
| 517 |
|
|
| 518 |
< |
void InteractionManager::doSelfCorrection(int *atid, RealType *eFrame, RealType *skippedCharge, RealType *pot, RealType *t){ |
| 518 |
> |
void InteractionManager::doSelfCorrection(SelfData sdat){ |
| 519 |
|
|
| 520 |
|
if (!initialized_) initialize(); |
| 521 |
|
|
| 522 |
< |
SelfCorrectionData scdat; |
| 640 |
< |
|
| 641 |
< |
scdat.atype = typeMap_[*atid]; |
| 642 |
< |
scdat.eFrame = Mat3x3d(eFrame); |
| 643 |
< |
scdat.skippedCharge = *skippedCharge; |
| 644 |
< |
for (int i = 0; i < 4; i++){ |
| 645 |
< |
scdat.pot[i] = pot[i]; |
| 646 |
< |
} |
| 647 |
< |
scdat.t = Vector3d(t); |
| 648 |
< |
|
| 649 |
< |
pair<AtomType*, AtomType*> key = make_pair(scdat.atype, scdat.atype); |
| 522 |
> |
pair<AtomType*, AtomType*> key = make_pair(sdat.atype, sdat.atype); |
| 523 |
|
set<NonBondedInteraction*>::iterator it; |
| 524 |
|
|
| 525 |
|
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it){ |
| 526 |
|
if ((*it)->getFamily() == ELECTROSTATIC_FAMILY) { |
| 527 |
< |
dynamic_cast<ElectrostaticInteraction*>(*it)->calcSelfCorrection(scdat); |
| 527 |
> |
dynamic_cast<ElectrostaticInteraction*>(*it)->calcSelfCorrection(sdat); |
| 528 |
|
} |
| 529 |
|
} |
| 530 |
|
|
| 658 |
– |
t[0] = scdat.t.x(); |
| 659 |
– |
t[1] = scdat.t.y(); |
| 660 |
– |
t[2] = scdat.t.z(); |
| 661 |
– |
|
| 531 |
|
return; |
| 532 |
|
} |
| 533 |
|
|
| 665 |
– |
|
| 534 |
|
RealType InteractionManager::getSuggestedCutoffRadius(int *atid) { |
| 535 |
|
if (!initialized_) initialize(); |
| 536 |
|
|
| 541 |
|
RealType cutoff = 0.0; |
| 542 |
|
|
| 543 |
|
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it) |
| 544 |
< |
cutoff = max(cutoff, (*it)->getSuggestedCutoffRadius(atype, atype)); |
| 544 |
> |
cutoff = max(cutoff, (*it)->getSuggestedCutoffRadius(key)); |
| 545 |
|
return cutoff; |
| 546 |
|
} |
| 547 |
|
|
| 553 |
|
RealType cutoff = 0.0; |
| 554 |
|
|
| 555 |
|
for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it) |
| 556 |
< |
cutoff = max(cutoff, (*it)->getSuggestedCutoffRadius(atype, atype)); |
| 556 |
> |
cutoff = max(cutoff, (*it)->getSuggestedCutoffRadius(key)); |
| 557 |
|
return cutoff; |
| 558 |
|
} |
| 691 |
– |
|
| 692 |
– |
void InteractionManager::getSwitch(RealType *r2, RealType *sw, RealType *dswdr, RealType *r, |
| 693 |
– |
int *in_switching_region) { |
| 694 |
– |
bool isr = switcher_->getSwitch( *r2 , *sw, *dswdr, *r); |
| 695 |
– |
*in_switching_region = (int)isr; |
| 696 |
– |
} |
| 697 |
– |
|
| 559 |
|
} //end namespace OpenMD |
| 699 |
– |
|
| 700 |
– |
extern "C" { |
| 701 |
– |
|
| 702 |
– |
#define fortranDoPrePair FC_FUNC(do_prepair, DO_PREPAIR) |
| 703 |
– |
#define fortranDoPreForce FC_FUNC(do_preforce, DO_PREFORCE) |
| 704 |
– |
#define fortranDoPair FC_FUNC(do_pair, DO_PAIR) |
| 705 |
– |
#define fortranDoSkipCorrection FC_FUNC(do_skip_correction, DO_SKIP_CORRECTION) |
| 706 |
– |
#define fortranDoSelfCorrection FC_FUNC(do_self_correction, DO_SELF_CORRECTION) |
| 707 |
– |
#define fortranGetCutoff FC_FUNC(get_cutoff, GET_CUTOFF) |
| 708 |
– |
#define fortranGetSwitch FC_FUNC(get_switch, GET_SWITCH) |
| 709 |
– |
|
| 710 |
– |
void fortranDoPrePair(int *atid1, int *atid2, RealType *rij, |
| 711 |
– |
RealType *rho_i_at_j, RealType *rho_j_at_i) { |
| 712 |
– |
|
| 713 |
– |
return OpenMD::InteractionManager::Instance()->doPrePair(atid1, atid2, rij, |
| 714 |
– |
rho_i_at_j, |
| 715 |
– |
rho_j_at_i); |
| 716 |
– |
} |
| 717 |
– |
void fortranDoPreForce(int *atid, RealType *rho, RealType *frho, |
| 718 |
– |
RealType *dfrhodrho) { |
| 719 |
– |
|
| 720 |
– |
return OpenMD::InteractionManager::Instance()->doPreForce(atid, rho, frho, |
| 721 |
– |
dfrhodrho); |
| 722 |
– |
} |
| 723 |
– |
|
| 724 |
– |
void fortranDoPair(int *atid1, int *atid2, RealType *d, RealType *r, RealType *r2, |
| 725 |
– |
RealType *sw, int *topoDist, RealType *A1, RealType *A2, |
| 726 |
– |
RealType *eFrame1, RealType *eFrame2, |
| 727 |
– |
RealType *vpair, RealType *pot, |
| 728 |
– |
RealType *f1, RealType *t1, RealType *t2, |
| 729 |
– |
RealType *rho1, RealType *rho2, RealType *dfrho1, RealType *dfrho2, |
| 730 |
– |
RealType *fshift1, RealType *fshift2){ |
| 731 |
– |
|
| 732 |
– |
return OpenMD::InteractionManager::Instance()->doPair(atid1, atid2, d, r, r2, |
| 733 |
– |
sw, topoDist, A1, A2, |
| 734 |
– |
eFrame1, eFrame2, |
| 735 |
– |
vpair, pot, |
| 736 |
– |
f1, t1, t2, |
| 737 |
– |
rho1, rho2, |
| 738 |
– |
dfrho1, dfrho2, |
| 739 |
– |
fshift1, fshift2); |
| 740 |
– |
} |
| 741 |
– |
|
| 742 |
– |
void fortranDoSkipCorrection(int *atid1, int *atid2, RealType *d, |
| 743 |
– |
RealType *r, RealType *skippedCharge1, |
| 744 |
– |
RealType *skippedCharge2, RealType *sw, |
| 745 |
– |
RealType *electroMult, RealType *pot, |
| 746 |
– |
RealType *vpair, RealType *f1, |
| 747 |
– |
RealType *eFrame1, RealType *eFrame2, |
| 748 |
– |
RealType *t1, RealType *t2){ |
| 749 |
– |
|
| 750 |
– |
return OpenMD::InteractionManager::Instance()->doSkipCorrection(atid1, |
| 751 |
– |
atid2, d, |
| 752 |
– |
r, |
| 753 |
– |
skippedCharge1, |
| 754 |
– |
skippedCharge2, |
| 755 |
– |
sw, electroMult, pot, |
| 756 |
– |
vpair, f1, eFrame1, |
| 757 |
– |
eFrame2, t1, t2); |
| 758 |
– |
} |
| 759 |
– |
|
| 760 |
– |
void fortranDoSelfCorrection(int *atid, RealType *eFrame, RealType *skippedCharge, |
| 761 |
– |
RealType *pot, RealType *t) { |
| 762 |
– |
|
| 763 |
– |
return OpenMD::InteractionManager::Instance()->doSelfCorrection(atid, |
| 764 |
– |
eFrame, |
| 765 |
– |
skippedCharge, |
| 766 |
– |
pot, t); |
| 767 |
– |
} |
| 768 |
– |
RealType fortranGetCutoff(int *atid) { |
| 769 |
– |
return OpenMD::InteractionManager::Instance()->getSuggestedCutoffRadius(atid); |
| 770 |
– |
} |
| 771 |
– |
|
| 772 |
– |
void fortranGetSwitch(RealType *r2, RealType *sw, RealType *dswdr, RealType *r, |
| 773 |
– |
int *in_switching_region) { |
| 774 |
– |
return OpenMD::InteractionManager::Instance()->getSwitch(r2, sw, dswdr, r, |
| 775 |
– |
in_switching_region); |
| 776 |
– |
} |
| 777 |
– |
} |