MultipoleMultipoleTerm.h

Go to the documentation of this file.

#ifndef NDDO_MULTIPOLE_MULTIPOLEMULTIPOLETERM_H
#define NDDO_MULTIPOLE_MULTIPOLEMULTIPOLETERM_H

#include <Utils/Math/AutomaticDifferentiation/AutomaticDifferentiationHelpers.h>

namespace Scine {
namespace Sparrow {

namespace nddo {

namespace multipole {

class MultipoleMultipoleTerm {
 public:
  MultipoleMultipoleTerm(double f, double x1, double y1, double z1, double x2, double y2, double z2)
    : f_(f), x1_(x1), x2_(x2), y1_(y1), y2_(y2), z1_(z1), z2_(z2) {
  }
  template<Utils::DerivativeOrder O>
  Utils::AutomaticDifferentiation::Value1DType<O> calculate(double R, double D1, double d2, double squaredRhos) const {
    double dx = d2 * x2_ - D1 * x1_;
    double dy = d2 * y2_ - D1 * y1_;
    double dz = R + d2 * z2_ - D1 * z1_;
    double invsqrt = 1.0 / std::sqrt(dx * dx + dy * dy + dz * dz + squaredRhos);
    return expr<O>(f_, dz, invsqrt);
  }
  template<Utils::DerivativeOrder O>
  Utils::AutomaticDifferentiation::Value1DType<O> expr(double f, double dz, double invsqrt) const;

 private:
  double f_;
  double x1_, x2_, y1_, y2_, z1_, z2_;
};

template<>
inline Utils::AutomaticDifferentiation::Value1DType<Utils::DerivativeOrder::Zero>
MultipoleMultipoleTerm::expr<Utils::DerivativeOrder::Zero>(double f, double /*dz*/, double invsqrt) const {
  return f * invsqrt;
}
template<>
inline Utils::AutomaticDifferentiation::Value1DType<Utils::DerivativeOrder::One>
MultipoleMultipoleTerm::expr<Utils::DerivativeOrder::One>(double f, double dz, double invsqrt) const {
  return {f * invsqrt, -dz * f * invsqrt * invsqrt * invsqrt};
}
template<>
inline Utils::AutomaticDifferentiation::Value1DType<Utils::DerivativeOrder::Two>
MultipoleMultipoleTerm::expr<Utils::DerivativeOrder::Two>(double f, double dz, double invsqrt) const {
  double inv3 = invsqrt * invsqrt * invsqrt;
  return {f * invsqrt, -dz * f * inv3, f * inv3 * (3 * dz * dz * invsqrt * invsqrt - 1)};
}

} // namespace multipole

} // namespace nddo

} // namespace Sparrow
} // namespace Scine
#endif // NDDO_MULTIPOLE_MULTIPOLEMULTIPOLETERM_H