52 |
|
namespace OpenMD { |
53 |
|
|
54 |
|
Electrostatic::Electrostatic(): name_("Electrostatic"), initialized_(false), |
55 |
< |
forceField_(NULL) {} |
55 |
> |
forceField_(NULL), info_(NULL), haveCutoffRadius_(false), |
56 |
> |
haveDampingAlpha_(false), haveDielectric_(false), |
57 |
> |
haveElectroSpline_(false) |
58 |
> |
{} |
59 |
|
|
60 |
|
void Electrostatic::initialize() { |
61 |
|
|
62 |
< |
Globals* simParams_; |
62 |
> |
Globals* simParams_ = info_->getSimParams(); |
63 |
|
|
64 |
|
summationMap_["HARD"] = esm_HARD; |
65 |
|
summationMap_["SWITCHING_FUNCTION"] = esm_SWITCHING_FUNCTION; |
100 |
|
screeningMethod_ = UNDAMPED; |
101 |
|
dielectric_ = 1.0; |
102 |
|
one_third_ = 1.0 / 3.0; |
100 |
– |
haveCutoffRadius_ = false; |
101 |
– |
haveDampingAlpha_ = false; |
102 |
– |
haveDielectric_ = false; |
103 |
– |
haveElectroSpline_ = false; |
103 |
|
|
104 |
|
// check the summation method: |
105 |
|
if (simParams_->haveElectrostaticSummationMethod()) { |
406 |
|
return; |
407 |
|
} |
408 |
|
|
409 |
< |
void Electrostatic::setElectrostaticCutoffRadius( RealType theECR, |
410 |
< |
RealType theRSW ) { |
412 |
< |
cutoffRadius_ = theECR; |
409 |
> |
void Electrostatic::setCutoffRadius( RealType rCut ) { |
410 |
> |
cutoffRadius_ = rCut; |
411 |
|
rrf_ = cutoffRadius_; |
414 |
– |
rt_ = theRSW; |
412 |
|
haveCutoffRadius_ = true; |
413 |
+ |
} |
414 |
+ |
|
415 |
+ |
void Electrostatic::setSwitchingRadius( RealType rSwitch ) { |
416 |
+ |
rt_ = rSwitch; |
417 |
|
} |
418 |
|
void Electrostatic::setElectrostaticSummationMethod( ElectrostaticSummationMethod esm ) { |
419 |
|
summationMethod_ = esm; |
1035 |
|
RealType mu1, preVal, chg1, self; |
1036 |
|
|
1037 |
|
if (!initialized_) initialize(); |
1038 |
< |
|
1038 |
> |
|
1039 |
|
ElectrostaticAtomData data = ElectrostaticMap[sdat.atype]; |
1040 |
|
|
1041 |
|
// logicals |
1041 |
– |
|
1042 |
|
bool i_is_Charge = data.is_Charge; |
1043 |
|
bool i_is_Dipole = data.is_Dipole; |
1044 |
|
|
1046 |
|
if (i_is_Dipole) { |
1047 |
|
mu1 = data.dipole_moment; |
1048 |
|
preVal = pre22_ * preRF2_ * mu1 * mu1; |
1049 |
< |
sdat.pot[2] -= 0.5 * preVal; |
1049 |
> |
(*(sdat.pot))[ELECTROSTATIC_FAMILY] -= 0.5 * preVal; |
1050 |
|
|
1051 |
|
// The self-correction term adds into the reaction field vector |
1052 |
|
Vector3d uz_i = sdat.eFrame->getColumn(2); |
1063 |
|
} else { |
1064 |
|
self = - 0.5 * rcuti_ * chg1 * (chg1 + *(sdat.skippedCharge)) * pre11_; |
1065 |
|
} |
1066 |
< |
sdat.pot[ELECTROSTATIC_FAMILY] += self; |
1066 |
> |
(*(sdat.pot))[ELECTROSTATIC_FAMILY] += self; |
1067 |
|
} |
1068 |
|
} |
1069 |
|
} |