transc.hpp 7.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232
  1. /* Boost interval/transc.hpp template implementation file
  2. *
  3. * Copyright 2000 Jens Maurer
  4. * Copyright 2002 Hervé Brönnimann, Guillaume Melquiond, Sylvain Pion
  5. *
  6. * Distributed under the Boost Software License, Version 1.0.
  7. * (See accompanying file LICENSE_1_0.txt or
  8. * copy at http://www.boost.org/LICENSE_1_0.txt)
  9. */
  10. #ifndef BOOST_NUMERIC_INTERVAL_TRANSC_HPP
  11. #define BOOST_NUMERIC_INTERVAL_TRANSC_HPP
  12. #include <boost/config.hpp>
  13. #include <boost/numeric/interval/detail/interval_prototype.hpp>
  14. #include <boost/numeric/interval/detail/bugs.hpp>
  15. #include <boost/numeric/interval/detail/test_input.hpp>
  16. #include <boost/numeric/interval/rounding.hpp>
  17. #include <boost/numeric/interval/constants.hpp>
  18. #include <boost/numeric/interval/arith.hpp>
  19. #include <boost/numeric/interval/arith2.hpp>
  20. #include <algorithm>
  21. namespace boost {
  22. namespace numeric {
  23. template<class T, class Policies> inline
  24. interval<T, Policies> exp(const interval<T, Policies>& x)
  25. {
  26. typedef interval<T, Policies> I;
  27. if (interval_lib::detail::test_input(x))
  28. return I::empty();
  29. typename Policies::rounding rnd;
  30. return I(rnd.exp_down(x.lower()), rnd.exp_up(x.upper()), true);
  31. }
  32. template<class T, class Policies> inline
  33. interval<T, Policies> log(const interval<T, Policies>& x)
  34. {
  35. typedef interval<T, Policies> I;
  36. if (interval_lib::detail::test_input(x) ||
  37. !interval_lib::user::is_pos(x.upper()))
  38. return I::empty();
  39. typename Policies::rounding rnd;
  40. typedef typename Policies::checking checking;
  41. T l = !interval_lib::user::is_pos(x.lower())
  42. ? checking::neg_inf() : rnd.log_down(x.lower());
  43. return I(l, rnd.log_up(x.upper()), true);
  44. }
  45. template<class T, class Policies> inline
  46. interval<T, Policies> cos(const interval<T, Policies>& x)
  47. {
  48. if (interval_lib::detail::test_input(x))
  49. return interval<T, Policies>::empty();
  50. typename Policies::rounding rnd;
  51. typedef interval<T, Policies> I;
  52. typedef typename interval_lib::unprotect<I>::type R;
  53. // get lower bound within [0, pi]
  54. const R pi2 = interval_lib::pi_twice<R>();
  55. R tmp = fmod((const R&)x, pi2);
  56. if (width(tmp) >= pi2.lower())
  57. return I(static_cast<T>(-1), static_cast<T>(1), true); // we are covering a full period
  58. if (tmp.lower() >= interval_lib::constants::pi_upper<T>())
  59. return -cos(tmp - interval_lib::pi<R>());
  60. T l = tmp.lower();
  61. T u = tmp.upper();
  62. BOOST_USING_STD_MIN();
  63. // separate into monotone subintervals
  64. if (u <= interval_lib::constants::pi_lower<T>())
  65. return I(rnd.cos_down(u), rnd.cos_up(l), true);
  66. else if (u <= pi2.lower())
  67. return I(static_cast<T>(-1), rnd.cos_up(min BOOST_PREVENT_MACRO_SUBSTITUTION(rnd.sub_down(pi2.lower(), u), l)), true);
  68. else
  69. return I(static_cast<T>(-1), static_cast<T>(1), true);
  70. }
  71. template<class T, class Policies> inline
  72. interval<T, Policies> sin(const interval<T, Policies>& x)
  73. {
  74. typedef interval<T, Policies> I;
  75. if (interval_lib::detail::test_input(x))
  76. return I::empty();
  77. typename Policies::rounding rnd;
  78. typedef typename interval_lib::unprotect<I>::type R;
  79. I r = cos((const R&)x - interval_lib::pi_half<R>());
  80. (void)&rnd;
  81. return r;
  82. }
  83. template<class T, class Policies> inline
  84. interval<T, Policies> tan(const interval<T, Policies>& x)
  85. {
  86. typedef interval<T, Policies> I;
  87. if (interval_lib::detail::test_input(x))
  88. return I::empty();
  89. typename Policies::rounding rnd;
  90. typedef typename interval_lib::unprotect<I>::type R;
  91. // get lower bound within [-pi/2, pi/2]
  92. const R pi = interval_lib::pi<R>();
  93. R tmp = fmod((const R&)x, pi);
  94. const T pi_half_d = interval_lib::constants::pi_half_lower<T>();
  95. if (tmp.lower() >= pi_half_d)
  96. tmp -= pi;
  97. if (tmp.lower() <= -pi_half_d || tmp.upper() >= pi_half_d)
  98. return I::whole();
  99. return I(rnd.tan_down(tmp.lower()), rnd.tan_up(tmp.upper()), true);
  100. }
  101. template<class T, class Policies> inline
  102. interval<T, Policies> asin(const interval<T, Policies>& x)
  103. {
  104. typedef interval<T, Policies> I;
  105. if (interval_lib::detail::test_input(x)
  106. || x.upper() < static_cast<T>(-1) || x.lower() > static_cast<T>(1))
  107. return I::empty();
  108. typename Policies::rounding rnd;
  109. T l = (x.lower() <= static_cast<T>(-1))
  110. ? -interval_lib::constants::pi_half_upper<T>()
  111. : rnd.asin_down(x.lower());
  112. T u = (x.upper() >= static_cast<T>(1) )
  113. ? interval_lib::constants::pi_half_upper<T>()
  114. : rnd.asin_up (x.upper());
  115. return I(l, u, true);
  116. }
  117. template<class T, class Policies> inline
  118. interval<T, Policies> acos(const interval<T, Policies>& x)
  119. {
  120. typedef interval<T, Policies> I;
  121. if (interval_lib::detail::test_input(x)
  122. || x.upper() < static_cast<T>(-1) || x.lower() > static_cast<T>(1))
  123. return I::empty();
  124. typename Policies::rounding rnd;
  125. T l = (x.upper() >= static_cast<T>(1) )
  126. ? static_cast<T>(0)
  127. : rnd.acos_down(x.upper());
  128. T u = (x.lower() <= static_cast<T>(-1))
  129. ? interval_lib::constants::pi_upper<T>()
  130. : rnd.acos_up (x.lower());
  131. return I(l, u, true);
  132. }
  133. template<class T, class Policies> inline
  134. interval<T, Policies> atan(const interval<T, Policies>& x)
  135. {
  136. typedef interval<T, Policies> I;
  137. if (interval_lib::detail::test_input(x))
  138. return I::empty();
  139. typename Policies::rounding rnd;
  140. return I(rnd.atan_down(x.lower()), rnd.atan_up(x.upper()), true);
  141. }
  142. template<class T, class Policies> inline
  143. interval<T, Policies> sinh(const interval<T, Policies>& x)
  144. {
  145. typedef interval<T, Policies> I;
  146. if (interval_lib::detail::test_input(x))
  147. return I::empty();
  148. typename Policies::rounding rnd;
  149. return I(rnd.sinh_down(x.lower()), rnd.sinh_up(x.upper()), true);
  150. }
  151. template<class T, class Policies> inline
  152. interval<T, Policies> cosh(const interval<T, Policies>& x)
  153. {
  154. typedef interval<T, Policies> I;
  155. if (interval_lib::detail::test_input(x))
  156. return I::empty();
  157. typename Policies::rounding rnd;
  158. if (interval_lib::user::is_neg(x.upper()))
  159. return I(rnd.cosh_down(x.upper()), rnd.cosh_up(x.lower()), true);
  160. else if (!interval_lib::user::is_neg(x.lower()))
  161. return I(rnd.cosh_down(x.lower()), rnd.cosh_up(x.upper()), true);
  162. else
  163. return I(static_cast<T>(1), rnd.cosh_up(-x.lower() > x.upper() ? x.lower() : x.upper()), true);
  164. }
  165. template<class T, class Policies> inline
  166. interval<T, Policies> tanh(const interval<T, Policies>& x)
  167. {
  168. typedef interval<T, Policies> I;
  169. if (interval_lib::detail::test_input(x))
  170. return I::empty();
  171. typename Policies::rounding rnd;
  172. return I(rnd.tanh_down(x.lower()), rnd.tanh_up(x.upper()), true);
  173. }
  174. template<class T, class Policies> inline
  175. interval<T, Policies> asinh(const interval<T, Policies>& x)
  176. {
  177. typedef interval<T, Policies> I;
  178. if (interval_lib::detail::test_input(x))
  179. return I::empty();
  180. typename Policies::rounding rnd;
  181. return I(rnd.asinh_down(x.lower()), rnd.asinh_up(x.upper()), true);
  182. }
  183. template<class T, class Policies> inline
  184. interval<T, Policies> acosh(const interval<T, Policies>& x)
  185. {
  186. typedef interval<T, Policies> I;
  187. if (interval_lib::detail::test_input(x) || x.upper() < static_cast<T>(1))
  188. return I::empty();
  189. typename Policies::rounding rnd;
  190. T l = x.lower() <= static_cast<T>(1) ? static_cast<T>(0) : rnd.acosh_down(x.lower());
  191. return I(l, rnd.acosh_up(x.upper()), true);
  192. }
  193. template<class T, class Policies> inline
  194. interval<T, Policies> atanh(const interval<T, Policies>& x)
  195. {
  196. typedef interval<T, Policies> I;
  197. if (interval_lib::detail::test_input(x)
  198. || x.upper() < static_cast<T>(-1) || x.lower() > static_cast<T>(1))
  199. return I::empty();
  200. typename Policies::rounding rnd;
  201. typedef typename Policies::checking checking;
  202. T l = (x.lower() <= static_cast<T>(-1))
  203. ? checking::neg_inf() : rnd.atanh_down(x.lower());
  204. T u = (x.upper() >= static_cast<T>(1) )
  205. ? checking::pos_inf() : rnd.atanh_up (x.upper());
  206. return I(l, u, true);
  207. }
  208. } // namespace numeric
  209. } // namespace boost
  210. #endif // BOOST_NUMERIC_INTERVAL_TRANSC_HPP