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 1502 by gezelter, Sat Oct 2 19:53:32 2010 UTC vs.
Revision 1571 by gezelter, Fri May 27 16:45:44 2011 UTC

# Line 302 | Line 302 | namespace OpenMD {
302      return;
303    }
304  
305 <  void SC::calcDensity(DensityData ddat) {
305 >  void SC::calcDensity(InteractionData &idat) {
306      
307      if (!initialized_) initialize();
308      
309 <    SCInteractionData mixer = MixingMap[make_pair(ddat.atype1, ddat.atype2)];
309 >    SCInteractionData mixer = MixingMap[ idat.atypes ];
310  
311      RealType rcij = mixer.rCut;
312  
313 <    if (ddat.rij < rcij) {
314 <      ddat.rho_i_at_j = mixer.phi->getValueAt(ddat.rij);
315 <      ddat.rho_j_at_i = ddat.rho_i_at_j;
316 <    } else {
317 <      ddat.rho_i_at_j = 0.0;
318 <      ddat.rho_j_at_i = 0.0;
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      }
320 <
320 >    
321      return;
322    }
323  
324 <  void SC::calcFunctional(FunctionalData fdat) {
324 >  void SC::calcFunctional(SelfData &sdat) {
325  
326      if (!initialized_) initialize();
327  
328 <    SCAtomData data1 = SCMap[fdat.atype];
328 >    SCAtomData data1 = SCMap[sdat.atype];
329      
330 <    fdat.frho = - data1.c * data1.epsilon * sqrt(fdat.rho);
331 <    fdat.dfrhodrho = 0.5 * fdat.frho / fdat.rho;
330 >    *(sdat.frho) = - data1.c * data1.epsilon * sqrt( *(sdat.rho) );
331 >    *(sdat.dfrhodrho) = 0.5 * *(sdat.frho) / *(sdat.rho);
332      
333      return;
334    }
335    
336  
337 <  void SC::calcForce(InteractionData idat) {
337 >  void SC::calcForce(InteractionData &idat) {
338      
339      if (!initialized_) initialize();
340      
341 <    SCAtomData data1 = SCMap[idat.atype1];
342 <    SCAtomData data2 = SCMap[idat.atype2];
341 >    SCAtomData data1 = SCMap[idat.atypes.first];
342 >    SCAtomData data2 = SCMap[idat.atypes.second];
343  
344 <    SCInteractionData mixer = MixingMap[make_pair(idat.atype1, idat.atype2)];
344 >    SCInteractionData mixer = MixingMap[idat.atypes];
345  
346      RealType rcij = mixer.rCut;
347  
348 <    if (idat.rij < rcij) {
348 >    if ( *(idat.rij)  < rcij) {
349        RealType vcij = mixer.vCut;
350        
351        pair<RealType, RealType> res;
352        
353 <      res = mixer.phi->getValueAndDerivativeAt(idat.rij);
353 >      res = mixer.phi->getValueAndDerivativeAt( *(idat.rij) );
354        RealType rhtmp = res.first;
355        RealType drhodr = res.second;
356        
357 <      res = mixer.V->getValueAndDerivativeAt(idat.rij);
357 >      res = mixer.V->getValueAndDerivativeAt( *(idat.rij) );
358        RealType vptmp = res.first;
359        RealType dvpdr = res.second;
360        
361        RealType pot_temp = vptmp - vcij;
362 <      idat.vpair += pot_temp;
362 >      *(idat.vpair) += pot_temp;
363        
364 <      RealType dudr = drhodr * (idat.dfrho1 + idat.dfrho2) + dvpdr;
364 >      RealType dudr = drhodr * ( *(idat.dfrho1) + *(idat.dfrho2) ) + dvpdr;
365        
366 <      idat.f1 += idat.d * dudr / idat.rij;
366 >      *(idat.f1) += *(idat.d) * dudr / *(idat.rij) ;
367          
368        // particle_pot is the difference between the full potential
369        // and the full potential without the presence of a particular
# Line 377 | Line 377 | namespace OpenMD {
377        // Most of the particle_pot heavy lifting comes from the
378        // pair interaction, and will be handled by vpair.
379        
380 <      idat.fshift1 = - data1.c * data1.epsilon * sqrt(idat.rho1 - rhtmp);
381 <      idat.fshift2 = - data2.c * data2.epsilon * sqrt(idat.rho2 - rhtmp);
380 >      *(idat.fshift1) = - data1.c * data1.epsilon * sqrt( *(idat.rho1) - rhtmp);
381 >      *(idat.fshift2) = - data2.c * data2.epsilon * sqrt( *(idat.rho2) - rhtmp);
382        
383 <      idat.pot += pot_temp;
383 >      idat.pot[METALLIC_FAMILY] += pot_temp;
384      }
385        
386      return;    
387    }
388 +
389 +  RealType SC::getSuggestedCutoffRadius(pair<AtomType*, AtomType*> atypes) {
390 +    if (!initialized_) initialize();  
391 +
392 +    map<pair<AtomType*, AtomType*>, SCInteractionData>::iterator it;
393 +    it = MixingMap.find(atypes);
394 +    if (it == MixingMap.end())
395 +      return 0.0;
396 +    else  {
397 +      SCInteractionData mixer = (*it).second;
398 +      return mixer.rCut;
399 +    }
400 +  }
401   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines