45 |
|
#include <cmath> |
46 |
|
#include "nonbonded/Morse.hpp" |
47 |
|
#include "utils/simError.h" |
48 |
< |
#include "types/NonBondedInteractionType.hpp" |
48 |
> |
#include "types/MorseInteractionType.hpp" |
49 |
|
|
50 |
|
using namespace std; |
51 |
|
|
55 |
|
|
56 |
|
void Morse::initialize() { |
57 |
|
|
58 |
– |
stringToEnumMap_["shiftedMorse"] = shiftedMorse; |
59 |
– |
stringToEnumMap_["repulsiveMorse"] = repulsiveMorse; |
60 |
– |
|
58 |
|
ForceField::NonBondedInteractionTypeContainer* nbiTypes = forceField_->getNonBondedInteractionTypes(); |
59 |
|
ForceField::NonBondedInteractionTypeContainer::MapTypeIterator j; |
60 |
|
NonBondedInteractionType* nbt; |
61 |
+ |
ForceField::NonBondedInteractionTypeContainer::KeyType keys; |
62 |
|
|
63 |
|
for (nbt = nbiTypes->beginType(j); nbt != NULL; |
64 |
|
nbt = nbiTypes->nextType(j)) { |
65 |
|
|
66 |
|
if (nbt->isMorse()) { |
67 |
< |
|
68 |
< |
pair<AtomType*, AtomType*> atypes = nbt->getAtomTypes(); |
69 |
< |
|
72 |
< |
GenericData* data = nbt->getPropertyByName("Morse"); |
73 |
< |
if (data == NULL) { |
74 |
< |
sprintf( painCave.errMsg, "Morse::initialize could not find\n" |
75 |
< |
"\tMorse parameters for %s - %s interaction.\n", |
76 |
< |
atypes.first->getName().c_str(), |
77 |
< |
atypes.second->getName().c_str()); |
78 |
< |
painCave.severity = OPENMD_ERROR; |
79 |
< |
painCave.isFatal = 1; |
80 |
< |
simError(); |
81 |
< |
} |
82 |
< |
|
83 |
< |
MorseData* morseData = dynamic_cast<MorseData*>(data); |
84 |
< |
if (morseData == NULL) { |
85 |
< |
sprintf( painCave.errMsg, |
86 |
< |
"Morse::initialize could not convert GenericData to\n" |
87 |
< |
"\tMorseData for %s - %s interaction.\n", |
88 |
< |
atypes.first->getName().c_str(), |
89 |
< |
atypes.second->getName().c_str()); |
90 |
< |
painCave.severity = OPENMD_ERROR; |
91 |
< |
painCave.isFatal = 1; |
92 |
< |
simError(); |
93 |
< |
} |
94 |
< |
|
95 |
< |
MorseParam morseParam = morseData->getData(); |
67 |
> |
keys = nbiTypes->getKeys(j); |
68 |
> |
AtomType* at1 = forceField_->getAtomType(keys[0]); |
69 |
> |
AtomType* at2 = forceField_->getAtomType(keys[1]); |
70 |
|
|
71 |
< |
RealType De = morseParam.De; |
98 |
< |
RealType Re = morseParam.Re; |
99 |
< |
RealType beta = morseParam.beta; |
100 |
< |
string interactionType = morseParam.interactionType; |
71 |
> |
MorseInteractionType* mit = dynamic_cast<MorseInteractionType*>(nbt); |
72 |
|
|
73 |
< |
toUpper(interactionType); |
103 |
< |
map<string, MorseInteractionType>::iterator i; |
104 |
< |
i = stringToEnumMap_.find(interactionType); |
105 |
< |
if (i != stringToEnumMap_.end()) { |
106 |
< |
addExplicitInteraction(atypes.first, atypes.second, |
107 |
< |
De, Re, beta, i->second ); |
108 |
< |
} else { |
73 |
> |
if (mit == NULL) { |
74 |
|
sprintf( painCave.errMsg, |
75 |
< |
"Morse::initialize found unknown Morse interaction type\n" |
76 |
< |
"\t(%s) for %s - %s interaction.\n", |
77 |
< |
morseParam.interactionType.c_str(), |
78 |
< |
atypes.first->getName().c_str(), |
114 |
< |
atypes.second->getName().c_str()); |
75 |
> |
"Morse::initialize could not convert NonBondedInteractionType\n" |
76 |
> |
"\tto MorseInteractionType for %s - %s interaction.\n", |
77 |
> |
at1->getName().c_str(), |
78 |
> |
at2->getName().c_str()); |
79 |
|
painCave.severity = OPENMD_ERROR; |
80 |
|
painCave.isFatal = 1; |
81 |
|
simError(); |
82 |
|
} |
83 |
+ |
|
84 |
+ |
RealType De = mit->getD(); |
85 |
+ |
RealType Re = mit->getR(); |
86 |
+ |
RealType beta = mit->getBeta(); |
87 |
+ |
|
88 |
+ |
MorseType variant = mit->getInteractionType(); |
89 |
+ |
addExplicitInteraction(at1, at2, De, Re, beta, variant ); |
90 |
|
} |
91 |
|
} |
92 |
|
initialized_ = true; |
94 |
|
|
95 |
|
void Morse::addExplicitInteraction(AtomType* atype1, AtomType* atype2, |
96 |
|
RealType De, RealType Re, RealType beta, |
97 |
< |
MorseInteractionType mit) { |
97 |
> |
MorseType mt) { |
98 |
|
|
99 |
|
MorseInteractionData mixer; |
100 |
|
mixer.De = De; |
101 |
|
mixer.Re = Re; |
102 |
|
mixer.beta = beta; |
103 |
< |
mixer.interactionType = mit; |
103 |
> |
mixer.variant = mt; |
104 |
|
|
105 |
|
pair<AtomType*, AtomType*> key1, key2; |
106 |
|
key1 = make_pair(atype1, atype2); |
129 |
|
RealType De = mixer.De; |
130 |
|
RealType Re = mixer.Re; |
131 |
|
RealType beta = mixer.beta; |
132 |
< |
MorseInteractionType interactionType = mixer.interactionType; |
132 |
> |
MorseType variant = mixer.variant; |
133 |
|
|
134 |
|
// V(r) = D_e exp(-a(r-re)(exp(-a(r-re))-2) |
135 |
|
|
148 |
|
} |
149 |
|
|
150 |
|
|
151 |
< |
switch(interactionType) { |
152 |
< |
case shiftedMorse : { |
151 |
> |
switch(variant) { |
152 |
> |
case mtShifted : { |
153 |
|
|
154 |
|
myPot = De * (expfnc2 - 2.0 * expfnc); |
155 |
|
myDeriv = 2.0 * De * beta * (expfnc - expfnc2); |
168 |
|
|
169 |
|
break; |
170 |
|
} |
171 |
< |
case repulsiveMorse : { |
171 |
> |
case mtRepulsive : { |
172 |
|
|
173 |
|
myPot = De * expfnc2; |
174 |
|
myDeriv = -2.0 * De * beta * expfnc2; |
187 |
|
|
188 |
|
break; |
189 |
|
} |
190 |
+ |
case mtUnknown: { |
191 |
+ |
// don't know what to do so don't do anything |
192 |
+ |
break; |
193 |
|
} |
194 |
+ |
} |
195 |
|
|
196 |
|
RealType pot_temp = *(idat.vdwMult) * (myPot - myPotC); |
197 |
|
*(idat.vpair) += pot_temp; |