6#ifndef GEMMI_FORMFACT_HPP_
7#define GEMMI_FORMFACT_HPP_
22 static_assert(std::numeric_limits<float>::is_iec559,
"float is not IEEE 754?");
24 float val = 12102203.1615614f * x + 1065353216.f;
26 std::int32_t
vali =
static_cast<std::int32_t
>(val);
27 std::int32_t
xu1 =
vali & 0x7F800000;
28 std::int32_t
xu2 = (
vali & 0x7FFFFF) | 0x3F800000;
30 std::memcpy(&a, &
xu1, 4);
31 std::memcpy(&b, &
xu2, 4);
32 return a * (0.509871020f + b * (0.312146713f + b * (0.166617139f + b *
33 (-2.190619930e-3f + b * 1.3555747234e-2f))));
37template<
int N,
typename Real>
43 for (
int i = 0;
i <
N; ++
i)
51 for (
int i = 0;
i <
N; ++
i) {
66 for (
int i = 0;
i <
N; ++
i)
68 for (
int i = 0;
i <
N; ++
i)
77 for (
int i = 0;
i <
N; ++
i)
78 tmp[
i] = std::max(
b[
i] * (
r *
r), -88.f);
79 for (
int i = 0;
i <
N; ++
i) {
89template<
int N,
typename Real>
96 for (
int i = 0;
i <
N; ++
i)
111 for (
int i = 0;
i <
N; ++
i)
112 tmp[
i] = std::max(
b[
i].r_u_r(
r), -88.f);
113 for (
int i = 0;
i <
N; ++
i)
119template<
typename Real>
123template<
int N,
int WithC,
typename Real>
137 for (
int i = 0;
i <
N; ++
i)
138 sf +=
a(
i) * std::exp(-
b(
i)*stol2);
146 for (
int i = 0;
i <
N; ++
i) {
157 for (
int i = 0;
i <
N; ++
i) {
174 std::exp(-4 *
pi2 *
B.inverse().r_u_r(
r));
175 for (
int i = 0;
i <
N; ++
i) {
178 std::exp(-4 *
pi2 *
Bb.inverse().r_u_r(
r));
188 for (
int i = 0;
i <
N; ++
i) {
constexpr float sq(float x)
float unsafe_expapprox(float x)
float calculate(const Vec3 &r_) const
Real calculate(const Vec3 &r) const
std::pair< float, float > calculate_with_derivative(float r) const
float calculate(float r2) const
Real calculate(Real r2) const
std::pair< Real, Real > calculate_with_derivative(Real r) const
void set_coefs(const std::array< Real, 2 *N+WithC > &c)
Real calculate_sf(Real stol2) const
ExpAnisoSum< N+WithC, Real > precalculate_density_aniso_b(const SMat33< Real > &B, Real addend=0) const
ExpSum< N+WithC, Real > precalculate_density_iso(Real B, Real addend=0) const
ExpAnisoSum< N+WithC, Real > precalculate_density_aniso_u(const SMat33< float > &U, Real addend=0) const
Real calculate_density_aniso(const Vec3 &r, const SMat33< float > &U) const
std::array< Real, 2 *N+WithC > coefs
Real calculate_density_iso(Real r2, Real B) const