00001
00007
00008
00009
00010
00011 {
00012 const BigReal coulomb = COULOMB * ComputeNonbondedUtil::scaling
00013 * ComputeNonbondedUtil::dielectric_1;
00014 BigReal *dp1 = data1;
00015 BigReal *rp1 = results1;
00016 int j_begin = 0;
00017 register BigReal electEnergy = 0.;
00018 register BigReal virial_xx = 0.;
00019 register BigReal virial_xy = 0.;
00020 register BigReal virial_xz = 0.;
00021 register BigReal virial_yy = 0.;
00022 register BigReal virial_yz = 0.;
00023 register BigReal virial_zz = 0.;
00024
00025 #ifdef FULLDIRECT_PERIODIC
00026 Vector a1 = lattice->a();
00027 Vector b1(0); if ( lattice->a_p() ) b1 = lattice->a_r();
00028 Vector a2 = lattice->b();
00029 Vector b2(0); if ( lattice->b_p() ) b2 = lattice->b_r();
00030 Vector a3 = lattice->c();
00031 Vector b3(0); if ( lattice->c_p() ) b3 = lattice->c_r();
00032 #endif
00033
00034 for(int i=0; i<n1; ++i)
00035 {
00036 register BigReal p_i_x = *(dp1++);
00037 register BigReal p_i_y = *(dp1++);
00038 register BigReal p_i_z = *(dp1++);
00039 register BigReal kq_i = coulomb * *(dp1++);
00040 register BigReal f_i_x = 0.;
00041 register BigReal f_i_y = 0.;
00042 register BigReal f_i_z = 0.;
00043 if ( selfmode )
00044 {
00045 ++j_begin; data2 += 4; results2 += 3;
00046 }
00047 register BigReal *dp2 = data2;
00048 register BigReal *rp2 = results2;
00049 register int n2c = n2;
00050 register int j;
00051 for( j = j_begin; j<n2c; ++j)
00052 {
00053 register BigReal p_ij_x = p_i_x - *(dp2++);
00054 register BigReal p_ij_y = p_i_y - *(dp2++);
00055 register BigReal p_ij_z = p_i_z - *(dp2++);
00056
00057 #ifdef FULLDIRECT_PERIODIC
00058 Vector p_ij(p_ij_x,p_ij_y,p_ij_z);
00059 p_ij -= ( a1*floor(0.5+b1*p_ij) + a2*floor(0.5+b2*p_ij) + a3*floor(0.5+b3*p_ij) );
00060 p_ij_x = p_ij.x;
00061 p_ij_y = p_ij.y;
00062 p_ij_z = p_ij.z;
00063 #endif
00064
00065 register BigReal r_1;
00066 r_1 = 1./sqrt(p_ij_x * p_ij_x + p_ij_y * p_ij_y + p_ij_z * p_ij_z);
00067 register BigReal f = *(dp2++) * kq_i * r_1;
00068 electEnergy += f;
00069 f *= r_1*r_1;
00070 virial_xx += f * p_ij_x * p_ij_x;
00071 virial_xy += f * p_ij_x * p_ij_y;
00072 virial_xz += f * p_ij_x * p_ij_z;
00073 p_ij_x *= f;
00074 virial_yy += f * p_ij_y * p_ij_y;
00075 virial_yz += f * p_ij_y * p_ij_z;
00076 p_ij_y *= f;
00077 virial_zz += f * p_ij_z * p_ij_z;
00078 p_ij_z *= f;
00079 f_i_x += p_ij_x;
00080 f_i_y += p_ij_y;
00081 f_i_z += p_ij_z;
00082 *(rp2++) -= p_ij_x;
00083 *(rp2++) -= p_ij_y;
00084 *(rp2++) -= p_ij_z;
00085 }
00086 *(rp1++) += f_i_x;
00087 *(rp1++) += f_i_y;
00088 *(rp1++) += f_i_z;
00089 }
00090
00091 virial.xx += virial_xx;
00092 virial.xy += virial_xy;
00093 virial.xz += virial_xz;
00094 virial.yx += virial_xy;
00095 virial.yy += virial_yy;
00096 virial.yz += virial_yz;
00097 virial.zx += virial_xz;
00098 virial.zy += virial_yz;
00099 virial.zz += virial_zz;
00100 return electEnergy;
00101 }
00102
