/*
 * ManyBodyPotential_Tersoff.hpp
 *
 *  Created on: Sep 26, 2012
 *      Author: heber
 */

#ifndef MANYBODYPOTENTIAL_TERSOFF_HPP_
#define MANYBODYPOTENTIAL_TERSOFF_HPP_

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <boost/function.hpp>
#include <cmath>

#include "Potentials/EmpiricalPotential.hpp"
#include "FunctionApproximation/FunctionModel.hpp"

/** This class is the implementation of the Tersoff potential function.
 *
 * \note The arguments_t argument list is here in the following order:
 * -# first \f$ r_{ij}$ \f$,
 * -# then all \f$ r_{ik}$ \f$ that are within the cutoff, i.e. \f$ r_{ik}$ < R + D\f$
 *
 */
class ManyBodyPotential_Tersoff : virtual public EmpiricalPotential, virtual public FunctionModel
{
  // some repeated typedefs to avoid ambiguities
  typedef FunctionModel::arguments_t arguments_t;
  typedef FunctionModel::result_t result_t;
  typedef FunctionModel::results_t results_t;
  typedef EmpiricalPotential::derivative_components_t derivative_components_t;
private:
  //!> cutoff_offset offset for cutoff parameter
  const double cutoff_offset;
  //!> cutoff_halfwidth half-width for cutoff parameter
  const double cutoff_halfwidth;
  //!> many body prefactor parameter for attractive part
  const double manybodyparameter_A;
  //!> many body prefactor parameter for attractive part
  const double manybodyparameter_B;
  //!> many body scale parameter for attractive part
  const double manybodyparameter_lambda1;
  //!> many body scale parameter for attractive part
  const double manybodyparameter_lambda2;
  //!> length scale for coordination influence
  const double manybodyparameter_lambda3;
  //!> many body parameter for attractive part
  const double manybodyparameter_alpha;
  //!> many body parameter for repulsive part
  const double manybodyparameter_beta;
  //!> many body type-dependent parameter giving the power
  const double manybodyparameter_n;
  //!> many-body type-dependent parameter in angular dependence
  const double manybodyparameter_c;
  //!> many-body type-dependent parameter in angular dependence
  const double manybodyparameter_d;
  //!> many-body type-dependent parameter in angular dependence
  const double manybodyparameter_h;
public:
  /** Constructor for class ManyBodyPotential_Tersoff.
   *
   * @param _triplefunction function that returns a list of triples (i.e. the
   *        two remaining distances) to a given pair of points (contained as
   *        indices within the argument)
   */
  ManyBodyPotential_Tersoff(
      const double _cutoff_offset,
      const double _cutoff_halfwidth,
      const double A,
      const double B,
      const double lambda1,
      const double lambda2,
      const double lambda3,
      const double alpha,
      const double beta,
      const double n,
      const double c,
      const double d,
      const double h,
      boost::function< std::vector<arguments_t>(const argument_t &, const double)> &_triplefunction) :
    cutoff_offset(_cutoff_offset),
    cutoff_halfwidth(_cutoff_offset),
    manybodyparameter_A(A),
    manybodyparameter_B(B),
    manybodyparameter_lambda1(lambda1),
    manybodyparameter_lambda2(lambda2),
    manybodyparameter_lambda3(lambda3),
    manybodyparameter_alpha(alpha),
    manybodyparameter_beta(beta),
    manybodyparameter_n(n),
    manybodyparameter_c(d),
    manybodyparameter_d(d),
    manybodyparameter_h(h),
    triplefunction(_triplefunction)
  {}
  /** Destructor of class ManyBodyPotential_Tersoff.
   *
   */
  virtual ~ManyBodyPotential_Tersoff() {}

  /** Evaluates the Tersoff potential for the given arguments.
   *
   * @param arguments single distance
   * @return value of the potential function
   */
  results_t operator()(const arguments_t &arguments) const;

  /** Evaluates the derivative of the Tersoff potential with respect to the
   * input variables.
   *
   * @param arguments single distance
   * @return vector with components of the derivative
   */
  derivative_components_t derivative(const arguments_t &arguments) const;

  /** Evaluates the derivative of the function with the given \a arguments
   * with respect to a specific parameter indicated by \a index.
   *
   * \param arguments set of arguments as input variables to the function
   * \param index derivative of which parameter
   * \return result vector containing the derivative with respect to the given
   *         input
   */
  virtual results_t paramter_derivative(const arguments_t &arguments, const size_t index) const=0;

private:
  /** This function represents the cutoff \f$ f_C \f$.
   *
   * @param distance variable of the function
   * @return a value in [0,1].
   */
  result_t function_cutoff(
      const double &distance
      ) const;
  /** This  function has the exponential feature from the Morse potential.
   *
   * @param distance variable of the function
   * @param prefactor prefactor parameter to exp function
   * @param lambda scale parameter of exp function's argument
   * @return
   */
  result_t function_smoother(
      const double &distance,
      const double &prefactor,
      const double &lambda
      ) const
  {
    return prefactor * exp(-lambda * distance);
  }

  /** This function represents \f$ (1 + \alpha^n \eta^n)^{-1/2n} \f$.
   *
   * @param alpha prefactor to eta function
   * @param r_ij distance argument
   * @param etafunction eta or zeta
   * @return \f$ (1 + \alpha^n \eta^n)^{-1/2n} \f$
   */
  result_t function_prefactor(
      const double &alpha,
      boost::function<result_t()> etafunction
        ) const;

  result_t
  function_eta(
      const argument_t &r_ij
    ) const;

  result_t
  function_zeta(
      const argument_t &r_ij
    ) const;

  result_t
  function_angle(
      const double &r_ij,
      const double &r_ik,
      const double &r_jk
    ) const;

private:
  //!> bound function that obtains the triples for the internal coordinationb summation.
  const boost::function< std::vector< arguments_t >(const argument_t &, const double)> &triplefunction;
};


#endif /* MANYBODYPOTENTIAL_TERSOFF_HPP_ */