--- branches/development/src/nonbonded/InteractionManager.cpp 2010/10/02 20:41:53 1504 +++ branches/development/src/nonbonded/InteractionManager.cpp 2010/12/29 19:59:21 1532 @@ -44,10 +44,22 @@ namespace OpenMD { namespace OpenMD { bool InteractionManager::initialized_ = false; + RealType InteractionManager::rCut_ = -1.0; + RealType InteractionManager::rSwitch_ = -1.0; ForceField* InteractionManager::forceField_ = NULL; InteractionManager* InteractionManager::_instance = NULL; map InteractionManager::typeMap_; map, set > InteractionManager::interactions_; + + LJ* InteractionManager::lj_ = new LJ(); + GB* InteractionManager::gb_ = new GB(); + Sticky* InteractionManager::sticky_ = new Sticky(); + Morse* InteractionManager::morse_ = new Morse(); + EAM* InteractionManager::eam_ = new EAM(); + SC* InteractionManager::sc_ = new SC(); + Electrostatic* InteractionManager::electrostatic_ = new Electrostatic(); + MAW* InteractionManager::maw_ = new MAW(); + SwitchingFunction* InteractionManager::switcher_ = new SwitchingFunction(); InteractionManager* InteractionManager::Instance() { if (!_instance) { @@ -58,14 +70,6 @@ namespace OpenMD { void InteractionManager::initialize() { - lj_ = new LJ(); - gb_ = new GB(); - sticky_ = new Sticky(); - eam_ = new EAM(); - sc_ = new SC(); - morse_ = new Morse(); - electrostatic_ = new Electrostatic(); - lj_->setForceField(forceField_); gb_->setForceField(forceField_); sticky_->setForceField(forceField_); @@ -73,7 +77,27 @@ namespace OpenMD { sc_->setForceField(forceField_); morse_->setForceField(forceField_); electrostatic_->setForceField(forceField_); + maw_->setForceField(forceField_); + ForceFieldOptions& fopts = forceField_->getForceFieldOptions(); + // Force fields can set options on how to scale van der Waals and electrostatic + // interactions for atoms connected via bonds, bends and torsions + // in this case the topological distance between atoms is: + // 0 = the atom itself + // 1 = bonded together + // 2 = connected via a bend + // 3 = connected via a torsion + + vdwScale_[0] = 0.0; + vdwScale_[1] = fopts.getvdw12scale(); + vdwScale_[2] = fopts.getvdw13scale(); + vdwScale_[3] = fopts.getvdw14scale(); + + electrostaticScale_[0] = 0.0; + electrostaticScale_[1] = fopts.getelectrostatic12scale(); + electrostaticScale_[2] = fopts.getelectrostatic13scale(); + electrostaticScale_[3] = fopts.getelectrostatic14scale(); + ForceField::AtomTypeContainer* atomTypes = forceField_->getAtomTypes(); ForceField::AtomTypeContainer::MapTypeIterator i1, i2; AtomType* atype1; @@ -209,6 +233,27 @@ namespace OpenMD { interactions_[key].insert(sc_); metExplicit = true; } + + if (nbiType->isMAW()) { + if (vdwExplicit) { + sprintf( painCave.errMsg, + "InteractionManager::initialize found more than one explicit\n" + "\tvan der Waals interaction for atom types %s - %s\n", + atype1->getName().c_str(), atype2->getName().c_str()); + painCave.severity = OPENMD_ERROR; + painCave.isFatal = 1; + simError(); + } + // We found an explicit MAW interaction. + // override all other vdw entries for this pair of atom types: + set::iterator it; + for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it) { + InteractionFamily ifam = (*it)->getFamily(); + if (ifam == VANDERWAALS_FAMILY) interactions_[key].erase(*it); + } + interactions_[key].insert(maw_); + vdwExplicit = true; + } } } @@ -232,7 +277,7 @@ namespace OpenMD { } } - void InteractionManager::do_prepair(int *atid1, int *atid2, RealType *rij, RealType *rho_i_at_j, RealType *rho_j_at_i){ + void InteractionManager::doPrePair(int *atid1, int *atid2, RealType *rij, RealType *rho_i_at_j, RealType *rho_j_at_i){ if (!initialized_) initialize(); @@ -256,7 +301,7 @@ namespace OpenMD { return; } - void InteractionManager::do_preforce(int *atid, RealType *rho, RealType *frho, RealType *dfrhodrho){ + void InteractionManager::doPreForce(int *atid, RealType *rho, RealType *frho, RealType *dfrhodrho){ if (!initialized_) initialize(); @@ -279,7 +324,7 @@ namespace OpenMD { return; } - void InteractionManager::do_pair(int *atid1, int *atid2, RealType *d, RealType *r, RealType *r2, RealType *rcut, RealType *sw, RealType *vdwMult,RealType *electroMult, RealType *pot, RealType *vpair, RealType *f1, RealType *eFrame1, RealType *eFrame2, RealType *A1, RealType *A2, RealType *t1, RealType *t2, RealType *rho1, RealType *rho2, RealType *dfrho1, RealType *dfrho2, RealType *fshift1, RealType *fshift2){ + void InteractionManager::doPair(int *atid1, int *atid2, RealType *d, RealType *r, RealType *r2, RealType *rcut, RealType *sw, int *topoDist, RealType *pot, RealType *vpair, RealType *f1, RealType *eFrame1, RealType *eFrame2, RealType *A1, RealType *A2, RealType *t1, RealType *t2, RealType *rho1, RealType *rho2, RealType *dfrho1, RealType *dfrho2, RealType *fshift1, RealType *fshift2){ if (!initialized_) initialize(); @@ -292,8 +337,8 @@ namespace OpenMD { idat.r2 = *r2; idat.rcut = *rcut; idat.sw = *sw; - idat.vdwMult = *vdwMult; - idat.electroMult = *electroMult; + idat.vdwMult = vdwScale_[*topoDist]; + idat.electroMult = electrostaticScale_[*topoDist]; idat.pot = *pot; idat.vpair = *vpair; idat.f1 = Vector3d(f1); @@ -331,7 +376,7 @@ namespace OpenMD { return; } - void InteractionManager::do_skip_correction(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){ + 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){ if (!initialized_) initialize(); @@ -377,7 +422,7 @@ namespace OpenMD { return; } - void InteractionManager::do_self_correction(int *atid, RealType *eFrame, RealType *skippedCharge, RealType *pot, RealType *t){ + void InteractionManager::doSelfCorrection(int *atid, RealType *eFrame, RealType *skippedCharge, RealType *pot, RealType *t){ if (!initialized_) initialize(); @@ -403,8 +448,46 @@ namespace OpenMD { t[2] = scdat.t.z(); return; + } + + + RealType InteractionManager::getSuggestedCutoffRadius(int *atid) { + if (!initialized_) initialize(); + + AtomType* atype = typeMap_[*atid]; + + pair key = make_pair(atype, atype); + set::iterator it; + RealType cutoff = 0.0; + + for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it) + cutoff = max(cutoff, (*it)->getSuggestedCutoffRadius(atype, atype)); + return cutoff; + } + + RealType InteractionManager::getSuggestedCutoffRadius(AtomType* atype) { + if (!initialized_) initialize(); + + pair key = make_pair(atype, atype); + set::iterator it; + RealType cutoff = 0.0; + + for (it = interactions_[key].begin(); it != interactions_[key].end(); ++it) + cutoff = max(cutoff, (*it)->getSuggestedCutoffRadius(atype, atype)); + return cutoff; + } + + + void InteractionManager::setSwitch(RealType *rIn, RealType *rOut) { + switcher_->setSwitch(*rIn, *rOut); } + void InteractionManager::getSwitch(RealType *r2, RealType *sw, RealType *dswdr, RealType *r, + int *in_switching_region) { + bool isr = switcher_->getSwitch(*r2, *sw, *dswdr, *r); + *in_switching_region = (int)isr; + } + } //end namespace OpenMD extern "C" { @@ -415,60 +498,79 @@ extern "C" { #define fortranDoSkipCorrection FC_FUNC(do_skip_correction, DO_SKIP_CORRECTION) #define fortranDoSelfCorrection FC_FUNC(do_self_correction, DO_SELF_CORRECTION) #define fortranGetCutoff FC_FUNC(get_cutoff, GET_CUTOFF) +#define fortranSetSwitch FC_FUNC(set_switch, SET_SWITCH) +#define fortranGetSwitch FC_FUNC(get_switch, GET_SWITCH) void fortranDoPrePair(int *atid1, int *atid2, RealType *rij, RealType *rho_i_at_j, RealType *rho_j_at_i) { - return OpenMD::InteractionManager::Instance()->do_prepair(atid1, atid2, rij, - rho_i_at_j, - rho_j_at_i); + return OpenMD::InteractionManager::Instance()->doPrePair(atid1, atid2, rij, + rho_i_at_j, + rho_j_at_i); } - void fortranDoPreforce(int *atid, RealType *rho, RealType *frho, + void fortranDoPreForce(int *atid, RealType *rho, RealType *frho, RealType *dfrhodrho) { - return OpenMD::InteractionManager::Instance()->do_preforce(atid, rho, frho, - dfrhodrho); + return OpenMD::InteractionManager::Instance()->doPreForce(atid, rho, frho, + dfrhodrho); } void fortranDoPair(int *atid1, int *atid2, RealType *d, RealType *r, - RealType *r2, RealType *rcut, RealType *sw, - RealType *vdwMult, RealType *electroMult, RealType *pot, - RealType *vpair, RealType *f1, RealType *eFrame1, + RealType *r2, RealType *rcut, RealType *sw, int *topoDist, + RealType *pot, RealType *vpair, RealType *f1, RealType *eFrame1, RealType *eFrame2, RealType *A1, RealType *A2, RealType *t1, RealType *t2, RealType *rho1, RealType *rho2, RealType *dfrho1, RealType *dfrho2, RealType *fshift1, RealType *fshift2){ - return OpenMD::InteractionManager::Instance()->do_pair(atid1, atid2, d, r, r2, rcut, - sw, vdwMult, electroMult, pot, - vpair, f1, eFrame1, eFrame2, - A1, A2, t1, t2, rho1, rho2, - dfrho1, dfrho2, fshift1, fshift2); + return OpenMD::InteractionManager::Instance()->doPair(atid1, atid2, d, r, + r2, rcut, sw, topoDist, + pot, vpair, f1, + eFrame1, eFrame2, + A1, A2, t1, t2, rho1, + rho2, dfrho1, dfrho2, + fshift1, fshift2); } - - void fortranDoSkipCorrection(int *atid1, int *atid2, RealType *d, RealType *r, - RealType *skippedCharge1, RealType *skippedCharge2, - RealType *sw, RealType *electroMult, RealType *pot, + + void fortranDoSkipCorrection(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){ - return OpenMD::InteractionManager::Instance()->do_skip_correction(atid1, atid2, d, r, - skippedCharge1, - skippedCharge2, - sw, electroMult, pot, - vpair, f1, eFrame1, - eFrame2, t1, t2); + return OpenMD::InteractionManager::Instance()->doSkipCorrection(atid1, + atid2, d, + r, + skippedCharge1, + skippedCharge2, + sw, electroMult, pot, + vpair, f1, eFrame1, + eFrame2, t1, t2); } - + void fortranDoSelfCorrection(int *atid, RealType *eFrame, RealType *skippedCharge, RealType *pot, RealType *t) { - return OpenMD::InteractionManager::Instance()->do_self_correction(atid, eFrame, - skippedCharge, - pot, t); + return OpenMD::InteractionManager::Instance()->doSelfCorrection(atid, + eFrame, + skippedCharge, + pot, t); } - RealType fortranGetCutoff() { - return OpenMD::InteractionManager::Instance()->getCutoff(); + RealType fortranGetCutoff(int *atid) { + return OpenMD::InteractionManager::Instance()->getSuggestedCutoffRadius(atid); } + + void fortranGetSwitch(RealType *r2, RealType *sw, RealType *dswdr, RealType *r, + int *in_switching_region) { + + return OpenMD::InteractionManager::Instance()->getSwitch(r2, sw, dswdr, r, + in_switching_region); + } + + void fortranSetSwitch(RealType *rIn, RealType *rOut) { + return OpenMD::InteractionManager::Instance()->setSwitch(rIn, rOut); + } + }