/usr/include/boost/math/special_functions
Edit: /usr/include/boost/math/special_functions/spherical_harmonic.hpp (6358B)
// (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
#include
#include
#include
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
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(pow(fabs(sin_theta), T(m))), pol);
T prefix = boost::math::tgamma_delta_ratio(static_cast(n - m + 1), static_cast(2 * m), pol);
prefix *= (2 * n + 1) / (4 * constants::pi());
prefix = sqrt(prefix);
return prefix * leg;
}
//
// Real Part:
//
template
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()));
if(mod < 0)
mod += 2 * constants::pi();
if(mod > constants::pi())
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
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()));
if(mod < 0)
mod += 2 * constants::pi();
if(mod > constants::pi())
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
std::complex 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()));
if(mod < 0)
mod += 2 * constants::pi();
if(mod > constants::pi())
{
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(policies::checked_narrowing_cast(r, function), policies::checked_narrowing_cast(i, function));
}
} // namespace detail
template
inline std::complex::type>
spherical_harmonic(unsigned n, int m, T1 theta, T2 phi, const Policy& pol)
{
typedef typename tools::promote_args::type result_type;
typedef typename policies::evaluation::type value_type;
return detail::spherical_harmonic(n, m, static_cast(theta), static_cast(phi), pol);
}
template
inline std::complex::type>
spherical_harmonic(unsigned n, int m, T1 theta, T2 phi)
{
return boost::math::spherical_harmonic(n, m, theta, phi, policies::policy<>());
}
template
inline typename tools::promote_args::type
spherical_harmonic_r(unsigned n, int m, T1 theta, T2 phi, const Policy& pol)
{
typedef typename tools::promote_args::type result_type;
typedef typename policies::evaluation::type value_type;
return policies::checked_narrowing_cast(detail::spherical_harmonic_r(n, m, static_cast(theta), static_cast(phi), pol), "bost::math::spherical_harmonic_r<%1%>(unsigned, int, %1%, %1%)");
}
template
inline typename tools::promote_args::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
inline typename tools::promote_args::type
spherical_harmonic_i(unsigned n, int m, T1 theta, T2 phi, const Policy& pol)
{
typedef typename tools::promote_args::type result_type;
typedef typename policies::evaluation::type value_type;
return policies::checked_narrowing_cast(detail::spherical_harmonic_i(n, m, static_cast(theta), static_cast(phi), pol), "boost::math::spherical_harmonic_i<%1%>(unsigned, int, %1%, %1%)");
}
template
inline typename tools::promote_args::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