ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/nonbonded/Electrostatic.cpp
(Generate patch)

Comparing branches/development/src/nonbonded/Electrostatic.cpp (file contents):
Revision 1502 by gezelter, Sat Oct 2 19:53:32 2010 UTC vs.
Revision 1505 by gezelter, Sun Oct 3 22:18:59 2010 UTC

# Line 283 | Line 283 | namespace OpenMD {
283            simError();                  
284          }
285          
286 +        // Quadrupoles in OpenMD are set as the diagonal elements
287 +        // of the diagonalized traceless quadrupole moment tensor.
288 +        // The column vectors of the unitary matrix that diagonalizes
289 +        // the quadrupole moment tensor become the eFrame (or the
290 +        // electrostatic version of the body-fixed frame.
291 +
292          Vector3dGenericData* v3dData = dynamic_cast<Vector3dGenericData*>(data);
293          if (v3dData == NULL) {
294            sprintf( painCave.errMsg,
# Line 883 | Line 889 | namespace OpenMD {
889  
890        if (i_is_Dipole) {
891          mu_i = data1.dipole_moment;
892 <        uz_i = skdat.eFrame1->getColumn(2);      
892 >        uz_i = skdat.eFrame1.getColumn(2);      
893          ct_i = dot(uz_i, rhat);
894          duduz_i = V3Zero;
895        }
896              
897        if (j_is_Dipole) {
898          mu_j = data2.dipole_moment;
899 <        uz_j = skdat.eFrame2->getColumn(2);      
899 >        uz_j = skdat.eFrame2.getColumn(2);      
900          ct_j = dot(uz_j, rhat);
901          duduz_j = V3Zero;
902        }
# Line 932 | Line 938 | namespace OpenMD {
938        skdat.f1 += dVdr;
939        
940        if (i_is_Dipole)
941 <        *(skdat.t1) -= cross(uz_i, duduz_i);
941 >        skdat.t1 -= cross(uz_i, duduz_i);
942        if (j_is_Dipole)
943 <        *(skdat.t2) -= cross(uz_j, duduz_j);
943 >        skdat.t2 -= cross(uz_j, duduz_j);
944      }
945    }
946      
# Line 957 | Line 963 | namespace OpenMD {
963          scdat.pot -= 0.5 * preVal;
964          
965          // The self-correction term adds into the reaction field vector
966 <        Vector3d uz_i = scdat.eFrame->getColumn(2);
966 >        Vector3d uz_i = scdat.eFrame.getColumn(2);
967          Vector3d ei = preVal * uz_i;
968  
969          // This looks very wrong.  A vector crossed with itself is zero.
970 <        *(scdat.t) -= cross(uz_i, ei);
970 >        scdat.t -= cross(uz_i, ei);
971        }
972      } else if (summationMethod_ == SHIFTED_FORCE || summationMethod_ == SHIFTED_POTENTIAL) {
973        if (i_is_Charge) {        
# Line 975 | Line 981 | namespace OpenMD {
981        }
982      }
983    }
984 +
985 +  RealType Electrostatic::getSuggestedCutoffRadius(AtomType* at1, AtomType* at2) {
986 +    // This seems to work moderately well as a default.  There's no
987 +    // inherent scale for 1/r interactions that we can standardize.
988 +    // 12 angstroms seems to be a reasonably good guess for most
989 +    // cases.
990 +    return 12.0;
991 +  }
992   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines