E.V.E  0.1-beta

◆ legendre

eve::legendre = {}
inlineconstexpr

Callable object evaluating the legendre functions.

Required header: #include <eve/function/legendre.hpp>

Members Functions

Member Effect
operator() the legendre operation

auto operator()( integral_value l, floating_value x) const noexcept;
auto operator()( integral_value l, integral_value m, floating_value x) const noexcept;
Definition: value.hpp:83
Definition: value.hpp:42

Parameters

l, m: real values. l and m must be [flint](eve::is_flint) or the result is U.B.

x: floating real value. x must be in the \([-1, 1]\) interval or the result is nan.

Return value

With two parameters returns elementwise the value of the first kind Legendre polynomial of order l at x.

With three parameters returns elementwise the value of the associated Legendre "polynomial" of order (l, m) at x.

The result type is of the index compatible type type of the parameters.

Warning
using float based inputs (instead of double) may introduce inaccuracies (peculiarly in the associated polynomials computation).

Supported decorators

  • eve::p_kind, eve::q_kind

    The expression p_kind(legendre)(n,x) is equivalent to legendre(n,x).

    The expression q_kind(legendre)(n,x) return the value at x of the second kind legendre function of order n.

  • eve::diff, eve::diff_1st, eve::diff_nth

    Required header: #include <eve/function/diff/legendre.hpp>

    The expression diff(legendre)(...,x) computes the derivative of the p_kind function relative to x.

  • eve::successor

    The expression successor(legendre)(l, x, ln, lnm1) (or successor(legendre)(l, m, x, ln, lnm1)) implements the three term recurrence relation for the (associated) Legendre polynomials, \(\displaystyle \mbox{P}^m_{l+1} = \left((2l+1)\mbox{P}^m_{l}(x)-l\mbox{P}^m_{l-1}(x)\right)/(l+m+1)\) This object function can be used to create a sequence of values evaluated at the same x and for rising l. ( \(m = 0\) and no \(m\) in call are equivalent here).

  • eve::condon_shortley

    The expression condon_shortley(legendre)(l, m, x) multiplies the associated legendre polynomial value by the Condon-Shortley phase \((-1)^m\) to match the definition given by Abramowitz and Stegun (8.6.6). This is currently the version implemented in boost::math.

Example

See it live on Compiler Explorer

#include <eve/function/legendre.hpp>
#include <eve/wide.hpp>
#include <eve/constant/inf.hpp>
#include <eve/constant/minf.hpp>
#include <eve/constant/nan.hpp>
#include <iostream>
int main()
{
wide_ft xd = {-0.1, -0.2, -0.3, -0.5, 0.0, 0.2, 0.3, 2.0};
wide_it n = {0, 1, 2, 3, 4, 5, 6, 7};
wide_ft x(0.5);
std::cout << "---- simd" << '\n'
<< "<- xd = " << xd << '\n'
<< "<- n = " << n << '\n'
<< "<- x = " << x << '\n'
<< "-> legendre(n, xd) = " << eve::legendre(n, xd) << '\n'
<< "-> legendre(3, xd) = " << eve::legendre(3, xd) << '\n'
<< "-> legendre(n, 0.5) = " << eve::legendre(n, 0.5) << '\n'
<< "-> legendre(n, x) = " << eve::legendre(n, x) << '\n';
double xs = 0.1;
std::cout << "---- scalar" << '\n'
<< "<- xs = " << xs << '\n'
<< "-> legendre(4, xs) = " << eve::legendre(4, xs) << '\n';
return 0;
}
constexpr callable_legendre_ legendre
Callable object evaluating the legendre functions.
Definition: legendre.hpp:89
Wrapper for SIMD registers.
Definition: wide.hpp:65