ComputeNonbondedFEP.C File Reference

#include "ComputeNonbondedInl.h"
#include "ComputeNonbondedAlch.h"
#include "ComputeNonbondedBase.h"

Go to the source code of this file.

Defines

#define FEPFLAG
#define CALCENERGY
#define NBTYPE   NBPAIR
#define FULLELECT
#define MERGEELECT
#define SLOWONLY
#define NBTYPE   NBSELF
#define FULLELECT
#define MERGEELECT
#define SLOWONLY

Functions

void fep_vdw_forceandenergies (BigReal A, BigReal B, BigReal r2, BigReal myVdwShift, BigReal myVdwShift2, BigReal switchdist2, BigReal cutoff2, BigReal switchfactor, Bool vdwForceSwitching, BigReal myVdwLambda, BigReal myVdwLambda2, Bool alchWCAOn, BigReal myRepLambda, BigReal myRepLambda2, BigReal *alch_vdw_energy, BigReal *alch_vdw_force, BigReal *alch_vdw_energy_2)

Define Documentation

#define CALCENERGY

Definition at line 170 of file ComputeNonbondedFEP.C.

#define FEPFLAG

Definition at line 169 of file ComputeNonbondedFEP.C.

#define FULLELECT

Definition at line 187 of file ComputeNonbondedFEP.C.

#define FULLELECT

Definition at line 187 of file ComputeNonbondedFEP.C.

#define MERGEELECT

Definition at line 189 of file ComputeNonbondedFEP.C.

#define MERGEELECT

Definition at line 189 of file ComputeNonbondedFEP.C.

#define NBTYPE   NBSELF

Definition at line 185 of file ComputeNonbondedFEP.C.

#define NBTYPE   NBPAIR

Definition at line 185 of file ComputeNonbondedFEP.C.

#define SLOWONLY

Definition at line 192 of file ComputeNonbondedFEP.C.

#define SLOWONLY

Definition at line 192 of file ComputeNonbondedFEP.C.


Function Documentation

void fep_vdw_forceandenergies ( BigReal  A,
BigReal  B,
BigReal  r2,
BigReal  myVdwShift,
BigReal  myVdwShift2,
BigReal  switchdist2,
BigReal  cutoff2,
BigReal  switchfactor,
Bool  vdwForceSwitching,
BigReal  myVdwLambda,
BigReal  myVdwLambda2,
Bool  alchWCAOn,
BigReal  myRepLambda,
BigReal  myRepLambda2,
BigReal alch_vdw_energy,
BigReal alch_vdw_force,
BigReal alch_vdw_energy_2 
) [inline]

Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved.

Definition at line 17 of file ComputeNonbondedFEP.C.

References vdw_energy(), vdw_forceandenergy(), vdw_fswitch_energy(), vdw_fswitch_forceandenergy(), vdw_fswitch_shift(), and vdw_switch().

