libhoomd/computes/BondTablePotential.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: phillicl
  51
  52 #include <boost/shared_ptr.hpp>
  53
  54 #include "ForceCompute.h"
  55 #include "Index1D.h"
  56 #include "GPUArray.h"
  57
  58 /*! \file BondTablePotential.h
  59     \brief Declares the BondTablePotential class
  60 */
  61
  62 #ifdef NVCC
  63 #error This header cannot be compiled by nvcc
  64 #endif
  65
  66 #ifndef __BONDTABLEPOTENTIAL_H__
  67 #define __BONDTABLEPOTENTIAL_H__
  68
  69 //! Computes the potential and force on bonds based on values given in a table
  70 /*! \b Overview
  71     Bond potentials and forces are evaluated for all bonded particle pairs in the system.
  72     Both the potentials and forces are provided the tables V(r) and F(r) at discreet \a r values between \a rmin and
  73     \a rmax. Evaluations are performed by simple linear interpolation, thus why F(r) must be explicitly specified to
  74     avoid large errors resulting from the numerical derivative. Note that F(r) should store - dV/dr.
  75
  76     \b Table memory layout
  77
  78     V(r) and F(r) are specified for each bond type.
  79
  80     Three parameters need to be stored for each bond potential: rmin, rmax, and dr, the minimum r, maximum r, and spacing
  81     between r values in the table respectively. For simple access on the GPU, these will be stored in a float4 where
  82     x is rmin, y is rmax, and z is dr.
  83
  84     V(0) is the value of V at r=rmin. V(i) is the value of V at r=rmin + dr * i where i is chosen such that r >= rmin
  85     and r <= rmax. V(r) for r < rmin and > rmax is 0. The same goes for F. Thus V and F are defined between the region
  86     [rmin,rmax], inclusive.
  87
  88     For ease of storing the data, all tables must be of the same number of points for all bonds.
  89
  90     \b Interpolation
  91     Values are interpolated linearly between two points straddling the given r. For a given r, the first point needed, i
  92     can be calculated via i = floorf((r - rmin) / dr). The fraction between ri and ri+1 can be calculated via
  93     f = (r - rmin) / dr - float(i). And the linear interpolation can then be performed via V(r) ~= Vi + f * (Vi+1 - Vi)
  94     \ingroup computes
  95 */
  96 class BondTablePotential : public ForceCompute
  97     {
  98     public:
  99         //! Constructs the compute
 100         BondTablePotential(boost::shared_ptr<SystemDefinition> sysdef,
 101                        unsigned int table_width,
 102                        const std::string& log_suffix="");
 103
 104         //! Destructor
 105         virtual ~BondTablePotential();
 106
 107         //! Set the table for a given type pair
 108         virtual void setTable(unsigned int type,
 109                               const std::vector<Scalar> &V,
 110                               const std::vector<Scalar> &F,
 111                               Scalar rmin,
 112                               Scalar rmax);
 113
 114         //! Returns a list of log quantities this compute calculates
 115         virtual std::vector< std::string > getProvidedLogQuantities();
 116
 117         //! Calculates the requested log value and returns it
 118         virtual Scalar getLogValue(const std::string& quantity, unsigned int timestep);
 119
 120         #ifdef ENABLE_MPI
 121         //! Get ghost particle fields requested by this pair potential
 122         /*! \param timestep Current time step
 123         */
 124         virtual CommFlags getRequestedCommFlags(unsigned int timestep)
 125             {
 126                 CommFlags flags = CommFlags(0);
 127                 flags[comm_flag::tag] = 1;
 128                 flags |= ForceCompute::getRequestedCommFlags(timestep);
 129                 return flags;
 130             }
 131         #endif
 132
 133     protected:
 134         boost::shared_ptr<BondData> m_bond_data;    //!< Bond data to use in computing bonds
 135         unsigned int m_table_width;                 //!< Width of the tables in memory
 136         GPUArray<Scalar2> m_tables;                  //!< Stored V and F tables
 137         GPUArray<Scalar4> m_params;                 //!< Parameters stored for each table
 138         Index2D m_table_value;                      //!< Index table helper
 139         std::string m_log_name;                     //!< Cached log name
 140
 141         //! Actually compute the forces
 142         virtual void computeForces(unsigned int timestep);
 143     };
 144
 145 //! Exports the TablePotential class to python
 146 void export_BondTablePotential();
 147
 148 #endif