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

Comparing branches/development/src/nonbonded/SC.cpp (file contents):
Revision 1574 by gezelter, Fri May 27 16:45:44 2011 UTC vs.
Revision 1575 by gezelter, Fri Jun 3 21:39:49 2011 UTC

# Line 311 | Line 311 | namespace OpenMD {
311      RealType rcij = mixer.rCut;
312  
313      if ( *(idat.rij)  < rcij) {
314 <      *(idat.rho_i_at_j) = mixer.phi->getValueAt( *(idat.rij) );
315 <      *(idat.rho_j_at_i) = *(idat.rho_i_at_j);
316 <    } else {
317 <      *(idat.rho_i_at_j) = 0.0;
318 <      *(idat.rho_j_at_i) = 0.0;
319 <    }
314 >      RealType rho = mixer.phi->getValueAt( *(idat.rij) );
315 >      *(idat.rho1) += rho;
316 >      *(idat.rho2) += rho;
317 >    }
318      
319      return;
320    }
# Line 326 | Line 324 | namespace OpenMD {
324      if (!initialized_) initialize();
325  
326      SCAtomData data1 = SCMap[sdat.atype];
327 <    
328 <    *(sdat.frho) = - data1.c * data1.epsilon * sqrt( *(sdat.rho) );
327 >  
328 >    RealType u = - data1.c * data1.epsilon * sqrt( *(sdat.rho) );
329 >    *(sdat.frho) = u;
330      *(sdat.dfrhodrho) = 0.5 * *(sdat.frho) / *(sdat.rho);
331 +
332 +    sdat.pot[METALLIC_FAMILY] += u;
333 +    *(sdat.particlePot) += u;
334      
335      return;
336    }
# Line 365 | Line 367 | namespace OpenMD {
367        
368        *(idat.f1) += *(idat.d) * dudr / *(idat.rij) ;
369          
370 <      // particle_pot is the difference between the full potential
371 <      // and the full potential without the presence of a particular
370 >      // particlePot is the difference between the full potential and
371 >      // the full potential without the presence of a particular
372        // particle (atom1).
373        //
374 <      // This reduces the density at other particle locations, so
375 <      // we need to recompute the density at atom2 assuming atom1
376 <      // didn't contribute.  This then requires recomputing the
377 <      // density functional for atom2 as well.
378 <      //
379 <      // Most of the particle_pot heavy lifting comes from the
380 <      // pair interaction, and will be handled by vpair.
374 >      // This reduces the density at other particle locations, so we
375 >      // need to recompute the density at atom2 assuming atom1 didn't
376 >      // contribute.  This then requires recomputing the density
377 >      // functional for atom2 as well.
378 >          
379 >      *(idat.particlePot1) -= data2.c * data2.epsilon *
380 >        sqrt( *(idat.rho2) - rhtmp) + *(idat.frho2);
381 >
382 >      *(idat.particlePot1) -= data1.c * data1.epsilon *
383 >        sqrt( *(idat.rho1) - rhtmp) + *(idat.frho1);
384        
380      *(idat.fshift1) = - data1.c * data1.epsilon * sqrt( *(idat.rho1) - rhtmp);
381      *(idat.fshift2) = - data2.c * data2.epsilon * sqrt( *(idat.rho2) - rhtmp);
382      
385        idat.pot[METALLIC_FAMILY] += pot_temp;
386      }
387        

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines