-
Notifications
You must be signed in to change notification settings - Fork 4
/
Copy pathtelectr.cpp
52 lines (43 loc) · 2.06 KB
/
telectr.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
#include "telectr.h"
//---------------------------------------------------------------------------------------------------
void TElectr::Scale14()
{
dQQ /= 1.2;
}
//---------------------------------------------------------------------------------------------------
double TElectr::dCalcEnergy() const
{
// Some polarizable medium. Does influence energy, but not gradients:
if(TConsts::USEBORN)
{
return dQQ/(sqrt(dLength_SQR())*TConsts::EpsilonRNA) +
-dQQ*(1./TConsts::EpsilonRNA - 1./TConsts::EpsilonWater)/dEffectiveBornDistance();
}
else return dQQ/(sqrt(dLength_SQR())*TConsts::EpsilonRNA); // Simplified approach.
}
//---------------------------------------------------------------------------------------------------
void TElectr::CalcGradients() const
{
TVector vec12 = vecMakeVector(patAtom1, patAtom2);
double dLen_SQR = vec12.dLength_SQR();
vec12 *= dQQ/(dLen_SQR*sqrt(dLen_SQR)*TConsts::EpsilonRNA);
patAtom1->vecEnergyGradient += vec12;
patAtom2->vecEnergyGradient -= vec12;
}
//---------------------------------------------------------------------------------------------------
double TElectr::dEffectiveBornDistance() const
{
double dLen_SQR = dLength_SQR(),
R1R2 = (patAtom1->dBornRadius)*(patAtom2->dBornRadius);
return sqrt(dLen_SQR + R1R2*exp(-dLen_SQR/(4.*R1R2)));
}
//---------------------------------------------------------------------------------------------------
TElectr elMakeElectrPair(TAtom& _atAtomA, TAtom& _atAtomB)
{
/**************************************************************************************************************
* It would be convenient to use some simple form of electrostatic energy, like const/r, so we define dQQ
* as the product of two interacting charges times Coulomb constant.
**************************************************************************************************************/
return TElectr(_atAtomA, _atAtomB, TConsts::COULOMB*_atAtomA.dCharge*_atAtomB.dCharge);
}
//---------------------------------------------------------------------------------------------------