sinhc.hpp 4.2 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142
  1. // boost sinhc.hpp header file
  2. // (C) Copyright Hubert Holin 2001.
  3. // Distributed under the Boost Software License, Version 1.0. (See
  4. // accompanying file LICENSE_1_0.txt or copy at
  5. // http://www.boost.org/LICENSE_1_0.txt)
  6. // See http://www.boost.org for updates, documentation, and revision history.
  7. #ifndef BOOST_SINHC_HPP
  8. #define BOOST_SINHC_HPP
  9. #ifdef _MSC_VER
  10. #pragma once
  11. #endif
  12. #include <boost/math/tools/precision.hpp>
  13. #include <boost/math/policies/error_handling.hpp>
  14. #include <boost/math/special_functions/math_fwd.hpp>
  15. #include <boost/math/special_functions/fpclassify.hpp>
  16. #include <limits>
  17. #include <string>
  18. #include <stdexcept>
  19. #include <cmath>
  20. // These are the the "Hyperbolic Sinus Cardinal" functions.
  21. namespace boost
  22. {
  23. namespace math
  24. {
  25. namespace detail
  26. {
  27. // This is the "Hyperbolic Sinus Cardinal" of index Pi.
  28. template<typename T, typename Policy>
  29. inline T sinhc_pi_imp(const T x, const Policy&)
  30. {
  31. using ::std::abs;
  32. using ::std::sinh;
  33. using ::std::sqrt;
  34. static T const taylor_0_bound = tools::epsilon<T>();
  35. static T const taylor_2_bound = sqrt(taylor_0_bound);
  36. static T const taylor_n_bound = sqrt(taylor_2_bound);
  37. if((boost::math::isinf)(x))
  38. {
  39. return policies::raise_overflow_error<T>("sinhc(%1%)", nullptr, Policy());
  40. }
  41. if (abs(x) >= taylor_n_bound)
  42. {
  43. return(sinh(x)/x);
  44. }
  45. else
  46. {
  47. // approximation by taylor series in x at 0 up to order 0
  48. T result = static_cast<T>(1);
  49. if (abs(x) >= taylor_0_bound)
  50. {
  51. T x2 = x*x;
  52. // approximation by taylor series in x at 0 up to order 2
  53. result += x2/static_cast<T>(6);
  54. if (abs(x) >= taylor_2_bound)
  55. {
  56. // approximation by taylor series in x at 0 up to order 4
  57. result += (x2*x2)/static_cast<T>(120);
  58. }
  59. }
  60. return(result);
  61. }
  62. }
  63. } // namespace detail
  64. template <class T, class Policy>
  65. inline typename tools::promote_args<T>::type sinhc_pi(T x, const Policy& pol)
  66. {
  67. typedef typename tools::promote_args<T>::type result_type;
  68. return policies::checked_narrowing_cast<T, Policy>(detail::sinhc_pi_imp(static_cast<result_type>(x), pol), "sinhc(%1%)");
  69. }
  70. template <class T>
  71. inline typename tools::promote_args<T>::type sinhc_pi(T x)
  72. {
  73. typedef typename tools::promote_args<T>::type result_type;
  74. return sinhc_pi(static_cast<result_type>(x), policies::policy<>());
  75. }
  76. template<typename T, template<typename> class U>
  77. inline U<T> sinhc_pi(const U<T> x)
  78. {
  79. using std::abs;
  80. using std::sinh;
  81. using std::sqrt;
  82. using ::std::numeric_limits;
  83. static T const taylor_0_bound = tools::epsilon<T>();
  84. static T const taylor_2_bound = sqrt(taylor_0_bound);
  85. static T const taylor_n_bound = sqrt(taylor_2_bound);
  86. if (abs(x) >= taylor_n_bound)
  87. {
  88. return(sinh(x)/x);
  89. }
  90. else
  91. {
  92. // approximation by taylor series in x at 0 up to order 0
  93. #ifdef __MWERKS__
  94. U<T> result = static_cast<U<T> >(1);
  95. #else
  96. U<T> result = U<T>(1);
  97. #endif
  98. if (abs(x) >= taylor_0_bound)
  99. {
  100. U<T> x2 = x*x;
  101. // approximation by taylor series in x at 0 up to order 2
  102. result += x2/static_cast<T>(6);
  103. if (abs(x) >= taylor_2_bound)
  104. {
  105. // approximation by taylor series in x at 0 up to order 4
  106. result += (x2*x2)/static_cast<T>(120);
  107. }
  108. }
  109. return(result);
  110. }
  111. }
  112. }
  113. }
  114. #endif /* BOOST_SINHC_HPP */