00022                                                          {
00023   BigReal U, F, U_2, dU, dU_2, switchmul, switchmul2;
00024   /*
00025    * The variable naming here is unavoidably awful. A number immediately after
00026    * a variable implies a distance to that power. The suffix "_2" indicates an
00027    * evaluation at alchLambda2. The prefix "inv" means the inverse (previously
00028    * this also used underscores - it was awful).
00029    */
00030   if (alchWCAOn) {
00031     // WCA-on, auxilliary, lambda-dependent cutoff based on Rmin
00032     //
00033     // Avoid divide by zero - corectly zeroes interaction below.
00034     const BigReal Rmin2 = (B <= 0.0 ? 0.0 : powf(2.0*A/B, 1.f/3));
00035     if (myRepLambda < 1.0) {
00036       // modified repulsive soft-core
00037       // All of the scaling is baked into the shift and cutoff values. There
00038       // are no linear coupling terms out front.
00039       const BigReal WCAshift = Rmin2*(1 - myRepLambda)*(1 - myRepLambda);
00040       if (r2 <= Rmin2 - WCAshift) {
00041         const BigReal epsilon = B*B/(4.0*A);
00042         vdw_forceandenergy(A, B, r2 + WCAshift, &U, &F);
00043         *alch_vdw_energy = U + epsilon;
00044         *alch_vdw_force = F;
00045       } else {
00046         *alch_vdw_energy = 0.0;
00047         *alch_vdw_force = 0.0;
00048       }
00049     } else { // myRepLambda == 1.0
00050       if (vdwForceSwitching) {
00051         // force and potential switching
00052         if (r2 <= Rmin2) {
00053           // normal LJ, potential w/two shifts
00054           const BigReal epsilon = B*B/(4.0*A);
00055           vdw_fswitch_shift(A, B, switchdist2, cutoff2, &dU);
00056           vdw_forceandenergy(A, B, r2, &U, &F);
00057           *alch_vdw_energy = U + (1 - myVdwLambda)*epsilon + myVdwLambda*dU;
00058           *alch_vdw_force = F;
00059         } else if (r2 <= switchdist2) {
00060           // normal LJ potential w/shift
00061           vdw_fswitch_shift(A, B, switchdist2, cutoff2, &dU);
00062           vdw_forceandenergy(A, B, r2, &U, &F);
00063           *alch_vdw_energy = myVdwLambda*(U + dU);
00064           *alch_vdw_force = myVdwLambda*F;
00065         } else { // r2 > switchdist
00066           // normal LJ potential with linear coupling   
00067           vdw_fswitch_forceandenergy(A, B, r2, switchdist2, cutoff2, &U, &F);
00068           *alch_vdw_energy = myVdwLambda*U;
00069           *alch_vdw_force = myVdwLambda*F;
00070         }
00071       } else {
00072         // potential switching - also correct for no switching
00073         if (r2 <= Rmin2) {
00074           // normal LJ potential w/shift
00075           const BigReal epsilon = B*B/(4.0*A);
00076           vdw_forceandenergy(A, B, r2, &U, &F);
00077           *alch_vdw_energy = U + (1 - myVdwLambda)*epsilon;
00078           *alch_vdw_force = F; 
00079         } else { // r2 > Rmin2
00080           // normal LJ potential with linear coupling
00081           vdw_switch(r2, switchdist2, cutoff2, switchfactor, &switchmul, \
00082               &switchmul2);
00083           vdw_forceandenergy(A, B, r2, &U, &F);
00084           *alch_vdw_energy = myVdwLambda*switchmul*U;
00085           *alch_vdw_force = myVdwLambda*(switchmul*F + switchmul2*U);
00086         }
00087       }
00088     }
00089     if (myRepLambda2 < 1.0) {
00090       // Same as above, but energy only. This is done separately bc the
00091       // cutoffs are lambda dependent.
00092       const BigReal WCAshift_2 = Rmin2*(1 - myRepLambda2)*(1 - myRepLambda2);
00093       if (r2 <= Rmin2 - WCAshift_2) {
00094         const BigReal epsilon = B*B/(4.0*A);
00095         vdw_energy(A, B, r2 + WCAshift_2, &U_2);
00096         *alch_vdw_energy_2 = U_2 + epsilon;
00097       } else {
00098         *alch_vdw_energy_2 = 0.0;
00099       }
00100     } else { // myRepLambda2 == 1.0
00101       if (vdwForceSwitching) {
00102         // force and potential switching
00103         if (r2 <= Rmin2) {
00104           // normal LJ, potential w/two shifts
00105           const BigReal epsilon = B*B/(4.0*A);
00106           vdw_fswitch_shift(A, B, switchdist2, cutoff2, &dU_2);
00107           vdw_energy(A, B, r2, &U_2);
00108           *alch_vdw_energy_2 = \
00109               U_2 + (1 - myVdwLambda2)*epsilon + myVdwLambda2*dU_2;
00110         } else if (r2 <= switchdist2) {
00111           // normal LJ potential w/shift
00112           vdw_fswitch_shift(A, B, switchdist2, cutoff2, &dU_2);
00113           vdw_energy(A, B, r2, &U_2);
00114           *alch_vdw_energy_2 = myVdwLambda2*(U_2 + dU_2);
00115         } else { // r2 > switchdist
00116           vdw_fswitch_energy(A, B, r2, switchdist2, cutoff2, &U_2);
00117           *alch_vdw_energy_2 = myVdwLambda2*U_2;
00118         }
00119       } else {
00120         // potential switching - also correct for no switching
00121         if (r2 <= Rmin2) {
00122           // normal LJ potential w/shift
00123           const BigReal epsilon = B*B/(4.0*A);
00124           vdw_energy(A, B, r2, &U_2);
00125           *alch_vdw_energy_2 = U_2 + (1 - myVdwLambda2)*epsilon;
00126         } else { // r2 > Rmin2
00127           // normal LJ potential with linear coupling
00128           vdw_switch(r2, switchdist2, cutoff2, switchfactor, &switchmul, \
00129               &switchmul2);
00130           vdw_energy(A, B, r2, &U_2);
00131           *alch_vdw_energy_2 = myVdwLambda2*switchmul*U_2;
00132         }
00133       }
00134     }
00135   } else { // WCA-off
00136     if (vdwForceSwitching) {
00137       // force and potential switching
00138       if (r2 <= switchdist2) {
00139         // normal LJ potential w/shift
00140         vdw_forceandenergy(A, B, r2 + myVdwShift, &U, &F);
00141         vdw_fswitch_shift(A, B, switchdist2 + myVdwShift, cutoff2, &dU);
00142         vdw_energy(A, B, r2 + myVdwShift2, &U_2);
00143         vdw_fswitch_shift(A, B, switchdist2 + myVdwShift2, cutoff2, &dU_2);
00144         *alch_vdw_energy = myVdwLambda*(U + dU);
00145         *alch_vdw_energy_2 = myVdwLambda2*(U_2 + dU_2);
00146         *alch_vdw_force = myVdwLambda*F;
00147       } else { // r2 > switchdist2
00148         vdw_fswitch_forceandenergy(A, B, r2 + myVdwShift, \
00149             switchdist2 + myVdwShift, cutoff2, &U, &F);
00150         vdw_fswitch_energy(A, B, r2 + myVdwShift2, switchdist2 + myVdwShift2, \
00151             cutoff2, &U_2);
00152         *alch_vdw_energy = myVdwLambda*U;
00153         *alch_vdw_energy_2 = myVdwLambda2*U_2;
00154         *alch_vdw_force = myVdwLambda*F;
00155       }
00156     } else {
00157       // potential switching - also correct for no switching
00158       vdw_switch(r2, switchdist2, cutoff2, switchfactor, &switchmul, \
00159           &switchmul2);
00160       vdw_forceandenergy(A, B, r2 + myVdwShift, &U, &F);
00161       vdw_energy(A, B, r2 + myVdwShift2, &U_2); 
00162       *alch_vdw_energy = myVdwLambda*switchmul*U;
00163       *alch_vdw_energy_2 = myVdwLambda2*switchmul*U_2;
00164       *alch_vdw_force = myVdwLambda*(switchmul*F + switchmul2*U);
00165     }
00166   }
00167 }


Generated on 21 Sep 2019 for NAMD by  doxygen 1.6.1