--- branches/development/src/nonbonded/InteractionManager.cpp 2011/07/08 20:25:32 1587 +++ branches/development/src/nonbonded/InteractionManager.cpp 2013/01/09 19:27:52 1825 @@ -36,7 +36,8 @@ * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). - * [4] Vardeman & Gezelter, in progress (2009). + * [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). + * [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). */ #include "nonbonded/InteractionManager.hpp" @@ -51,6 +52,7 @@ namespace OpenMD { gb_ = new GB(); sticky_ = new Sticky(); morse_ = new Morse(); + repulsivePower_ = new RepulsivePower(); eam_ = new EAM(); sc_ = new SC(); electrostatic_ = new Electrostatic(); @@ -58,7 +60,9 @@ namespace OpenMD { } void InteractionManager::initialize() { - + + if (initialized_) return; + ForceField* forceField_ = info_->getForceField(); lj_->setForceField(forceField_); @@ -70,6 +74,7 @@ namespace OpenMD { electrostatic_->setSimInfo(info_); electrostatic_->setForceField(forceField_); maw_->setForceField(forceField_); + repulsivePower_->setForceField(forceField_); ForceField::AtomTypeContainer* atomTypes = forceField_->getAtomTypes(); ForceField::AtomTypeContainer::MapTypeIterator i1, i2; @@ -82,14 +87,13 @@ namespace OpenMD { atype1 = atomTypes->nextType(i1)) { // add it to the map: - AtomTypeProperties atp = atype1->getATP(); pair::iterator,bool> ret; - ret = typeMap_.insert( pair(atp.ident, atype1) ); + ret = typeMap_.insert( pair(atype1->getIdent(), atype1) ); if (ret.second == false) { sprintf( painCave.errMsg, "InteractionManager already had a previous entry with ident %d\n", - atp.ident); + atype1->getIdent()); painCave.severity = OPENMD_INFO; painCave.isFatal = 0; simError(); @@ -107,7 +111,7 @@ namespace OpenMD { bool vdwExplicit = false; bool metExplicit = false; - bool hbExplicit = false; + // bool hbExplicit = false; key = make_pair(atype1, atype2); @@ -176,9 +180,33 @@ namespace OpenMD { if (ifam == VANDERWAALS_FAMILY) interactions_[key].erase(*it); } interactions_[key].insert(morse_); + vdwExplicit = true; + } + + if (nbiType->isRepulsivePower()) { + 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 RepulsivePower 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(repulsivePower_); vdwExplicit = true; } + if (nbiType->isEAM()) { // We found an explicit EAM interaction. // override all other metallic entries for this pair of atom types: @@ -276,10 +304,6 @@ namespace OpenMD { eam_->setCutoffRadius(rcut); } - void InteractionManager::setSwitchingRadius(RealType rswitch) { - electrostatic_->setSwitchingRadius(rswitch); - } - void InteractionManager::doPrePair(InteractionData idat){ if (!initialized_) initialize(); @@ -353,7 +377,7 @@ namespace OpenMD { RealType InteractionManager::getSuggestedCutoffRadius(int *atid) { if (!initialized_) initialize(); - + AtomType* atype = typeMap_[*atid]; pair key = make_pair(atype, atype);