libhoomd/potentials/EvaluatorPairYukawa.h

   1 /*
   2 Highly Optimized Object-oriented Many-particle Dynamics -- Blue Edition
   3 (HOOMD-blue) Open Source Software License Copyright 2009-2014 The Regents of
   4 the University of Michigan All rights reserved.
   5
   6 HOOMD-blue may contain modifications ("Contributions") provided, and to which
   7 copyright is held, by various Contributors who have granted The Regents of the
   8 University of Michigan the right to modify and/or distribute such Contributions.
   9
  10 You may redistribute, use, and create derivate works of HOOMD-blue, in source
  11 and binary forms, provided you abide by the following conditions:
  12
  13 * Redistributions of source code must retain the above copyright notice, this
  14 list of conditions, and the following disclaimer both in the code and
  15 prominently in any materials provided with the distribution.
  16
  17 * Redistributions in binary form must reproduce the above copyright notice, this
  18 list of conditions, and the following disclaimer in the documentation and/or
  19 other materials provided with the distribution.
  20
  21 * All publications and presentations based on HOOMD-blue, including any reports
  22 or published results obtained, in whole or in part, with HOOMD-blue, will
  23 acknowledge its use according to the terms posted at the time of submission on:
  24 http://codeblue.umich.edu/hoomd-blue/citations.html
  25
  26 * Any electronic documents citing HOOMD-Blue will link to the HOOMD-Blue website:
  27 http://codeblue.umich.edu/hoomd-blue/
  28
  29 * Apart from the above required attributions, neither the name of the copyright
  30 holder nor the names of HOOMD-blue's contributors may be used to endorse or
  31 promote products derived from this software without specific prior written
  32 permission.
  33
  34 Disclaimer
  35
  36 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDER AND CONTRIBUTORS ``AS IS'' AND
  37 ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
  38 WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, AND/OR ANY
  39 WARRANTIES THAT THIS SOFTWARE IS FREE OF INFRINGEMENT ARE DISCLAIMED.
  40
  41 IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
  42 INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
  43 BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
  44 DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
  45 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
  46 OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
  47 ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  48 */
  49
  50 // Maintainer: joaander
  51
  52 #ifndef __PAIR_EVALUATOR_YUKAWA_H__
  53 #define __PAIR_EVALUATOR_YUKAWA_H__
  54
  55 #ifndef NVCC
  56 #include <string>
  57 #endif
  58
  59 #include "HOOMDMath.h"
  60
  61 /*! \file EvaluatorPairYukawa.h
  62     \brief Defines the pair evaluator class for Yukawa potentials
  63 */
  64
  65 // need to declare these class methods with __device__ qualifiers when building in nvcc
  66 // DEVICE is __host__ __device__ when included in nvcc and blank when included into the host compiler
  67 #ifdef NVCC
  68 #define DEVICE __device__
  69 #else
  70 #define DEVICE
  71 #endif
  72
  73 //! Class for evaluating the Yukawa pair potential
  74 /*! <b>General Overview</b>
  75
  76     See EvaluatorPairLJ
  77
  78     <b>Yukawa specifics</b>
  79
  80     EvaluatorPairYukawa evaluates the function:
  81     \f[ V_{\mathrm{yukawa}}(r) = \varepsilon \frac{ \exp \left( -\kappa r \right) }{r} \f]
  82
  83     The Yukawa potential does not need diameter or charge. Two parameters are specified and stored in a Scalar2.
  84     \a epsilon is placed in \a params.x and \a kappa is in \a params.y.
  85
  86     These are related to the standard lj parameters sigma and epsilon by:
  87     - \a epsilon = \f$ \varepsilon \f$
  88     - \a kappa = \f$ \kappa \f$
  89
  90 */
  91 class EvaluatorPairYukawa
  92     {
  93     public:
  94         //! Define the parameter type used by this pair potential evaluator
  95         typedef Scalar2 param_type;
  96
  97         //! Constructs the pair potential evaluator
  98         /*! \param _rsq Squared distance beteen the particles
  99             \param _rcutsq Sqauared distance at which the potential goes to 0
 100             \param _params Per type pair parameters of this potential
 101         */
 102         DEVICE EvaluatorPairYukawa(Scalar _rsq, Scalar _rcutsq, const param_type& _params)
 103             : rsq(_rsq), rcutsq(_rcutsq), epsilon(_params.x), kappa(_params.y)
 104             {
 105             }
 106
 107         //! Yukawa doesn't use diameter
 108         DEVICE static bool needsDiameter() { return false; }
 109         //! Accept the optional diameter values
 110         /*! \param di Diameter of particle i
 111             \param dj Diameter of particle j
 112         */
 113         DEVICE void setDiameter(Scalar di, Scalar dj) { }
 114
 115         //! Yukawa doesn't use charge
 116         DEVICE static bool needsCharge() { return false; }
 117         //! Accept the optional diameter values
 118         /*! \param qi Charge of particle i
 119             \param qj Charge of particle j
 120         */
 121         DEVICE void setCharge(Scalar qi, Scalar qj) { }
 122
 123         //! Evaluate the force and energy
 124         /*! \param force_divr Output parameter to write the computed force divided by r.
 125             \param pair_eng Output parameter to write the computed pair energy
 126             \param energy_shift If true, the potential must be shifted so that V(r) is continuous at the cutoff
 127             \note There is no need to check if rsq < rcutsq in this method. Cutoff tests are performed
 128                   in PotentialPair.
 129
 130             \return True if they are evaluated or false if they are not because we are beyond the cuttoff
 131         */
 132         DEVICE bool evalForceAndEnergy(Scalar& force_divr, Scalar& pair_eng, bool energy_shift)
 133             {
 134             // compute the force divided by r in force_divr
 135             if (rsq < rcutsq && epsilon != 0)
 136                 {
 137                 Scalar rinv = fast::rsqrt(rsq);
 138                 Scalar r = Scalar(1.0) / rinv;
 139                 Scalar r2inv = Scalar(1.0) / rsq;
 140
 141                 Scalar exp_val = fast::exp(-kappa * r);
 142
 143                 force_divr = epsilon * exp_val * r2inv * (rinv + kappa);
 144                 pair_eng = epsilon * exp_val * rinv;
 145
 146                 if (energy_shift)
 147                     {
 148                     Scalar rcutinv = fast::rsqrt(rcutsq);
 149                     Scalar rcut = Scalar(1.0) / rcutinv;
 150                     pair_eng -= epsilon * fast::exp(-kappa * rcut) * rcutinv;
 151                     }
 152                 return true;
 153                 }
 154             else
 155                 return false;
 156             }
 157
 158         #ifndef NVCC
 159         //! Get the name of this potential
 160         /*! \returns The potential name. Must be short and all lowercase, as this is the name energies will be logged as
 161             via analyze.log.
 162         */
 163         static std::string getName()
 164             {
 165             return std::string("yukawa");
 166             }
 167         #endif
 168
 169     protected:
 170         Scalar rsq;     //!< Stored rsq from the constructor
 171         Scalar rcutsq;  //!< Stored rcutsq from the constructor
 172         Scalar epsilon; //!< epsilon parameter extracted from the params passed to the constructor
 173         Scalar kappa;   //!< kappa parameter extracted from the params passed to the constructor
 174     };
 175
 176
 177 #endif // __PAIR_EVALUATOR_YUKAWA_H__