123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205 |
- // (C) Copyright John Maddock 2006.
- // Use, modification and distribution are subject to the
- // Boost Software License, Version 1.0. (See accompanying file
- // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
- #ifndef BOOST_MATH_SPECIAL_SPHERICAL_HARMONIC_HPP
- #define BOOST_MATH_SPECIAL_SPHERICAL_HARMONIC_HPP
- #ifdef _MSC_VER
- #pragma once
- #endif
- #include <boost/math/special_functions/math_fwd.hpp>
- #include <boost/math/special_functions/legendre.hpp>
- #include <boost/math/tools/workaround.hpp>
- #include <complex>
- namespace boost{
- namespace math{
- namespace detail{
- //
- // Calculates the prefix term that's common to the real
- // and imaginary parts. Does *not* fix up the sign of the result
- // though.
- //
- template <class T, class Policy>
- inline T spherical_harmonic_prefix(unsigned n, unsigned m, T theta, const Policy& pol)
- {
- BOOST_MATH_STD_USING
- if(m > n)
- return 0;
- T sin_theta = sin(theta);
- T x = cos(theta);
- T leg = detail::legendre_p_imp(n, m, x, static_cast<T>(pow(fabs(sin_theta), T(m))), pol);
-
- T prefix = boost::math::tgamma_delta_ratio(static_cast<T>(n - m + 1), static_cast<T>(2 * m), pol);
- prefix *= (2 * n + 1) / (4 * constants::pi<T>());
- prefix = sqrt(prefix);
- return prefix * leg;
- }
- //
- // Real Part:
- //
- template <class T, class Policy>
- T spherical_harmonic_r(unsigned n, int m, T theta, T phi, const Policy& pol)
- {
- BOOST_MATH_STD_USING // ADL of std functions
- bool sign = false;
- if(m < 0)
- {
- // Reflect and adjust sign if m < 0:
- sign = m&1;
- m = abs(m);
- }
- if(m&1)
- {
- // Check phase if theta is outside [0, PI]:
- T mod = boost::math::tools::fmod_workaround(theta, T(2 * constants::pi<T>()));
- if(mod < 0)
- mod += 2 * constants::pi<T>();
- if(mod > constants::pi<T>())
- sign = !sign;
- }
- // Get the value and adjust sign as required:
- T prefix = spherical_harmonic_prefix(n, m, theta, pol);
- prefix *= cos(m * phi);
- return sign ? T(-prefix) : prefix;
- }
- template <class T, class Policy>
- T spherical_harmonic_i(unsigned n, int m, T theta, T phi, const Policy& pol)
- {
- BOOST_MATH_STD_USING // ADL of std functions
- bool sign = false;
- if(m < 0)
- {
- // Reflect and adjust sign if m < 0:
- sign = !(m&1);
- m = abs(m);
- }
- if(m&1)
- {
- // Check phase if theta is outside [0, PI]:
- T mod = boost::math::tools::fmod_workaround(theta, T(2 * constants::pi<T>()));
- if(mod < 0)
- mod += 2 * constants::pi<T>();
- if(mod > constants::pi<T>())
- sign = !sign;
- }
- // Get the value and adjust sign as required:
- T prefix = spherical_harmonic_prefix(n, m, theta, pol);
- prefix *= sin(m * phi);
- return sign ? T(-prefix) : prefix;
- }
- template <class T, class U, class Policy>
- std::complex<T> spherical_harmonic(unsigned n, int m, U theta, U phi, const Policy& pol)
- {
- BOOST_MATH_STD_USING
- //
- // Sort out the signs:
- //
- bool r_sign = false;
- bool i_sign = false;
- if(m < 0)
- {
- // Reflect and adjust sign if m < 0:
- r_sign = m&1;
- i_sign = !(m&1);
- m = abs(m);
- }
- if(m&1)
- {
- // Check phase if theta is outside [0, PI]:
- U mod = boost::math::tools::fmod_workaround(theta, U(2 * constants::pi<U>()));
- if(mod < 0)
- mod += 2 * constants::pi<U>();
- if(mod > constants::pi<U>())
- {
- r_sign = !r_sign;
- i_sign = !i_sign;
- }
- }
- //
- // Calculate the value:
- //
- U prefix = spherical_harmonic_prefix(n, m, theta, pol);
- U r = prefix * cos(m * phi);
- U i = prefix * sin(m * phi);
- //
- // Add in the signs:
- //
- if(r_sign)
- r = -r;
- if(i_sign)
- i = -i;
- static const char* function = "boost::math::spherical_harmonic<%1%>(int, int, %1%, %1%)";
- return std::complex<T>(policies::checked_narrowing_cast<T, Policy>(r, function), policies::checked_narrowing_cast<T, Policy>(i, function));
- }
- } // namespace detail
- template <class T1, class T2, class Policy>
- inline std::complex<typename tools::promote_args<T1, T2>::type>
- spherical_harmonic(unsigned n, int m, T1 theta, T2 phi, const Policy& pol)
- {
- typedef typename tools::promote_args<T1, T2>::type result_type;
- typedef typename policies::evaluation<result_type, Policy>::type value_type;
- return detail::spherical_harmonic<result_type, value_type>(n, m, static_cast<value_type>(theta), static_cast<value_type>(phi), pol);
- }
- template <class T1, class T2>
- inline std::complex<typename tools::promote_args<T1, T2>::type>
- spherical_harmonic(unsigned n, int m, T1 theta, T2 phi)
- {
- return boost::math::spherical_harmonic(n, m, theta, phi, policies::policy<>());
- }
- template <class T1, class T2, class Policy>
- inline typename tools::promote_args<T1, T2>::type
- spherical_harmonic_r(unsigned n, int m, T1 theta, T2 phi, const Policy& pol)
- {
- typedef typename tools::promote_args<T1, T2>::type result_type;
- typedef typename policies::evaluation<result_type, Policy>::type value_type;
- return policies::checked_narrowing_cast<result_type, Policy>(detail::spherical_harmonic_r(n, m, static_cast<value_type>(theta), static_cast<value_type>(phi), pol), "bost::math::spherical_harmonic_r<%1%>(unsigned, int, %1%, %1%)");
- }
- template <class T1, class T2>
- inline typename tools::promote_args<T1, T2>::type
- spherical_harmonic_r(unsigned n, int m, T1 theta, T2 phi)
- {
- return boost::math::spherical_harmonic_r(n, m, theta, phi, policies::policy<>());
- }
- template <class T1, class T2, class Policy>
- inline typename tools::promote_args<T1, T2>::type
- spherical_harmonic_i(unsigned n, int m, T1 theta, T2 phi, const Policy& pol)
- {
- typedef typename tools::promote_args<T1, T2>::type result_type;
- typedef typename policies::evaluation<result_type, Policy>::type value_type;
- return policies::checked_narrowing_cast<result_type, Policy>(detail::spherical_harmonic_i(n, m, static_cast<value_type>(theta), static_cast<value_type>(phi), pol), "boost::math::spherical_harmonic_i<%1%>(unsigned, int, %1%, %1%)");
- }
- template <class T1, class T2>
- inline typename tools::promote_args<T1, T2>::type
- spherical_harmonic_i(unsigned n, int m, T1 theta, T2 phi)
- {
- return boost::math::spherical_harmonic_i(n, m, theta, phi, policies::policy<>());
- }
- } // namespace math
- } // namespace boost
- #endif // BOOST_MATH_SPECIAL_SPHERICAL_HARMONIC_HPP
|