ComputeNonbondedAlch.h File Reference

Go to the source code of this file.

Functions

void vdw_forceandenergy (const BigReal A, const BigReal B, const BigReal r2, BigReal *U, BigReal *F)
void vdw_energy (const BigReal A, const BigReal B, const BigReal r2, BigReal *U)
void vdw_switch (const BigReal r2, const BigReal switchdist2, const BigReal cutoff2, const BigReal switchfactor, BigReal *switchmul, BigReal *switchmul2)
void vdw_fswitch_forceandenergy (const BigReal A, const BigReal B, const BigReal r2, const BigReal switchdist2, const BigReal cutoff2, BigReal *U, BigReal *F)
void vdw_fswitch_energy (const BigReal A, const BigReal B, const BigReal r2, const BigReal switchdist2, const BigReal cutoff2, BigReal *U)
void vdw_fswitch_shift (const BigReal A, const BigReal B, const BigReal switchdist2, const BigReal cutoff2, BigReal *dU)

Function Documentation

void vdw_energy ( const BigReal  A,
const BigReal  B,
const BigReal  r2,
BigReal U 
) [inline]

Definition at line 25 of file ComputeNonbondedAlch.h.

Referenced by fep_vdw_forceandenergies().

00026                                   {
00027   const BigReal inv_r2 = 1. / r2;
00028   const BigReal inv_r6 = inv_r2*inv_r2*inv_r2;
00029   *U = inv_r6*(A*inv_r6 - B);
00030 }

void vdw_forceandenergy ( const BigReal  A,
const BigReal  B,
const BigReal  r2,
BigReal U,
BigReal F 
) [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 16 of file ComputeNonbondedAlch.h.

Referenced by fep_vdw_forceandenergies(), and ti_vdw_force_energy_dUdl().

00017                                               {
00018   const BigReal inv_r2 = 1. / r2;
00019   const BigReal inv_r6 = inv_r2*inv_r2*inv_r2;
00020   *U = inv_r6*(A*inv_r6 - B);
00021   *F = 6*inv_r2*(2*(*U) + B*inv_r6);
00022 }

void vdw_fswitch_energy ( const BigReal  A,
const BigReal  B,
const BigReal  r2,
const BigReal  switchdist2,
const BigReal  cutoff2,
BigReal U 
) [inline]

Definition at line 75 of file ComputeNonbondedAlch.h.

Referenced by fep_vdw_forceandenergies().

00077                 {  
00078   // TODO: Some rearrangement could be done to use rsqrt() instead?
00079   const BigReal inv_r2 = 1. / r2;
00080   const BigReal inv_r6 = inv_r2*inv_r2*inv_r2;
00081   const BigReal inv_r3 = sqrt(inv_r6);
00082   const BigReal inv_cutoff6 = 1. / (cutoff2*cutoff2*cutoff2);
00083   const BigReal inv_cutoff3 = sqrt(inv_cutoff6);
00084   const BigReal switchdist6 = switchdist2*switchdist2*switchdist2;
00085   const BigReal k_vdwa = A / (1 - switchdist6*inv_cutoff6);
00086   const BigReal k_vdwb = B / (1 - sqrt(switchdist6*inv_cutoff6));
00087   const BigReal tmpa = inv_r6 - inv_cutoff6;
00088   const BigReal tmpb = inv_r3 - inv_cutoff3;
00089   *U = k_vdwa*tmpa*tmpa - k_vdwb*tmpb*tmpb;
00090 }

void vdw_fswitch_forceandenergy ( const BigReal  A,
const BigReal  B,
const BigReal  r2,
const BigReal  switchdist2,
const BigReal  cutoff2,
BigReal U,
BigReal F 
) [inline]

Definition at line 56 of file ComputeNonbondedAlch.h.

Referenced by fep_vdw_forceandenergies(), and ti_vdw_force_energy_dUdl().

00058                             {
00059   // TODO: Some rearrangement could be done to use rsqrt() instead?
00060   const BigReal inv_r2 = 1. / r2;
00061   const BigReal inv_r6 = inv_r2*inv_r2*inv_r2;
00062   const BigReal inv_r3 = sqrt(inv_r6);
00063   const BigReal inv_cutoff6 = 1. / (cutoff2*cutoff2*cutoff2);
00064   const BigReal inv_cutoff3 = sqrt(inv_cutoff6);
00065   const BigReal switchdist6 = switchdist2*switchdist2*switchdist2;
00066   const BigReal k_vdwa = A / (1 - switchdist6*inv_cutoff6);
00067   const BigReal k_vdwb = B / (1 - sqrt(switchdist6*inv_cutoff6));
00068   const BigReal tmpa = inv_r6 - inv_cutoff6;
00069   const BigReal tmpb = inv_r3 - inv_cutoff3;
00070   *U = k_vdwa*tmpa*tmpa - k_vdwb*tmpb*tmpb;
00071   *F = 6*inv_r2*(2*k_vdwa*tmpa*inv_r6 - k_vdwb*tmpb*inv_r3);
00072 }

void vdw_fswitch_shift ( const BigReal  A,
const BigReal  B,
const BigReal  switchdist2,
const BigReal  cutoff2,
BigReal dU 
) [inline]

Definition at line 97 of file ComputeNonbondedAlch.h.

Referenced by fep_vdw_forceandenergies(), and ti_vdw_force_energy_dUdl().

00098                                                                    {
00099   // TODO: Some rearrangement could be done to use rsqrt() instead.
00100   const BigReal cutoff6 = cutoff2*cutoff2*cutoff2;
00101   const BigReal switchdist6 = switchdist2*switchdist2*switchdist2;
00102   const BigReal v_vdwa = -A / (cutoff6*switchdist6);
00103   const BigReal v_vdwb = -B / sqrt(cutoff6*switchdist6);
00104   *dU = v_vdwa - v_vdwb; //deltaV2 from Steinbach & Brooks
00105 }

void vdw_switch ( const BigReal  r2,
const BigReal  switchdist2,
const BigReal  cutoff2,
const BigReal  switchfactor,
BigReal switchmul,
BigReal switchmul2 
) [inline]

Definition at line 35 of file ComputeNonbondedAlch.h.

Referenced by fep_vdw_forceandenergies(), and ti_vdw_force_energy_dUdl().

00037                          {
00038   if (r2 <= switchdist2) {
00039     // "switching off" will always fall through here.
00040     *switchmul = 1.0;
00041     *switchmul2 = 0.0;
00042   } else {
00043     const BigReal cdiff2 = (cutoff2 - r2);
00044     const BigReal sdiff2 = (r2 - switchdist2);
00045     *switchmul = switchfactor*cdiff2*cdiff2*(cdiff2 + 3*sdiff2);
00046     *switchmul2 = 12*switchfactor*cdiff2*sdiff2;
00047   }
00048 }


Generated on 13 Dec 2019 for NAMD by  doxygen 1.6.1