pme_util.h

// -*-c++-*-
#ifndef STORMM_PME_UTIL_H
#define STORMM_PME_UTIL_H
 
#include <cmath>
#include "copyright.h"
#include "Accelerator/hybrid.h"
#include "Constants/symbol_values.h"
#include "DataTypes/common_types.h"
#include "Math/math_enumerators.h"
#include "Math/radial_derivatives.h"
 
namespace stormm {
namespace energy {
 
using card::Hybrid;
using data_types::isSignedIntegralScalarType;
using stmath::FunctionLevel;
using stmath::radialFirstDerivative;
using stmath::radialSecondDerivative;
using stmath::radialThirdDerivative;
using stmath::radialPartialDerivative;
  
constexpr int default_charge_mapping_order = 5;
constexpr double default_dsum_tol = 1.0e-5;
constexpr double max_dsum_tol = 1.0e-4;
constexpr double default_pme_cutoff = 8.0;
constexpr double default_pme_grid_spacing_target = 1.0;

constexpr double minimum_ewald_coefficient = 0.0625;

constexpr double maximum_ewald_coefficient = 1.0;

double ewaldCoefficient(double cutoff, double direct_sum_tol);

double pmeGaussianSpread(double cutoff, double direct_sum_tol);

double recoverDirectSumTolerance(double cutoff, double ewald_coefficient);

template <typename Tcalc>
Tcalc elecPMEDirectSpace(double ew_coeff, Tcalc kcoul, Tcalc r, Tcalc r2, int order);
 
template <typename Tcalc, typename T3>
Tcalc elecPMEDirectSpace(const T3 pa, const T3 pb, double ew_coeff, Tcalc kcoul,
                         FunctionLevel order = FunctionLevel::VALUE);
 
template <typename Tcalc, typename T3>
std::vector<Tcalc> elecPMEDirectSpace(const T3 pa, const T3 pb, double ew_coeff, Tcalc kcoul,
                                      const std::vector<FunctionLevel> &orders);

double pmeGammaSum(int m, int mesh_length, int ordr);

std::vector<double> pmeLoadBPrefactor(int ordr, int mesh_length);

std::vector<double> pmeLoadMVec(int mesh_length);

std::vector<double> pmeLoadMVecShift(int mesh_length);
  
} // namespace energy
} // namespace stormm
 
#include "pme_util.tpp"
 
#endif