[BE] Replace M_PI with c10::pi constexpr variable (#50819)

Summary:
Also, get rid of MSVC specific `_USE_MATH_DEFINES`

Test at compile time that c10::pi<double> == M_PI

Pull Request resolved: https://github.com/pytorch/pytorch/pull/50819

Reviewed By: albanD

Differential Revision: D25976330

Pulled By: malfet

fbshipit-source-id: 8f3ddfd58a5aa4bd382da64ad6ecc679706d1284
This commit is contained in:
Nikita Shulga 2021-01-22 15:12:51 -08:00 committed by Facebook GitHub Bot
parent a66851a2ad
commit de8cd6b201
21 changed files with 84 additions and 76 deletions

View file

@ -2,6 +2,7 @@
#include <ATen/Utils.h>
#include <ATen/core/MT19937RNGEngine.h>
#include <c10/util/C++17.h>
#include <c10/util/MathConstants.h>
#include <algorithm>
namespace at {
@ -153,7 +154,7 @@ void CPUGeneratorImpl::set_state(const c10::TensorImpl& new_state) {
// intermediate values.
if (legacy_pod->normal_is_valid) {
auto r = legacy_pod->normal_rho;
auto theta = 2.0 * M_PI * legacy_pod->normal_x;
auto theta = 2.0 * c10::pi<double> * legacy_pod->normal_x;
// we return the sin version of the normal sample when in caching mode
double_normal_sample = c10::optional<double>(r * ::sin(theta));
}

View file

@ -1,17 +1,10 @@
#pragma once
// define constants like M_PI and C keywords for MSVC
#ifdef _MSC_VER
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES
#endif
#include <math.h>
#endif
#include <ATen/core/Array.h>
#include <ATen/core/TransformationHelper.h>
#include <c10/util/Half.h>
#include <c10/util/BFloat16.h>
#include <c10/util/MathConstants.h>
#include <c10/util/Optional.h>
#include <c10/macros/Macros.h>
@ -220,7 +213,7 @@ struct normal_distribution {
const dist_acctype<T> u1 = uniform(generator);
const dist_acctype<T> u2 = uniform(generator);
const dist_acctype<T> r = ::sqrt(static_cast<T>(-2.0) * ::log(static_cast<T>(1.0)-u2));
const dist_acctype<T> theta = static_cast<T>(2.0) * static_cast<T>(M_PI) * u1;
const dist_acctype<T> theta = static_cast<T>(2.0) * c10::pi<T> * u1;
if (std::is_same<T, double>::value) {
maybe_set_next_double_normal_sample(generator, r * ::sin(theta));
} else {

View file

@ -1,6 +1,7 @@
#include <c10/macros/Macros.h>
#include <c10/util/Half.h>
#include <c10/util/BFloat16.h>
#include <c10/util/MathConstants.h>
#include <ATen/NumericUtils.h>
#include <limits>
#include <cstdint>
@ -101,7 +102,7 @@ C10_HOST_DEVICE inline T normal(T val, T mean, T std) {
template <typename T>
C10_HOST_DEVICE inline T cauchy(T val, T median, T sigma) {
// https://en.wikipedia.org/wiki/Cauchy_distribution#Cumulative_distribution_function
return median + sigma * at::tan(static_cast<T>(M_PI) * (val - static_cast<T>(0.5)));
return median + sigma * at::tan(c10::pi<T> * (val - static_cast<T>(0.5)));
}
/**

View file

@ -210,7 +210,7 @@ public:
const auto nan_vec = _mm256_set1_ps(NAN);
const auto not_nan_mask = _mm256_cmp_ps(values, values, _CMP_EQ_OQ);
const auto nan_mask = _mm256_cmp_ps(not_nan_mask, zero_vec, _CMP_EQ_OQ);
const auto pi = _mm256_set1_ps(M_PI);
const auto pi = _mm256_set1_ps(c10::pi<float>);
const auto neg_mask = _mm256_cmp_ps(values, zero_vec, _CMP_LT_OQ);
auto angle = _mm256_blendv_ps(zero_vec, pi, neg_mask);

View file

@ -204,7 +204,8 @@ public:
}
Vec256<c10::complex<double>> acos() const {
// acos(x) = pi/2 - asin(x)
const __m256d pi_2 = _mm256_setr_pd(M_PI/2, 0.0, M_PI/2, 0.0);
constexpr auto pi_2d = c10::pi<double> / 2;
const __m256d pi_2 = _mm256_setr_pd(pi_2d, 0.0, pi_2d, 0.0);
return _mm256_sub_pd(pi_2, asin());
}
Vec256<c10::complex<double>> atan() const;

View file

@ -242,7 +242,8 @@ public:
}
Vec256<c10::complex<float>> acos() const {
// acos(x) = pi/2 - asin(x)
const __m256 pi_2 = _mm256_setr_ps(M_PI/2, 0.0, M_PI/2, 0.0, M_PI/2, 0.0, M_PI/2, 0.0);
constexpr float pi_2f = c10::pi<float> / 2;
const __m256 pi_2 = _mm256_setr_ps(pi_2f, 0.0, pi_2f, 0.0, pi_2f, 0.0, pi_2f, 0.0);
return _mm256_sub_ps(pi_2, asin());
}
Vec256<c10::complex<float>> atan() const;

View file

@ -112,7 +112,7 @@ public:
const auto nan_vec = _mm256_set1_pd(NAN);
const auto not_nan_mask = _mm256_cmp_pd(values, values, _CMP_EQ_OQ);
const auto nan_mask = _mm256_cmp_pd(not_nan_mask, zero_vec, _CMP_EQ_OQ);
const auto pi = _mm256_set1_pd(M_PI);
const auto pi = _mm256_set1_pd(c10::pi<double>);
const auto neg_mask = _mm256_cmp_pd(values, zero_vec, _CMP_LT_OQ);
auto angle = _mm256_blendv_pd(zero_vec, pi, neg_mask);

View file

@ -119,7 +119,7 @@ public:
const auto nan_vec = _mm256_set1_ps(NAN);
const auto not_nan_mask = _mm256_cmp_ps(values, values, _CMP_EQ_OQ);
const auto nan_mask = _mm256_cmp_ps(not_nan_mask, zero_vec, _CMP_EQ_OQ);
const auto pi = _mm256_set1_ps(M_PI);
const auto pi = _mm256_set1_ps(c10::pi<float>);
const auto neg_mask = _mm256_cmp_ps(values, zero_vec, _CMP_LT_OQ);
auto angle = _mm256_blendv_ps(zero_vec, pi, neg_mask);

View file

@ -3,6 +3,7 @@
#include <ATen/ATen.h>
#include <ATen/ExpandUtils.h>
#include <c10/macros/Macros.h>
#include <c10/util/MathConstants.h>
// ROCM hcc doesn't work well with using std:: in kernel functions
#if defined(__CUDA_ARCH__)
@ -276,8 +277,8 @@ C10_DEVICE static inline scalar_t digamma_one(scalar_t x) {
}
// it is more standard to write this as recursion, but
// nvcc does not like that
additional_summand = -static_cast<accscalar_t>(M_PI) /
compat_tan(static_cast<accscalar_t>(M_PI) * x);
additional_summand = -c10::pi<scalar_t> /
compat_tan(c10::pi<scalar_t> * x);
x = 1 - x;
}

View file

@ -1,10 +1,3 @@
// define constants like M_PI and C keywords for MSVC
#ifdef _MSC_VER
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES
#endif
#include <math.h>
#endif
#include <ATen/ATen.h>
#include <ATen/NativeFunctions.h>
#include <ATen/Dispatch.h>
@ -247,7 +240,7 @@ Tensor poisson_nll_loss(const Tensor& input, const Tensor& target, const bool lo
}
if (full) {
auto stirling_term = target * at::log(target) - target + 0.5 * at::log(2 * M_PI * target);
auto stirling_term = target * at::log(target) - target + 0.5 * at::log(2 * c10::pi<double> * target);
loss += stirling_term.masked_fill(target <= 1, 0);
}

View file

@ -6,11 +6,9 @@
#include <cfloat>
#include <limits>
#include <type_traits>
#include <c10/util/MathConstants.h>
#include <c10/util/math_compat.h>
#ifndef M_PIf
#define M_PIf 3.1415926535f
#endif // M_PIf
/* The next function is taken from https://github.com/antelopeusersgroup/antelope_contrib/blob/master/lib/location/libgenloc/erfinv.c.
Below is the copyright.
@ -105,8 +103,8 @@ Date: February 1996
#endif
}
/* Two steps of Newton-Raphson correction */
x = x - (std::erf(x) - y) / ((static_cast<T>(2.0)/static_cast<T>(std::sqrt(M_PI)))*std::exp(-x*x));
x = x - (std::erf(x) - y) / ((static_cast<T>(2.0)/static_cast<T>(std::sqrt(M_PI)))*std::exp(-x*x));
x = x - (std::erf(x) - y) / ((static_cast<T>(2.0)/static_cast<T>(std::sqrt(c10::pi<double>)))*std::exp(-x*x));
x = x - (std::erf(x) - y) / ((static_cast<T>(2.0)/static_cast<T>(std::sqrt(c10::pi<double>)))*std::exp(-x*x));
return(x);
}
@ -241,8 +239,8 @@ static inline double trigamma(double x) {
double result = 0;
if (x < 0.5) {
sign = -1;
const double sin_pi_x = sin(M_PI * x);
result -= (M_PI * M_PI) / (sin_pi_x * sin_pi_x);
const double sin_pi_x = sin(c10::pi<double> * x);
result -= (c10::pi<double> * c10::pi<double>) / (sin_pi_x * sin_pi_x);
x = 1 - x;
}
for (int i = 0; i < 6; ++i) {
@ -259,8 +257,8 @@ static inline float trigamma(float x) {
float result = 0;
if (x < 0.5f) {
sign = -1;
const float sin_pi_x = sinf(M_PIf * x);
result -= (M_PIf * M_PIf) / (sin_pi_x * sin_pi_x);
const float sin_pi_x = sinf(c10::pi<float> * x);
result -= (c10::pi<float> * c10::pi<float>) / (sin_pi_x * sin_pi_x);
x = 1 - x;
}
for (int i = 0; i < 6; ++i) {
@ -292,7 +290,7 @@ static inline double calc_digamma(double x) {
// If the argument is a negative integer, NaN is returned
return std::numeric_limits<double>::quiet_NaN();
}
return calc_digamma(1 - x) - M_PI / tan(M_PI * x);
return calc_digamma(1 - x) - c10::pi<double> / tan(c10::pi<double> * x);
}
// Push x to be >= 10
@ -346,7 +344,7 @@ static inline float calc_digamma(float x) {
}
// Avoid rounding errors for `tan`'s input.
// Those make a big difference at extreme values.
float pi_over_tan_pi_x = (float)(M_PI / tan(M_PI * (double)x));
float pi_over_tan_pi_x = (float)(c10::pi<double> / tan(c10::pi<double> * (double)x));
return calc_digamma(1 - x) - pi_over_tan_pi_x;
}
@ -883,7 +881,7 @@ static scalar_t _igam_helper_asymptotic_series(scalar_t a, scalar_t x, bool igam
absoldterm = absterm;
afac /= a;
}
res += sgn * std::exp(-0.5 * a * eta * eta) * sum / std::sqrt(2 * M_PIf * a);
res += sgn * std::exp(-0.5 * a * eta * eta) * sum / std::sqrt(2 * c10::pi<float> * a);
return res;
}

View file

@ -1,11 +1,3 @@
// define constants like M_PI and C keywords for MSVC
#ifdef _MSC_VER
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES
#endif
#include <math.h>
#endif
#include <ATen/ATen.h>
#include <ATen/Config.h>
#include <ATen/NativeFunctions.h>

View file

@ -1,11 +1,3 @@
// define constants like M_PI and C keywords for MSVC
#ifdef _MSC_VER
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES
#endif
#include <math.h>
#endif
#include <ATen/ATen.h>
#include <ATen/CPUGeneratorImpl.h>
#include <ATen/Utils.h>
@ -14,6 +6,7 @@
#include <ATen/TracerMode.h>
#include <c10/core/ScalarType.h>
#include <c10/util/Deprecated.h>
#include <ATen/native/Math.h>
#include <ATen/native/Resize.h>
#include <ATen/native/TensorFactories.h>
#include <c10/core/TensorOptions.h>
@ -915,7 +908,7 @@ Tensor blackman_window(
window_length += 1;
}
// from https://en.wikipedia.org/wiki/Window_function#Blackman_window
auto window = native::arange(window_length, options).mul_(M_PI / static_cast<double>(window_length - 1));
auto window = native::arange(window_length, options).mul_(c10::pi<double> / static_cast<double>(window_length - 1));
window = window.mul(4).cos_().mul_(0.08) - window.mul(2).cos_().mul_(0.5) + 0.42;
return periodic ? window.narrow(0, 0, window_length - 1) : window;
}
@ -960,7 +953,7 @@ Tensor hamming_window(
window_length += 1;
}
auto window = native::arange(window_length, options);
window.mul_(M_PI * 2. / static_cast<double>(window_length - 1)).cos_().mul_(-beta).add_(alpha);
window.mul_(c10::pi<double> * 2. / static_cast<double>(window_length - 1)).cos_().mul_(-beta).add_(alpha);
return periodic ? window.narrow(0, 0, window_length - 1) : window;
}

View file

@ -1,11 +1,3 @@
// define constants like M_PI and C keywords for MSVC
#ifdef _MSC_VER
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES
#endif
#include <math.h>
#endif
#include <ATen/ATen.h>
#include <ATen/Dispatch.h>
#include <ATen/ExpandUtils.h>
@ -16,6 +8,7 @@
#include <ATen/CPUApplyUtils.h>
#include <ATen/Parallel.h>
#include <ATen/native/Math.h>
#include <ATen/native/UnaryOps.h>
#include <ATen/native/TensorIterator.h>
#include <ATen/NamedTensorUtils.h>
@ -647,14 +640,14 @@ Tensor mvlgamma(const Tensor& self, int64_t p) {
mvlgamma_check(self, p);
Tensor args = native::arange(-p / 2. + 0.5, 0.5, 0.5, self.options());
args = args.add(self.unsqueeze(-1));
return args.lgamma_().sum(-1).add_(p * (p - 1) * std::log(M_PI) / 4.);
return args.lgamma_().sum(-1).add_(p * (p - 1) * std::log(c10::pi<double>) / 4.);
}
Tensor& mvlgamma_(Tensor& self, int64_t p) {
mvlgamma_check(self, p);
Tensor args = native::arange(-p / 2. + 0.5, 0.5, 0.5, self.options());
args = args.add(self.unsqueeze(-1));
return self.copy_(args.lgamma_().sum(-1).add_(p * (p - 1) * std::log(M_PI) / 4.));
return self.copy_(args.lgamma_().sum(-1).add_(p * (p - 1) * std::log(c10::pi<double>) / 4.));
}
// NB: If you use this macro, you may also need to add a CUDA forwarding

View file

@ -112,7 +112,7 @@ void normal_fill_AVX2(Tensor& self, const float mean, const float std, RNG gener
at::uniform_real_distribution<float> uniform(0, 1);
data[i] = uniform(generator);
}
const __m256 two_pi = _mm256_set1_ps(2.0f * M_PI);
const __m256 two_pi = _mm256_set1_ps(2.0f * c10::pi<double>);
const __m256 one = _mm256_set1_ps(1.0f);
const __m256 minus_two = _mm256_set1_ps(-2.0f);
const __m256 mean_v = _mm256_set1_ps(mean);
@ -140,7 +140,7 @@ static void normal_fill_16(scalar_t *data, const scalar_t mean, const scalar_t s
const scalar_t u1 = 1 - data[j]; // [0, 1) -> (0, 1] for log.
const scalar_t u2 = data[j + 8];
const scalar_t radius = std::sqrt(-2 * std::log(u1));
const scalar_t theta = 2.0f * M_PI * u2;
const scalar_t theta = 2.0f * c10::pi<double> * u2;
data[j] = radius * std::cos(theta) * std + mean;
data[j + 8] = radius * std::sin(theta) * std + mean;
}

View file

@ -16,11 +16,11 @@
#include <ATen/cpu/vml.h>
#include <ATen/native/Distributions.h>
#include <ATen/native/TensorFactories.h>
#include <ATen/native/Math.h>
#include <ATen/native/TensorIterator.h>
#include <ATen/native/cpu/DistributionTemplates.h>
#include <ATen/native/cpu/Loops.h>
#include <ATen/native/cpu/zmath.h>
#include <c10/util/MathConstants.h>
#if AT_MKL_ENABLED()
#include <mkl.h>
@ -310,7 +310,7 @@ static void sinc_kernel(TensorIterator& iter) {
if (a == scalar_t(0)) {
return scalar_t(1);
} else {
scalar_t product = scalar_t(M_PI) * a;
scalar_t product = c10::pi<scalar_t> * a;
return std::sin(product) / product;
}
});

View file

@ -3,6 +3,7 @@
// Complex number math operations that act as no-ops for other dtypes.
#include <c10/util/complex.h>
#include <c10/util/math_compat.h>
#include <c10/util/MathConstants.h>
#include<ATen/NumericUtils.h>
namespace at { namespace native {
@ -44,7 +45,7 @@ inline VALUE_TYPE angle_impl (SCALAR_TYPE z) {
if (at::_isnan(z)) {
return z;
}
return z < 0 ? M_PI : 0;
return z < 0 ? c10::pi<double> : 0;
}
template<>

View file

@ -6,6 +6,7 @@
#include <ATen/Context.h>
#include <ATen/Dispatch.h>
#include <ATen/native/DispatchStub.h>
#include <ATen/native/Math.h>
#include <ATen/native/TensorFactories.h>
#include <ATen/native/TensorIterator.h>
#include <ATen/native/cuda/Loops.cuh>
@ -109,7 +110,8 @@ void sinc_kernel_cuda(TensorIterator& iter) {
if (a == scalar_t(0)) {
return scalar_t(1);
} else {
scalar_t product = scalar_t(M_PI) * a;
// NVCC says constexpr var is not accessible from device
scalar_t product = c10::detail::pi<scalar_t>() * a;
return std::sin(product) / product;
}
});

View file

@ -69,10 +69,6 @@
# define TH_UNUSED
#endif
#ifndef M_PI
# define M_PI 3.14159265358979323846
#endif
typedef void (*THErrorHandlerFunction)(const char *msg, void *data);
typedef void (*THArgErrorHandlerFunction)(int argNumber, const char *msg, void *data);

View file

@ -0,0 +1,12 @@
// define constants like M_PI and C keywords for MSVC
#ifdef _MSC_VER
#ifndef _USE_MATH_DEFINES
#define _USE_MATH_DEFINES
#endif
#endif
#include <c10/util/MathConstants.h>
#include <math.h>
static_assert(M_PI == c10::pi<double>, "c10::pi<double> must be equal to M_PI");

30
c10/util/MathConstants.h Normal file
View file

@ -0,0 +1,30 @@
#pragma once
#include <c10/macros/Macros.h>
#include <c10/util/BFloat16.h>
#include <c10/util/Half.h>
namespace c10 {
// TODO: Replace me with inline constexpr variable when C++17 becomes available
namespace detail {
template <typename T>
C10_HOST_DEVICE inline constexpr T pi() {
return static_cast<T>(3.14159265358979323846L);
}
template<>
C10_HOST_DEVICE inline constexpr BFloat16 pi<BFloat16>() {
// According to https://en.wikipedia.org/wiki/Bfloat16_floating-point_format#Special_values pi is encoded as 4049
return BFloat16(0x4049, BFloat16::from_bits());
}
template<>
C10_HOST_DEVICE inline constexpr Half pi<Half>() {
return Half(0x4248, Half::from_bits());
}
} // namespace detail
// TODO: Replace me with std::numbers::pi when C++20 is there
template<typename T>
constexpr T pi = c10::detail::pi<T>();
} // namespace c10