--- branches/development/src/nonbonded/LJ.cpp	2010/07/19 18:59:59	1471
+++ branches/development/src/nonbonded/LJ.cpp	2010/07/21 18:31:05	1476
@@ -107,6 +107,7 @@ namespace OpenMD {
   }
 
   RealType LJ::getSigma(int atid) { 
+    if (!initialized_) initialize();
     std::map<int, AtomType*> :: const_iterator it;
     it = LJMap.find(atid);
     if (it == LJMap.end()) {
@@ -142,6 +143,7 @@ namespace OpenMD {
   }
 
   RealType LJ::getEpsilon(int atid) {    
+    if (!initialized_) initialize();
     std::map<int, AtomType*> :: const_iterator it;
     it = LJMap.find(atid);
     if (it == LJMap.end()) {
@@ -164,23 +166,23 @@ namespace OpenMD {
   }
 
   void LJ::initialize() {    
-    ForceField::AtomTypeContainer atomTypes = forceField_->getAtomTypes();
+    ForceField::AtomTypeContainer* atomTypes = forceField_->getAtomTypes();
     ForceField::AtomTypeContainer::MapTypeIterator i;
     AtomType* at;
 
-    for (at = atomTypes.beginType(i); at != NULL; 
-         at = atomTypes.nextType(i)) {
+    for (at = atomTypes->beginType(i); at != NULL; 
+         at = atomTypes->nextType(i)) {
       
       if (at->isLennardJones())
         addType(at);
     }
 
-    ForceField::NonBondedInteractionTypeContainer nbiTypes = forceField_->getNonBondedInteractionTypes();
+    ForceField::NonBondedInteractionTypeContainer* nbiTypes = forceField_->getNonBondedInteractionTypes();
     ForceField::NonBondedInteractionTypeContainer::MapTypeIterator j;
     NonBondedInteractionType* nbt;
 
-    for (nbt = nbiTypes.beginType(j); nbt != NULL; 
-         nbt = nbiTypes.nextType(j)) {
+    for (nbt = nbiTypes->beginType(j); nbt != NULL; 
+         nbt = nbiTypes->nextType(j)) {
       
       if (nbt->isLennardJones()) {
         
@@ -228,7 +230,8 @@ namespace OpenMD {
 
     // add it to the map:
     AtomTypeProperties atp = atomType->getATP();    
-    std::pair<std::map<int,AtomType*>::iterator,bool> ret;
+
+    std::pair<std::map<int,AtomType*>::iterator,bool> ret;    
     ret = LJMap.insert( std::pair<int, AtomType*>(atp.ident, atomType) );
     if (ret.second == false) {
       sprintf( painCave.errMsg,
@@ -285,9 +288,10 @@ namespace OpenMD {
     }    
   }
  
-  void LJ::calcForce(AtomType* at1, AtomType* at2, Vector3d d, RealType rij, 
-                     RealType r2, RealType rcut, RealType sw, RealType vdwMult,
-                     RealType vpair, RealType pot, Vector3d f1) {
+  void LJ::calcForce(AtomType* at1, AtomType* at2, Vector3d d, 
+                     RealType rij, RealType r2, RealType rcut, RealType sw, 
+                     RealType vdwMult, RealType &vpair, RealType &pot, 
+                     Vector3d &f1) {
 
     if (!initialized_) initialize();
     
@@ -304,6 +308,7 @@ namespace OpenMD {
     RealType sigmai = mixer.sigmai;
     RealType epsilon = mixer.epsilon;
     
+
     ros = rij * sigmai;
 
     getLJfunc(ros, myPot, myDeriv);
@@ -330,6 +335,7 @@ namespace OpenMD {
     f1 = d * dudr / rij;
 
     return;
+
   }
 
   void LJ::do_lj_pair(int *atid1, int *atid2, RealType *d, RealType *rij, 
@@ -347,23 +353,30 @@ namespace OpenMD {
     
     calcForce(atype1, atype2, disp, *rij, *r2, *rcut, *sw, *vdwMult, *vpair, 
               *pot, frc);
+      
+    f1[0] = frc.x();
+    f1[1] = frc.y();
+    f1[2] = frc.z();
+
     return;    
   }
   
-  void LJ::getLJfunc(const RealType r, RealType pot, RealType deriv) {
+  void LJ::getLJfunc(RealType r, RealType &pot, RealType &deriv) {
+
     RealType ri = 1.0 / r;
     RealType ri2 = ri * ri;
     RealType ri6 = ri2 * ri2 * ri2;
     RealType ri7 = ri6 * ri;
     RealType ri12 = ri6 * ri6;
     RealType ri13 = ri12 * ri;
-
+    
     pot = 4.0 * (ri12 - ri6);
     deriv = 24.0 * (ri7 - 2.0 * ri13);
+
     return;
   }
+  
 
-
   void LJ::setLJDefaultCutoff(RealType *thisRcut, int *sP, int *sF) {
     shiftedPot_ = (bool)(*sP);
     shiftedFrc_ = (bool)(*sF);