| 766 |
|
// Excluded potential that is still computed for fluctuating charges |
| 767 |
|
excluded_Pot= 0.0; |
| 768 |
|
|
| 769 |
– |
|
| 769 |
|
// some variables we'll need independent of electrostatic type: |
| 770 |
|
|
| 771 |
|
ri = 1.0 / *(idat.rij); |
| 886 |
|
Ea += pre14_ * (trQb * rhat * dv21 + 2.0 * Qbr * v22or |
| 887 |
|
+ rdQbr * rhat * (dv22 - 2.0*v22or)); |
| 888 |
|
} |
| 889 |
< |
|
| 890 |
< |
if ((a_is_Fluctuating || b_is_Fluctuating) |
| 891 |
< |
&& (idat.excluded || idat.sameRegion)) { |
| 889 |
> |
|
| 890 |
> |
|
| 891 |
> |
if ((a_is_Fluctuating || b_is_Fluctuating) && idat.excluded) { |
| 892 |
|
J = Jij[FQtids[idat.atid1]][FQtids[idat.atid2]]; |
| 893 |
|
} |
| 894 |
< |
|
| 894 |
> |
|
| 895 |
|
if (a_is_Charge) { |
| 896 |
|
|
| 897 |
|
if (b_is_Charge) { |
| 898 |
|
pref = pre11_ * *(idat.electroMult); |
| 899 |
|
U += C_a * C_b * pref * v01; |
| 900 |
|
F += C_a * C_b * pref * dv01 * rhat; |
| 901 |
< |
|
| 901 |
> |
|
| 902 |
|
// If this is an excluded pair, there are still indirect |
| 903 |
|
// interactions via the reaction field we must worry about: |
| 904 |
|
|
| 907 |
|
indirect_Pot += rfContrib; |
| 908 |
|
indirect_F += rfContrib * 2.0 * ri * rhat; |
| 909 |
|
} |
| 910 |
< |
|
| 910 |
> |
|
| 911 |
|
// Fluctuating charge forces are handled via Coulomb integrals |
| 912 |
< |
// for excluded pairs (i.e. those connected via bonds), atoms |
| 913 |
< |
// within the same region (for metals) and with the standard |
| 915 |
< |
// charge-charge interaction otherwise. |
| 912 |
> |
// for excluded pairs (i.e. those connected via bonds) and |
| 913 |
> |
// with the standard charge-charge interaction otherwise. |
| 914 |
|
|
| 915 |
< |
if (idat.excluded || idat.sameRegion) { |
| 915 |
> |
if (idat.excluded) { |
| 916 |
|
if (a_is_Fluctuating || b_is_Fluctuating) { |
| 917 |
|
coulInt = J->getValueAt( *(idat.rij) ); |
| 918 |
< |
if (a_is_Fluctuating) dUdCa += coulInt * C_b; |
| 919 |
< |
if (b_is_Fluctuating) dUdCb += coulInt * C_a; |
| 920 |
< |
if (idat.excluded) |
| 923 |
< |
excluded_Pot += C_a * C_b * coulInt; |
| 924 |
< |
} |
| 918 |
> |
if (a_is_Fluctuating) dUdCa += C_b * coulInt; |
| 919 |
> |
if (b_is_Fluctuating) dUdCb += C_a * coulInt; |
| 920 |
> |
} |
| 921 |
|
} else { |
| 922 |
|
if (a_is_Fluctuating) dUdCa += C_b * pref * v01; |
| 923 |
|
if (b_is_Fluctuating) dUdCb += C_a * pref * v01; |
| 924 |
< |
} |
| 924 |
> |
} |
| 925 |
|
} |
| 926 |
|
|
| 927 |
|
if (b_is_Dipole) { |
| 987 |
|
F -= pref * (rdDa * rdDb) * (dv22 - 2.0*v22or) * rhat; |
| 988 |
|
Ta += pref * ( v21 * DaxDb - v22 * rdDb * rxDa); |
| 989 |
|
Tb += pref * (-v21 * DaxDb - v22 * rdDa * rxDb); |
| 994 |
– |
|
| 990 |
|
// Even if we excluded this pair from direct interactions, we |
| 991 |
|
// still have the reaction-field-mediated dipole-dipole |
| 992 |
|
// interaction: |
| 1046 |
|
trQaQb = QaQb.trace(); |
| 1047 |
|
rQaQb = rhat * QaQb; |
| 1048 |
|
QaQbr = QaQb * rhat; |
| 1049 |
< |
QaxQb = cross(Q_a, Q_b); |
| 1049 |
> |
QaxQb = mCross(Q_a, Q_b); |
| 1050 |
|
rQaQbr = dot(rQa, Qbr); |
| 1051 |
|
rQaxQbr = cross(rQa, Qbr); |
| 1052 |
|
|
| 1077 |
|
// + 4.0 * cross(rhat, QbQar) |
| 1078 |
|
|
| 1079 |
|
Tb += pref * 2.0 * cross(rhat,Qbr) * rdQar * v43; |
| 1085 |
– |
|
| 1080 |
|
} |
| 1081 |
|
} |
| 1082 |
|
|
| 1437 |
|
dks[i] = dk * skr[i]; |
| 1438 |
|
} |
| 1439 |
|
if (data.is_Quadrupole) { |
| 1440 |
< |
Q = atom->getQuadrupole(); |
| 1441 |
< |
Q *= mPoleConverter; |
| 1448 |
< |
Qk = Q * kVec; |
| 1440 |
> |
Q = atom->getQuadrupole() * mPoleConverter; |
| 1441 |
> |
Qk = Q * kVec; |
| 1442 |
|
qk = dot(kVec, Qk); |
| 1443 |
|
qxk[i] = cross(kVec, Qk); |
| 1444 |
|
qkc[i] = qk * ckr[i]; |
| 1499 |
|
RealType qtrq1 = AK[kk]*(skr[i]*(ckcs-dkss-qkcs) |
| 1500 |
|
-ckr[i]*(ckss+dkcs-qkss)); |
| 1501 |
|
RealType qtrq2 = 2.0*AK[kk]*(ckr[i]*(ckcs-dkss-qkcs) |
| 1502 |
< |
+skr[i]*(ckss+dkcs-qkss)); |
| 1502 |
> |
+skr[i]*(ckss+dkcs-qkss)); |
| 1503 |
|
|
| 1504 |
|
atom->addFrc( 4.0 * rvol * qfrc * kVec ); |
| 1505 |
|
|
