boost/math/special_functions/sinc.hpp

   1 //  boost sinc.hpp header file
   2
   3 //  (C) Copyright Hubert Holin 2001.
   4 //  Distributed under the Boost Software License, Version 1.0. (See
   5 //  accompanying file LICENSE_1_0.txt or copy at
   6 //  http://www.boost.org/LICENSE_1_0.txt)
   7
   8 // See http://www.boost.org for updates, documentation, and revision history.
   9
  10 #ifndef BOOST_SINC_HPP
  11 #define BOOST_SINC_HPP
  12
  13
  14 #ifdef _MSC_VER
  15 #pragma once
  16 #endif
  17
  18 #include <boost/math/tools/config.hpp>
  19 #include <boost/math/tools/precision.hpp>
  20 #include <boost/math/policies/policy.hpp>
  21 #include <boost/math/special_functions/math_fwd.hpp>
  22 #include <boost/config/no_tr1/cmath.hpp>
  23 #include <boost/limits.hpp>
  24 #include <string>
  25 #include <stdexcept>
  26
  27
  28 #include <boost/config.hpp>
  29
  30
  31 // These are the the "Sinus Cardinal" functions.
  32
  33 namespace boost
  34 {
  35     namespace math
  36     {
  37        namespace detail
  38        {
  39         // This is the "Sinus Cardinal" of index Pi.
  40
  41         template<typename T>
  42         inline T    sinc_pi_imp(const T x)
  43         {
  44             BOOST_MATH_STD_USING
  45
  46             T const    taylor_0_bound = tools::epsilon<T>();
  47             T const    taylor_2_bound = tools::root_epsilon<T>();
  48             T const    taylor_n_bound = tools::forth_root_epsilon<T>();
  49
  50             if    (abs(x) >= taylor_n_bound)
  51             {
  52                 return(sin(x)/x);
  53             }
  54             else
  55             {
  56                 // approximation by taylor series in x at 0 up to order 0
  57                 T    result = static_cast<T>(1);
  58
  59                 if    (abs(x) >= taylor_0_bound)
  60                 {
  61                     T    x2 = x*x;
  62
  63                     // approximation by taylor series in x at 0 up to order 2
  64                     result -= x2/static_cast<T>(6);
  65
  66                     if    (abs(x) >= taylor_2_bound)
  67                     {
  68                         // approximation by taylor series in x at 0 up to order 4
  69                         result += (x2*x2)/static_cast<T>(120);
  70                     }
  71                 }
  72
  73                 return(result);
  74             }
  75         }
  76
  77        } // namespace detail
  78
  79        template <class T>
  80        inline typename tools::promote_args<T>::type sinc_pi(T x)
  81        {
  82           typedef typename tools::promote_args<T>::type result_type;
  83           return detail::sinc_pi_imp(static_cast<result_type>(x));
  84        }
  85
  86        template <class T, class Policy>
  87        inline typename tools::promote_args<T>::type sinc_pi(T x, const Policy&)
  88        {
  89           typedef typename tools::promote_args<T>::type result_type;
  90           return detail::sinc_pi_imp(static_cast<result_type>(x));
  91        }
  92
  93 #ifndef    BOOST_NO_TEMPLATE_TEMPLATES
  94         template<typename T, template<typename> class U>
  95         inline U<T>    sinc_pi(const U<T> x)
  96         {
  97             BOOST_MATH_STD_USING
  98             using    ::std::numeric_limits;
  99
 100             T const    taylor_0_bound = tools::epsilon<T>();
 101             T const    taylor_2_bound = tools::root_epsilon<T>();
 102             T const    taylor_n_bound = tools::forth_root_epsilon<T>();
 103
 104             if    (abs(x) >= taylor_n_bound)
 105             {
 106                 return(sin(x)/x);
 107             }
 108             else
 109             {
 110                 // approximation by taylor series in x at 0 up to order 0
 111 #ifdef __MWERKS__
 112                 U<T>    result = static_cast<U<T> >(1);
 113 #else
 114                 U<T>    result = U<T>(1);
 115 #endif
 116
 117                 if    (abs(x) >= taylor_0_bound)
 118                 {
 119                     U<T>    x2 = x*x;
 120
 121                     // approximation by taylor series in x at 0 up to order 2
 122                     result -= x2/static_cast<T>(6);
 123
 124                     if    (abs(x) >= taylor_2_bound)
 125                     {
 126                         // approximation by taylor series in x at 0 up to order 4
 127                         result += (x2*x2)/static_cast<T>(120);
 128                     }
 129                 }
 130
 131                 return(result);
 132             }
 133         }
 134
 135         template<typename T, template<typename> class U, class Policy>
 136         inline U<T>    sinc_pi(const U<T> x, const Policy&)
 137         {
 138            return sinc_pi(x);
 139         }
 140 #endif    /* BOOST_NO_TEMPLATE_TEMPLATES */
 141     }
 142 }
 143
 144 #endif /* BOOST_SINC_HPP */
 145