diff --git a/libcxx/include/cmath b/libcxx/include/cmath --- a/libcxx/include/cmath +++ b/libcxx/include/cmath @@ -296,6 +296,16 @@ float truncf(float x); long double truncl(long double x); +// C++17 Mathematical special functions [sf.cmath] + +floating_point assoc_laguerre (unsigned n, unsigned m, arithmetic x); +float assoc_laguerref(unsigned n, unsigned m, float x); +long double assoc_laguerrel(unsigned n, unsigned m, long double x); + +floating_point laguerre (unsigned n, arithmetic x); +float laguerref(unsigned n, float x); +long double laguerrel(unsigned n, long double x); + } // std */ @@ -303,6 +313,7 @@ #include <__config> #include #include +#include #if !defined(_LIBCPP_HAS_NO_PRAGMA_SYSTEM_HEADER) #pragma GCC system_header @@ -606,6 +617,124 @@ return isfinite(__lcpp_x); } +#if _LIBCPP_STD_VER > 14 + +template struct __promote_floating_point {}; +template <> struct __promote_floating_point { typedef double type; }; +template <> struct __promote_floating_point { typedef long double type; }; +template <> struct __promote_floating_point { typedef long double type; }; + +/// \brief Evaluates the generalized/associated Laguerre polynomial +/// \f$ L_{n}^{(\alpha)}(x) \f$. +/// +/// This method uses the three-term recurrence relation +/// \f[ +/// nL_{n}^{(\alpha)}(x) = (-x + 2n + \alpha - 1) L_{n-1}^{(\alpha)}(x) - +/// (n + \alpha - 1) L_{n-2}^{(\alpha)}(x) +/// \f] +/// in forward direction. It is not ill conditioned in both directions. +/// See Sec. 2.2 in Gil, A., Segura, J., & Temme, N. M. (2017). +/// Efficient computation of Laguerre polynomials. +/// Computer Physics Communications, 210, 124-131. +/// https://doi.org/10.1016/j.cpc.2016.09.002 +/// \tparam _Real one of float, double, or long double +/// \param [in] __n corresponds to \f$ n \f$. +/// \param [in] __alpha corresponds to \f$ \alpha \f$. +/// \param [in] __x corresponds to the argument \f$ x \f$. +/// \return the generalized Laguerre polynomial \f$ L_{n}^{(\alpha)}(x) \f$. +template +inline _LIBCPP_INLINE_VISIBILITY +_Real __libcpp_generalized_laguerre_recurrence(unsigned __n, _Real __alpha, + _Real __x) _NOEXCEPT { + if (__n == 0u) + return _Real(1); + + _Real __lm2 = _Real(1); + _Real __lm1 = _Real(1) + __alpha - __x; + const _Real __alpham1 = __alpha - _Real(1); + const _Real __alpham1x = __alpham1 - __x; + for (unsigned __i = 2; __i <= __n; ++__i) { + const _Real __li = ((2 * __i + __alpham1x) * __lm1 - (__i + __alpham1) * __lm2) / __i; + __lm2 = __lm1; + __lm1 = __li; + } + return __lm1; +} + +template +inline _LIBCPP_INLINE_VISIBILITY +_Real __libcpp_assoc_laguerre(unsigned __n, unsigned __m, _Real __x) { + if (std::isnan(__x)) + return std::numeric_limits<_Real>::quiet_NaN(); + + if (__x < _Real(0)) + __throw_domain_error("Argument of assoc_laguerre function is out of range"); + + return __libcpp_generalized_laguerre_recurrence(__n, _Real(__m), __x); +} + +template +inline _LIBCPP_INLINE_VISIBILITY +_Real __libcpp_laguerre(unsigned __n, _Real __x) { + if (std::isnan(__x)) + return std::numeric_limits<_Real>::quiet_NaN(); + + if (__x < _Real(0)) + __throw_domain_error("Argument of laguerre function is out of range"); + + return __libcpp_generalized_laguerre_recurrence(__n, _Real(0), __x); +} + + +template +inline _LIBCPP_INLINE_VISIBILITY +typename __lazy_enable_if::value, __promote<_A1>>::type +assoc_laguerre(unsigned __lcpp_n, unsigned __lcpp_m, _A1 __lcpp_x) { + typedef typename __promote<_A1>::type __result_type; + typedef typename __promote_floating_point<__result_type>::type __internal_type; + + return (__result_type) __libcpp_assoc_laguerre<__internal_type>( + __lcpp_n, __lcpp_m, (__internal_type) __lcpp_x); +} + +inline _LIBCPP_INLINE_VISIBILITY +float assoc_laguerref(unsigned __lcpp_n, unsigned __lcpp_m, + float __lcpp_x) { + return assoc_laguerre(__lcpp_n, __lcpp_m, __lcpp_x); +} + +inline _LIBCPP_INLINE_VISIBILITY +long double assoc_laguerrel(unsigned __lcpp_n, unsigned __lcpp_m, + long double __lcpp_x) { + return assoc_laguerre(__lcpp_n, __lcpp_m, __lcpp_x); +} + +template +inline _LIBCPP_INLINE_VISIBILITY +typename __lazy_enable_if::value, __promote<_A1>>::type +laguerre(unsigned __lcpp_n, _A1 __lcpp_x) { + typedef typename __promote<_A1>::type __result_type; + typedef typename __promote_floating_point<__result_type>::type __internal_type; + return (__result_type)__libcpp_laguerre<__internal_type>(__lcpp_n, + (__internal_type)__lcpp_x); +} + +inline _LIBCPP_INLINE_VISIBILITY +float laguerref(unsigned __lcpp_n, float __lcpp_x) { + return laguerre(__lcpp_n, __lcpp_x); +} + +inline _LIBCPP_INLINE_VISIBILITY +long double laguerrel(unsigned __lcpp_n, long double __lcpp_x) { + return laguerre(__lcpp_n, __lcpp_x); +} + + + +#endif + + + _LIBCPP_END_NAMESPACE_STD #endif // _LIBCPP_CMATH diff --git a/libcxx/test/std/numerics/c.math/cmath.pass.cpp b/libcxx/test/std/numerics/c.math/cmath.pass.cpp --- a/libcxx/test/std/numerics/c.math/cmath.pass.cpp +++ b/libcxx/test/std/numerics/c.math/cmath.pass.cpp @@ -100,6 +100,11 @@ Ambiguous tgamma(Ambiguous){ return Ambiguous(); } Ambiguous trunc(Ambiguous){ return Ambiguous(); } +#if _LIBCPP_STD_VER > 14 +Ambiguous assoc_laguerre(unsigned, unsigned, Ambiguous){ return Ambiguous(); } +Ambiguous laguerre(unsigned, Ambiguous){ return Ambiguous(); } +#endif + void test_abs() { static_assert((std::is_same::value), ""); @@ -1514,6 +1519,46 @@ assert(std::trunc(1) == 1); } +#if _LIBCPP_STD_VER > 14 + +void test_assoc_laguerre() +{ + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); +} + +void test_laguerre() +{ + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); + static_assert((std::is_same::value), ""); +} + +#endif + int main(int, char**) { test_abs(); @@ -1587,5 +1632,11 @@ test_tgamma(); test_trunc(); + +#if _LIBCPP_STD_VER > 14 + test_assoc_laguerre(); + test_laguerre(); +#endif + return 0; } diff --git a/libcxx/test/std/numerics/laguerre/assoc_laguerre_basic.pass.cpp b/libcxx/test/std/numerics/laguerre/assoc_laguerre_basic.pass.cpp new file mode 100644 --- /dev/null +++ b/libcxx/test/std/numerics/laguerre/assoc_laguerre_basic.pass.cpp @@ -0,0 +1,571 @@ +//===----------------------------------------------------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +// UNSUPPORTED: c++98, c++03, c++11, c++14 + +// + +// floating_point assoc_laguerre (unsigned n, unsigned m, arithmetic x); +// float assoc_laguerref(unsigned n, unsigned m, float x); +// long double assoc_laguerrel(unsigned n, unsigned m, long double x); + +#include +#include +#include +#include +#include +#include + +/// \brief std::assoc_laguerre(n,m,x) must silently return NaN if x is NaN. +/// +/// \tparam T one of float, double, or long double +template +void testNaNPropagation() { + for (unsigned n = 0; n < 128; ++n) { + for (unsigned m = 0; m < 128; ++m) { + const T QuietNaN = std::numeric_limits::quiet_NaN(); + const T ExpectNaN = std::assoc_laguerre(n, m, QuietNaN); + assert(std::isnan(ExpectNaN)); + } + } +} + +/// See Table 22.4 on p. 777 of Abramowitz, M., & Stegun, I. A. (1965). +/// Handbook of mathematical functions with formulas, graphs, and mathematical table. +/// In US Department of Commerce. National Bureau of Standards Applied Mathematics series 55. +template +void testAssocLaguerreAtZero(unsigned n, unsigned m, T RelTolerance, + T AbsTolerance) { + const long double LogExpected = + std::lgammal(n + m + 1u) - std::lgammal(n + 1u) - std::lgammal(m + 1u); + const T Expected = static_cast(std::exp(LogExpected)); + if (std::isinf(Expected)) + return; + + const T Result = std::assoc_laguerre(n, m, T(0)); + const T Tolerance = RelTolerance * std::abs(Expected) + AbsTolerance; + const T Error = std::abs(Result - Expected); + assert(Error <= Tolerance); +} + +template +void testAssocLaguerreAtZero(T RelTolerance, T AbsTolerance) { + for (unsigned n = 0u; n < 128u; ++n) { + for (unsigned m = 0u; m < 128u; ++m) { + testAssocLaguerreAtZero(n, m, RelTolerance, AbsTolerance); + } + } +} + +/// \brief Test if an analytic upper bound holds. +/// +/// For \f$ x \geq 0 \f$ and \f$ \alpha \geq 0 \f$, there is a simple +/// upper limit for all generalized Legendre polynomials, +/// \f[ +/// \left| L_{n}^{(\alpha)}(x) \right| \leq \frac{\Gamma\!\left( +/// \alpha + n + 1 \right)}{n!\Gamma\!\left( \alpha + 1 \right)} +/// \exp\!\left(\tfrac{1}{2}x\right). +/// \f] +/// See 22.14.13 on p. 786 of Abramowitz, M., & Stegun, I. A. (1965). +/// Handbook of mathematical functions with formulas, graphs, and mathematical table. +/// In US Department of Commerce. National Bureau of Standards Applied Mathematics series 55. +/// \param [in] x argument of assoc_laguerre +/// \param [in] RelTolerance relative tolerance for the upper limit of assoc_laguerre(x) +/// \param [in] AbsTolerance absolute tolerance for the upper limit of assoc_laguerre(x) +template +void testAssocLaguerreUpperLimit(unsigned n, unsigned m, T x, T RelTolerance, + T AbsTolerance) { + const T LogValueAtZero = + std::lgammal(n + m + 1u) - std::lgammal(n + 1u) - std::lgammal(m + 1u); + const T UpperLimit = static_cast(std::exp(0.5L * x + LogValueAtZero)); + if (std::isinf(UpperLimit)) + return; + const T UpperLimitWithTolerance = + UpperLimit + UpperLimit * RelTolerance + AbsTolerance; + const T Result = std::abs(std::laguerre(n, x)); + assert(Result <= UpperLimitWithTolerance); +} + +template +void testAssocLaguerreUpperLimit(T RelTolerance, T AbsTolerance) { + const T Samples[] = {T(0), T(0.5), T(1.0), T(10.0), T(100)}; + for (unsigned n = 0u; n < 128u; ++n) { + for (unsigned m = 0u; m < 128u; ++m) { + for (T x : Samples) + testAssocLaguerreUpperLimit(n, m, x, RelTolerance, AbsTolerance); + } + } +} + +/// The values are taken from: +/// Rabinowitz, P., & Weiss, G. (1959). +/// Tables of Abscissas and Weights for Numerical Evaluation of Integrals +/// of the Form \f$\int \limits_{0}^{\infty} e^{-x}x^n f (x)\, dx \f$. +/// Mathematical Tables and Other Aids to Computation, 285-294. +template +void testRoots(T RelTolerance, T AbsTolerance) { + std::map, std::vector > Table; + // clang-format off + Table[{4, 0}] = { + 3.22547689619392312e-1L, + 1.74576110115834658e0L, + 4.53662029692112798e0L, + 9.39507091230113313e0L + }; + + Table[{8, 0}] = { + 1.70279632305101000e-1L, + 9.03701776799379912e-1L, + 2.25108662986613069e0L, + 4.26670017028765879e0L, + 7.04590540239346570e0L, + 1.07585160101809952e1L, + 1.57406786412780046e1L, + 2.28631317368892641e1L + }; + + Table[{12, 0}] = { + 1.15722117358020675e-1L, + 6.11757484515130665e-1L, + 1.51261026977641879e0L, + 2.83375133774350723e0L, + 4.59922763941834848e0L, + 6.84452545311517735e0L, + 9.62131684245686704e0L, + 1.30060549933063477e1L, + 1.71168551874622557e1L, + 2.21510903793970057e1L, + 2.84879672509840003e1L, + 3.70991210444669203e1L + }; + + Table[{16, 0}] = { + 8.76494104789278403e-2L, + 4.62696328915080832e-1L, + 1.14105777483122686e0L, + 2.12928364509838062e0L, + 3.43708663389320665e0L, + 5.07801861454976791e0L, + 7.07033853504823413e0L, + 9.43831433639193878e0L, + 1.22142233688661587e1L, + 1.54415273687816171e1L, + 1.91801568567531349e1L, + 2.35159056939919085e1L, + 2.85787297428821404e1L, + 3.45833987022866258e1L, + 4.19404526476883326e1L, + 5.17011603395433184e1L + }; + + Table[{20, 0}] = { + 7.05398896919887534e-2L, + 3.72126818001611444e-1L, + 9.16582102483273565e-1L, + 1.70730653102834388e0L, + 2.74919925530943213e0L, + 4.04892531385088692e0L, + 5.61517497086161651e0L, + 7.45901745367106331e0L, + 9.59439286958109677e0L, + 1.20388025469643163e1L, + 1.48142934426307400e1L, + 1.79488955205193760e1L, + 2.14787882402850110e1L, + 2.54517027931869055e1L, + 2.99325546317006120e1L, + 3.50134342404790000e1L, + 4.08330570567285711e1L, + 4.76199940473465021e1L, + 5.58107957500638989e1L, + 6.65244165256157538e1L + }; + + Table[{24, 0}] = { + 5.90198521815079770e-2L, + 3.11239146198483727e-1L, + 7.66096905545936646e-1L, + 1.42559759080361309e0L, + 2.29256205863219029e0L, + 3.37077426420899772e0L, + 4.66508370346717079e0L, + 6.18153511873676541e0L, + 7.92753924717215218e0L, + 9.91209801507770602e0L, + 1.21461027117297656e1L, + 1.46427322895966743e1L, + 1.74179926465089787e1L, + 2.04914600826164247e1L, + 2.38873298481697332e1L, + 2.76359371743327174e1L, + 3.17760413523747233e1L, + 3.63584058016516217e1L, + 4.14517204848707670e1L, + 4.71531064451563230e1L, + 5.36085745446950698e1L, + 6.10585314472187616e1L, + 6.99622400351050304e1L, + 8.14982792339488854e1L + }; + + Table[{28, 0}] = { + 5.07346248498738876e-2L, + 2.67487268640741084e-1L, + 6.58136628354791519e-1L, + 1.22397180838490772e0L, + 1.96676761247377770e0L, + 2.88888332603032189e0L, + 3.99331165925011414e0L, + 5.28373606284344256e0L, + 6.76460340424350515e0L, + 8.44121632827132449e0L, + 1.03198504629932601e1L, + 1.24079034144606717e1L, + 1.47140851641357488e1L, + 1.72486634156080563e1L, + 2.00237833299517127e1L, + 2.30538901350302960e1L, + 2.63562973744013176e1L, + 2.99519668335961821e1L, + 3.38666055165844592e1L, + 3.81322544101946468e1L, + 4.27896723707725763e1L, + 4.78920716336227437e1L, + 5.35112979596642942e1L, + 5.97487960846412408e1L, + 6.67569772839064696e1L, + 7.47867781523391618e1L, + 8.43178371072270431e1L, + 9.65824206275273191e1L + }; + + Table[{32, 0}] = { + 4.44893658332670184e-2L, + 2.34526109519618537e-1L, + 5.76884629301886426e-1L, + 1.07244875381781763e0L, + 1.72240877644464544e0L, + 2.52833670642579488e0L, + 3.49221327302199449e0L, + 4.61645676974976739e0L, + 5.90395850417424395e0L, + 7.35812673318624111e0L, + 8.98294092421259610e0L, + 1.07830186325399721e1L, + 1.27636979867427251e1L, + 1.49311397555225573e1L, + 1.72924543367153148e1L, + 1.98558609403360547e1L, + 2.26308890131967745e1L, + 2.56286360224592478e1L, + 2.88621018163234747e1L, + 3.23466291539647370e1L, + 3.61004948057519738e1L, + 4.01457197715394415e1L, + 4.45092079957549380e1L, + 4.92243949873086392e1L, + 5.43337213333969073e1L, + 5.98925091621340182e1L, + 6.59753772879350528e1L, + 7.26876280906627086e1L, + 8.01874469779135231e1L, + 8.87353404178923987e1L, + 9.88295428682839726e1L, + 1.11751398097937695e2L + }; + + Table[{4, 1}] = { + 7.43291927981431435e-1L, + 2.57163500764627847e0L, + 5.73117875168909963e0L, + 1.09538943126831905e1L + }; + + Table[{8, 1}] = { + 4.09383573203185153e-1L, + 1.38496318480313988e0L, + 2.95625455616886207e0L, + 5.18194310104007138e0L, + 8.16170968814581733e0L, + 1.20700551268371548e1L, + 1.72497355261489875e1L, + 2.45859552436527819e1L + }; + + Table[{12, 1}] = { + 2.82858348239914112e-1L, + 9.52326041364613503e-1L, + 2.01649213857768073e0L, + 3.49235406977802504e0L, + 5.40549102001572284e0L, + 7.79281394041241498e0L, + 1.07073886889890908e1L, + 1.42271523637899807e1L, + 1.84719966342299991e1L, + 2.36417837524008953e1L, + 3.01200586261063443e1L, + 3.88892843760953185e1L + }; + + Table[{16, 1}] = { + 2.16140305239452255e-1L, + 7.26388243251803954e-1L, + 1.53359316037354132e0L, + 2.64497099861191096e0L, + 4.07097816088019065e0L, + 5.82585551510560455e0L, + 7.92850418530666709e0L, + 1.04038082899510393e1L, + 1.32846610707070382e1L, + 1.66151732168666126e1L, + 2.04560060200272169e1L, + 2.48938470253519108e1L, + 3.00598629202025761e1L, + 3.61706945436791780e1L, + 4.36403651841768370e1L, + 5.35291511602684204e1L + }; + + Table[{4, 2}] = { + 1.22676326350030207e0L, + 3.41250735869694597e0L, + 6.90269260585161340e0L, + 1.24580367719511386e1L + }; + + Table[{8, 2}] = { + 6.99330392297772446e-1L, + 1.89881649533754637e0L, + 3.67761476834163509e0L, + 6.09929454816086345e0L, + 9.26742581328238619e0L, + 1.33607382722601106e1L, + 1.87281386688430816e1L, + 2.62686410414766043e1L + }; + + Table[{12, 2}] = { + 4.90239109177454061e-1L, + 1.32377645578306233e0L, + 2.54213223523224643e0L, + 4.16451935324397066e0L, + 6.21800163077006985e0L, + 8.74078952492664456e0L, + 1.17872181154281359e1L, + 1.54365948600615973e1L, + 1.98104815433429476e1L, + 2.51111194768538770e1L, + 3.17262495017079404e1L, + 4.06488781934720538e1L + }; + + Table[{16, 2}] = { + 3.77613508344740734e-1L, + 1.01749195760256570e0L, + 1.94775802042424383e0L, + 3.17692724488986868e0L, + 4.71624006979179562e0L, + 6.58058826577491250e0L, + 8.78946527064707413e0L, + 1.13683230828333189e1L, + 1.43506267274370444e1L, + 1.77810957248416460e1L, + 2.17210847965713089e1L, + 2.62581386751110680e1L, + 3.15245960042758187e1L, + 3.77389210025289391e1L, + 4.53185461100898426e1L, + 5.53325835388358122e1L + }; + + Table[{4, 3}] = { + 1.75552164718549145e0L, + 4.26560586565682345e0L, + 8.05794068313800185e0L, + 1.39209318040196832e1L + }; + + Table[{8, 3}] = { + 1.02996168735087822e0L, + 2.43991423401425973e0L, + 4.41318676383851689e0L, + 7.01921044285466378e0L, + 1.03653585615970610e1L, + 1.46343284911153273e1L, + 2.01808472740399451e1L, + 2.79171925451893479e1L + }; + + Table[{12, 3}] = { + 7.31333453524153573e-1L, + 1.72200587771537790e0L, + 3.08711820886747629e0L, + 4.84907183643667852e0L, + 7.03654192274923744e0L, + 9.68900419917576858e0L, + 1.28619916279440914e1L, + 1.66361029842912300e1L, + 2.11344882654905403e1L, + 2.65616764973530448e1L, + 3.33094793553748032e1L, + 4.23811857710775980e1L + }; + + Table[{16, 3}] = { + 5.67443458991574145e-1L, + 1.33290773275989334e0L, + 2.38148248007005849e0L, + 3.72382664209342697e0L, + 5.37212395216187100e0L, + 7.34193662826135234e0L, + 9.65333213726122905e0L, + 1.23323014070182379e1L, + 1.54128500654077155e1L, + 1.89402755758602723e1L, + 2.29765957156018692e1L, + 2.76101814474260749e1L, + 3.29745092032584955e1L, + 3.92898232533641479e1L, + 4.69768962767103126e1L, + 5.71135140237534688e1L + }; + + Table[{4, 4}] = { + 2.31915524835569955e0L, + 5.12867199360117916e0L, + 9.20089134897243093e0L, + 1.53512814090706903e1L + }; + + Table[{8, 4}] = { + 1.39445874535841327e0L, + 3.00412262031589098e0L, + 5.16118127238140954e0L, + 7.94175644134279969e0L, + 1.14570496332457916e1L, + 1.58935676714650236e1L, + 2.16116099678724373e1L, + 2.95362536480182341e1L + }; + + Table[{12, 4}] = { + 1.00157108495383152e0L, + 2.14370613511225703e0L, + 3.64933321432680176e0L, + 5.54486813974946731e0L, + 7.86076259241096344e0L, + 1.06377654649934495e1L, + 1.39325683951155931e1L, + 1.78270063100842237e1L, + 2.24457538203894399e1L, + 2.79955513560437942e1L, + 3.48721707324949280e1L, + 4.40889427543252506e1L + }; + + Table[{16, 4}] = { + 7.82339164085635910e-1L, + 1.67007183671874001e0L, + 2.83292171947257737e0L, + 4.28441836819021891e0L, + 6.03786871096889479e0L, + 8.10954933413974021e0L, + 1.05201017048502218e1L, + 1.32960399123563184e1L, + 1.64718887086961996e1L, + 2.00934999730530394e1L, + 2.42235299428062549e1L, + 2.89511495806519383e1L, + 3.44109430408020216e1L, + 4.08248952713820481e1L, + 4.86170549734729113e1L, + 5.88737277583532394e1L + }; + + Table[{4, 5}] = { + 2.91078636936770048e0L, + 6.00000000000000000e0L, + 1.03341103467862578e1L, + 1.67551032838460417e1L + }; + + Table[{8, 5}] = { + 1.78786047157082745e0L, + 3.58823885764338781e0L, + 5.92009002612269490e0L, + 8.86684640118214593e0L, + 1.25435889444281568e1L, + 1.71405288833306908e1L, + 2.30233426297284257e1L, + 3.11295037859936706e1L + }; + + Table[{12, 5}] = { + 1.29748784000685508e0L, + 2.58620943634406719e0L, + 4.22697995386950259e0L, + 6.25087547037452539e0L, + 8.69028125595059872e0L, + 1.15872437894981120e1L, + 1.49995922575396278e1L, + 1.90103582323823253e1L, + 2.37456904925278601e1L, + 2.94144761288863084e1L, + 3.64163482053466644e1L, + 4.57744569372735530e1L + }; + + Table[{16, 5}] = { + 1.01975621006729237e0L, + 2.02686733891905119e0L, + 3.30047946506096635e0L, + 4.85757645713809794e0L, + 6.71275559761595590e0L, + 8.88305877614446853e0L, + 1.13897106617410895e1L, + 1.42597409258898484e1L, + 1.75281787111426743e1L, + 2.12414122811650920e1L, + 2.54627162726269377e1L, + 3.02820392701258409e1L, + 3.58350465381012185e1L, + 4.23454277633176611e1L, + 5.02404476674638389e1L, + 6.06147860634799663e1L + }; + + // clang-format on + + for (const auto& Entry : Table) { + const unsigned n = Entry.first.first; + const unsigned m = Entry.first.second; + for (long double Root : Entry.second) { + const T HalfWidthInterval = T(0.5) * (RelTolerance * Root + AbsTolerance); + const T LeftOfRoot = T(Root) - HalfWidthInterval; + const T RightOfRoot = T(Root) + HalfWidthInterval; + const T ValueLeft = std::assoc_laguerre(n, m, LeftOfRoot); + const T ValueRight = std::assoc_laguerre(n, m, RightOfRoot); + assert(std::signbit(ValueLeft) ^ std::signbit(ValueRight)); + } + } +} + +int main(int, char**) { + testNaNPropagation(); + testNaNPropagation(); + testNaNPropagation(); + + testAssocLaguerreAtZero(1e-6f, 0.0f); + testAssocLaguerreAtZero(1e-10, 0.0); + + testAssocLaguerreUpperLimit(1e-6f, 0.0f); + testAssocLaguerreUpperLimit(1e-14, 0.0); + + testRoots(1e-6f, 0.0f); + testRoots(1e-12, 0.0); +} diff --git a/libcxx/test/std/numerics/laguerre/assoc_laguerre_comparisons.pass.cpp b/libcxx/test/std/numerics/laguerre/assoc_laguerre_comparisons.pass.cpp new file mode 100644 --- /dev/null +++ b/libcxx/test/std/numerics/laguerre/assoc_laguerre_comparisons.pass.cpp @@ -0,0 +1,468 @@ +//===----------------------------------------------------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +// UNSUPPORTED: c++98, c++03, c++11, c++14 + +// + +// floating_point assoc_laguerre (unsigned n, unsigned m, arithmetic x); +// float assoc_laguerref(unsigned n, unsigned m, float x); +// long double assoc_laguerrel(unsigned n, unsigned m, long double x); + +#include +#include +#include +#include + +template +void compareWithTable( + const std::map, + std::map >& Table, + T RelTolerance, T AbsTolerance) { + for (const auto& Entry : Table) { + const unsigned l = Entry.first.first; + const unsigned m = Entry.first.second; + for (const auto& xy : Entry.second) { + const T ExpectedValue = static_cast(xy.second); + if (std::isinf(ExpectedValue)) + continue; + const T Value = std::assoc_laguerre(l, m, static_cast(xy.first)); + const T Tolerance = std::abs(ExpectedValue) * RelTolerance + AbsTolerance; + const T Error = std::abs(Value - ExpectedValue); + assert(Error <= Tolerance); + } + } +} + +template +void compareWithBoost(T RelTolerance, T AbsTolerance) { + // Generated with gcc 8.3.1 and boost 1.66.0 + std::map, std::map > + Table; + + Table[{10, 10}] = {{0.0L, 1.847560000000000000000e+05L}, + {0.5L, 1.150042533944161228305e+05L}, + {1.0L, 6.835330763365299824130e+04L}, + {1.5L, 3.819497840595790318474e+04L}, + {2.0L, 1.952449269841269840242e+04L}, + {2.5L, 8.617790029165816267032e+03L}, + {3.0L, 2.761132790178571428941e+03L}, + {3.5L, 2.647443095454463275404e+01L}, + {4.0L, -9.133717107583774248991e+02L}, + {4.5L, -9.351607177516392301331e+02L}, + {5.0L, -5.861136119378306876815e+02L}, + {5.5L, -1.796631517410951106584e+02L}, + {6.0L, 1.300342857142857142944e+02L}, + {6.5L, 2.913334329296293712663e+02L}, + {7.0L, 3.129948823302469135821e+02L}, + {7.5L, 2.340251775469098772520e+02L}, + {8.0L, 1.030916402116402116329e+02L}, + {8.5L, -3.454145314884353989981e+01L}, + {9.0L, -1.443173883928571428481e+02L}, + {9.5L, -2.054278306854086578043e+02L}}; + Table[{20, 10}] = {{0.0L, 3.004501500000000000000e+07L}, + {0.5L, 1.138930927314853252938e+07L}, + {1.0L, 3.651201251651593484667e+06L}, + {1.5L, 8.473698098486793293773e+05L}, + {2.0L, 4.566251493573249573643e+04L}, + {2.5L, -7.338348878522165050953e+04L}, + {3.0L, -3.208538228843214608688e+04L}, + {3.5L, 5.614549693109722527939e+03L}, + {4.0L, 1.414790961935669154048e+04L}, + {4.5L, 6.838363861624681725893e+03L}, + {5.0L, -1.795486420939170529643e+03L}, + {5.5L, -5.316437926431718982379e+03L}, + {6.0L, -3.887263891795776093119e+03L}, + {6.5L, -4.294692178238768064946e+02L}, + {7.0L, 2.300650817841901841776e+03L}, + {7.5L, 3.004218471787209729262e+03L}, + {8.0L, 1.857655658194353205381e+03L}, + {8.5L, -1.126184178317543087972e+02L}, + {9.0L, -1.763110242143649016411e+03L}, + {9.5L, -2.355059370043133450734e+03L}}; + Table[{20, 20}] = {{0.0L, 1.378465288200000000000e+11L}, + {0.5L, 8.468295361673511394858e+10L}, + {1.0L, 5.077912460175697552040e+10L}, + {1.5L, 2.960576939629002842866e+10L}, + {2.0L, 1.669670827385554017127e+10L}, + {2.5L, 9.043858989194316968322e+09L}, + {3.0L, 4.656130531018169694580e+09L}, + {3.5L, 2.241401912717946667224e+09L}, + {4.0L, 9.800663873313004742377e+08L}, + {4.5L, 3.659874821411927664303e+08L}, + {5.0L, 9.651825546133988800284e+07L}, + {5.5L, -2.223016973336849265706e+06L}, + {6.0L, -2.510182459882302286678e+07L}, + {6.5L, -2.030857764686113734570e+07L}, + {7.0L, -9.790386453177882433920e+06L}, + {7.5L, -1.710558937907162579336e+06L}, + {8.0L, 2.291830176433364631521e+06L}, + {8.5L, 3.127649252979586820175e+06L}, + {9.0L, 2.276440826891129583373e+06L}, + {9.5L, 9.499740595601278424738e+05L}}; + Table[{50, 10}] = {{0.0L, 7.539402756600000000000e+10L}, + {0.5L, 5.365943119336201372556e+09L}, + {1.0L, -6.058846804287141340683e+07L}, + {1.5L, -5.428935628098645252066e+06L}, + {2.0L, 6.491875483166804218854e+06L}, + {2.5L, -3.125723663735501904739e+06L}, + {3.0L, 3.215290936254617409418e+05L}, + {3.5L, 7.737700630158116541111e+05L}, + {4.0L, -3.911037482683599916982e+05L}, + {4.5L, -2.527733547889387023275e+05L}, + {5.0L, 2.294630650721308950750e+05L}, + {5.5L, 1.446393186292154792909e+05L}, + {6.0L, -1.282699408320068905027e+05L}, + {6.5L, -1.258756273224852412014e+05L}, + {7.0L, 5.353708357092831458601e+04L}, + {7.5L, 1.185268439074285104198e+05L}, + {8.0L, 1.536260133508163853566e+04L}, + {8.5L, -9.111765379659439557258e+04L}, + {9.0L, -7.288173794915159600549e+04L}, + {9.5L, 2.992603738587986626563e+04L}}; + Table[{50, 20}] = {{0.0L, 1.618846036626578760000e+17L}, + {0.5L, 4.682199741184386789062e+16L}, + {1.0L, 1.197682287075144469531e+16L}, + {1.5L, 2.562563860642318859131e+15L}, + {2.0L, 4.025512140617272661438e+14L}, + {2.5L, 2.558241994957465184402e+13L}, + {3.0L, -7.107135583922841755390e+12L}, + {3.5L, -1.527944352509347697258e+12L}, + {4.0L, 4.962109792862451581657e+11L}, + {4.5L, 2.207360501475066798627e+11L}, + {5.0L, -5.495990181724055286497e+10L}, + {5.5L, -5.395164965408203375712e+10L}, + {6.0L, 2.160824469618739394005e+09L}, + {6.5L, 1.625278916791839989182e+10L}, + {7.0L, 4.561836324614418234210e+09L}, + {7.5L, -4.311154213993686890230e+09L}, + {8.0L, -3.772710079809093437856e+09L}, + {8.5L, 6.289714539106999114301e+07L}, + {9.0L, 1.862043401563646103023e+09L}, + {9.5L, 1.093508422577525381697e+09L}}; + Table[{50, 50}] = {{0.0L, 1.008913445455641933873e+29L}, + {0.5L, 6.150335228601113623285e+28L}, + {1.0L, 3.712351462622653283938e+28L}, + {1.5L, 2.217375806127585064487e+28L}, + {2.0L, 1.309702308652025846160e+28L}, + {2.5L, 7.643899665374721918712e+27L}, + {3.0L, 4.404413646826269799461e+27L}, + {3.5L, 2.502994762342186406679e+27L}, + {4.0L, 1.401305429217811545136e+27L}, + {4.5L, 7.718403057066373194334e+26L}, + {5.0L, 4.176012189047808123250e+26L}, + {5.5L, 2.215242482188571356776e+26L}, + {6.0L, 1.149522800633857380555e+26L}, + {6.5L, 5.818720365519067235851e+25L}, + {7.0L, 2.862862686765511838217e+25L}, + {7.5L, 1.362737172361160253296e+25L}, + {8.0L, 6.236142696693560635818e+24L}, + {8.5L, 2.718917011811881114665e+24L}, + {9.0L, 1.113959548616688995271e+24L}, + {9.5L, 4.190270485617718532833e+23L}}; + Table[{100, 10}] = {{0.0L, 4.689763662398100000000e+13L}, + {0.5L, -4.294878427688370894874e+09L}, + {1.0L, 3.944092408530536249978e+09L}, + {1.5L, -2.013797229225323652645e+08L}, + {2.0L, -2.157299751648581376503e+07L}, + {2.5L, 1.459430547397130502941e+07L}, + {3.0L, -1.016058418748330159019e+06L}, + {3.5L, -6.908288622716346500511e+06L}, + {4.0L, 9.787558305163868595628e+06L}, + {4.5L, -8.176816207612928389608e+06L}, + {5.0L, 3.377298430833765111174e+06L}, + {5.5L, 1.735890291779104317243e+06L}, + {6.0L, -3.807192070466917534077e+06L}, + {6.5L, 1.734983554763113708077e+06L}, + {7.0L, 1.706476023943294066271e+06L}, + {7.5L, -2.355218427716766161211e+06L}, + {8.0L, -3.058164836881507359863e+05L}, + {8.5L, 2.149540815590451997195e+06L}, + {9.0L, -3.345374513274238412919e+05L}, + {9.5L, -1.920373498643499776108e+06L}}; + Table[{100, 20}] = {{0.0L, 2.946222729117663572173e+22L}, + {0.5L, 2.261308869389457750144e+21L}, + {1.0L, 9.475799170461032115200e+19L}, + {1.5L, -3.637696286243510428125e+17L}, + {2.0L, -2.432420166605434042578e+16L}, + {2.5L, 9.128749654970749480469e+15L}, + {3.0L, -2.745178778987279971191e+15L}, + {3.5L, 7.321257828016063543701e+14L}, + {4.0L, -1.137164349493819836121e+14L}, + {4.5L, -3.827471412315585189056e+13L}, + {5.0L, 4.041283789325254424667e+13L}, + {5.5L, -1.225276122801313934898e+13L}, + {6.0L, -4.723168970277397794247e+12L}, + {6.5L, 5.788209026781585838795e+12L}, + {7.0L, -6.214739115929799389839e+11L}, + {7.5L, -2.016398224166447072029e+12L}, + {8.0L, 8.280834859221748547554e+11L}, + {8.5L, 7.409824494487057602406e+11L}, + {9.0L, -5.314481262655516004860e+11L}, + {9.5L, -3.340793205142770996094e+11L}}; + Table[{100, 50}] = {{0.0L, 2.012866090973193230660e+40L}, + {0.5L, 7.443534596401896623943e+39L}, + {1.0L, 2.669138788507098282845e+39L}, + {1.5L, 9.249176234411826075623e+38L}, + {2.0L, 3.084251014932004951249e+38L}, + {2.5L, 9.845130400028869439863e+37L}, + {3.0L, 2.987893170071304389611e+37L}, + {3.5L, 8.543446586522976376137e+36L}, + {4.0L, 2.272446170635370974820e+36L}, + {4.5L, 5.516531148637918796375e+35L}, + {5.0L, 1.184363719341424931225e+35L}, + {5.5L, 2.115939411877618798396e+34L}, + {6.0L, 2.679064308539579242504e+33L}, + {6.5L, 6.930372893313020755552e+31L}, + {7.0L, -7.006151656594741370244e+31L}, + {7.5L, -1.654854053550525560767e+31L}, + {8.0L, 5.124436639236123969254e+29L}, + {8.5L, 1.212547004925553862607e+30L}, + {9.0L, 2.130085958501174472883e+29L}, + {9.5L, -6.983885241918135655551e+28L}}; + Table[{100, 100}] = {{0.0L, 9.054851465610328115415e+58L}, + {0.5L, 5.505859466697853237783e+58L}, + {1.0L, 3.331262988031657886627e+58L}, + {1.5L, 2.005242577321074531189e+58L}, + {2.0L, 1.200686309789629298134e+58L}, + {2.5L, 7.150282668151546652635e+57L}, + {3.0L, 4.234184715733205324058e+57L}, + {3.5L, 2.492794079123697664219e+57L}, + {4.0L, 1.458763178926480349637e+57L}, + {4.5L, 8.483418739565801742550e+56L}, + {5.0L, 4.901672141399097796573e+56L}, + {5.5L, 2.813178471492599513520e+56L}, + {6.0L, 1.603299095004078548328e+56L}, + {6.5L, 9.071348729442847391696e+55L}, + {7.0L, 5.093722441315075924182e+55L}, + {7.5L, 2.837656817120185206977e+55L}, + {8.0L, 1.567791758256770327412e+55L}, + {8.5L, 8.587159575060838647900e+54L}, + {9.0L, 4.660753669629245405396e+54L}, + {9.5L, 2.505548374372461101626e+54L}}; + + compareWithTable(Table, RelTolerance, AbsTolerance); +} + +template +void compareWithGCC(T RelTolerance, T AbsTolerance) { + // Generated with gcc 8.3.1 + std::map, std::map > + Table; + + // Generated with gcc 8.3.1 + Table[{10, 10}] = {{0.0L, 1.847559999999999999858e+05L}, + {0.5L, 1.150042533944161228305e+05L}, + {1.0L, 6.835330763365299823420e+04L}, + {1.5L, 3.819497840595790315632e+04L}, + {2.0L, 1.952449269841269841130e+04L}, + {2.5L, 8.617790029165816266143e+03L}, + {3.0L, 2.761132790178571431605e+03L}, + {3.5L, 2.647443095454463224403e+01L}, + {4.0L, -9.133717107583774249546e+02L}, + {4.5L, -9.351607177516392298000e+02L}, + {5.0L, -5.861136119378306877925e+02L}, + {5.5L, -1.796631517410951106584e+02L}, + {6.0L, 1.300342857142857142666e+02L}, + {6.5L, 2.913334329296293712663e+02L}, + {7.0L, 3.129948823302469135266e+02L}, + {7.5L, 2.340251775469098772520e+02L}, + {8.0L, 1.030916402116402116051e+02L}, + {8.5L, -3.454145314884353987206e+01L}, + {9.0L, -1.443173883928571428481e+02L}, + {9.5L, -2.054278306854086578181e+02L}}; + Table[{20, 10}] = {{0.0L, 3.004501499999999999454e+07L}, + {0.5L, 1.138930927314853253574e+07L}, + {1.0L, 3.651201251651593483757e+06L}, + {1.5L, 8.473698098486793282405e+05L}, + {2.0L, 4.566251493573249561564e+04L}, + {2.5L, -7.338348878522165048111e+04L}, + {3.0L, -3.208538228843214611175e+04L}, + {3.5L, 5.614549693109722516837e+03L}, + {4.0L, 1.414790961935669154048e+04L}, + {4.5L, 6.838363861624681717899e+03L}, + {5.0L, -1.795486420939170525424e+03L}, + {5.5L, -5.316437926431718985931e+03L}, + {6.0L, -3.887263891795776096671e+03L}, + {6.5L, -4.294692178238768076604e+02L}, + {7.0L, 2.300650817841901841554e+03L}, + {7.5L, 3.004218471787209729262e+03L}, + {8.0L, 1.857655658194353205159e+03L}, + {8.5L, -1.126184178317543080894e+02L}, + {9.0L, -1.763110242143649016189e+03L}, + {9.5L, -2.355059370043133450956e+03L}}; + Table[{20, 20}] = {{0.0L, 1.378465288200000000000e+11L}, + {0.5L, 8.468295361673511391878e+10L}, + {1.0L, 5.077912460175697550923e+10L}, + {1.5L, 2.960576939629002842307e+10L}, + {2.0L, 1.669670827385554017220e+10L}, + {2.5L, 9.043858989194316969253e+09L}, + {3.0L, 4.656130531018169695046e+09L}, + {3.5L, 2.241401912717946666060e+09L}, + {4.0L, 9.800663873313004756928e+08L}, + {4.5L, 3.659874821411927664594e+08L}, + {5.0L, 9.651825546133988801739e+07L}, + {5.5L, -2.223016973336849307088e+06L}, + {6.0L, -2.510182459882302287770e+07L}, + {6.5L, -2.030857764686113736207e+07L}, + {7.0L, -9.790386453177882432101e+06L}, + {7.5L, -1.710558937907162581723e+06L}, + {8.0L, 2.291830176433364629020e+06L}, + {8.5L, 3.127649252979586819947e+06L}, + {9.0L, 2.276440826891129583601e+06L}, + {9.5L, 9.499740595601278420190e+05L}}; + Table[{50, 10}] = {{0.0L, 7.539402756600000000000e+10L}, + {0.5L, 5.365943119336201387458e+09L}, + {1.0L, -6.058846804287141332679e+07L}, + {1.5L, -5.428935628098645301179e+06L}, + {2.0L, 6.491875483166804202028e+06L}, + {2.5L, -3.125723663735501907013e+06L}, + {3.0L, 3.215290936254617415102e+05L}, + {3.5L, 7.737700630158116542816e+05L}, + {4.0L, -3.911037482683599921813e+05L}, + {4.5L, -2.527733547889387019154e+05L}, + {5.0L, 2.294630650721308950040e+05L}, + {5.5L, 1.446393186292154793762e+05L}, + {6.0L, -1.282699408320068904743e+05L}, + {6.5L, -1.258756273224852412369e+05L}, + {7.0L, 5.353708357092831457180e+04L}, + {7.5L, 1.185268439074285104695e+05L}, + {8.0L, 1.536260133508163865201e+04L}, + {8.5L, -9.111765379659439550153e+04L}, + {9.0L, -7.288173794915159593444e+04L}, + {9.5L, 2.992603738587986625852e+04L}}; + Table[{50, 20}] = {{0.0L, 1.618846036626578760000e+17L}, + {0.5L, 4.682199741184386786719e+16L}, + {1.0L, 1.197682287075144469336e+16L}, + {1.5L, 2.562563860642318861816e+15L}, + {2.0L, 4.025512140617272669678e+14L}, + {2.5L, 2.558241994957465185547e+13L}, + {3.0L, -7.107135583922841745853e+12L}, + {3.5L, -1.527944352509347693205e+12L}, + {4.0L, 4.962109792862451580167e+11L}, + {4.5L, 2.207360501475066800117e+11L}, + {5.0L, -5.495990181724055255949e+10L}, + {5.5L, -5.395164965408203381300e+10L}, + {6.0L, 2.160824469618739462458e+09L}, + {6.5L, 1.625278916791839991789e+10L}, + {7.0L, 4.561836324614418232813e+09L}, + {7.5L, -4.311154213993686886970e+09L}, + {8.0L, -3.772710079809093443211e+09L}, + {8.5L, 6.289714539106999291107e+07L}, + {9.0L, 1.862043401563646102557e+09L}, + {9.5L, 1.093508422577525383327e+09L}}; + Table[{50, 50}] = {{0.0L, 1.008913445455641933615e+29L}, + {0.5L, 6.150335228601113628010e+28L}, + {1.0L, 3.712351462622653281791e+28L}, + {1.5L, 2.217375806127585061480e+28L}, + {2.0L, 1.309702308652025845945e+28L}, + {2.5L, 7.643899665374721938039e+27L}, + {3.0L, 4.404413646826269802683e+27L}, + {3.5L, 2.502994762342186410974e+27L}, + {4.0L, 1.401305429217811544465e+27L}, + {4.5L, 7.718403057066373177556e+26L}, + {5.0L, 4.176012189047808123585e+26L}, + {5.5L, 2.215242482188571357447e+26L}, + {6.0L, 1.149522800633857379884e+26L}, + {6.5L, 5.818720365519067237948e+25L}, + {7.0L, 2.862862686765511838636e+25L}, + {7.5L, 1.362737172361160253191e+25L}, + {8.0L, 6.236142696693560632148e+24L}, + {8.5L, 2.718917011811881114927e+24L}, + {9.0L, 1.113959548616688994288e+24L}, + {9.5L, 4.190270485617718538404e+23L}}; + Table[{100, 10}] = {{0.0L, 4.689763662398099999619e+13L}, + {0.5L, -4.294878427688372131437e+09L}, + {1.0L, 3.944092408530536258360e+09L}, + {1.5L, -2.013797229225323626306e+08L}, + {2.0L, -2.157299751648581400877e+07L}, + {2.5L, 1.459430547397130521313e+07L}, + {3.0L, -1.016058418748330104791e+06L}, + {3.5L, -6.908288622716346549169e+06L}, + {4.0L, 9.787558305163868584714e+06L}, + {4.5L, -8.176816207612928385061e+06L}, + {5.0L, 3.377298430833765122770e+06L}, + {5.5L, 1.735890291779104324632e+06L}, + {6.0L, -3.807192070466917537033e+06L}, + {6.5L, 1.734983554763113707850e+06L}, + {7.0L, 1.706476023943294063884e+06L}, + {7.5L, -2.355218427716766165304e+06L}, + {8.0L, -3.058164836881507324051e+05L}, + {8.5L, 2.149540815590451998105e+06L}, + {9.0L, -3.345374513274238433382e+05L}, + {9.5L, -1.920373498643499780655e+06L}}; + Table[{100, 20}] = {{0.0L, 2.946222729117663572378e+22L}, + {0.5L, 2.261308869389457764992e+21L}, + {1.0L, 9.475799170461032092800e+19L}, + {1.5L, -3.637696286243510460000e+17L}, + {2.0L, -2.432420166605433971875e+16L}, + {2.5L, 9.128749654970749492188e+15L}, + {3.0L, -2.745178778987279970215e+15L}, + {3.5L, 7.321257828016063547363e+14L}, + {4.0L, -1.137164349493819837799e+14L}, + {4.5L, -3.827471412315585184479e+13L}, + {5.0L, 4.041283789325254420090e+13L}, + {5.5L, -1.225276122801313933945e+13L}, + {6.0L, -4.723168970277397769451e+12L}, + {6.5L, 5.788209026781585855007e+12L}, + {7.0L, -6.214739115929799386263e+11L}, + {7.5L, -2.016398224166447076082e+12L}, + {8.0L, 8.280834859221748533249e+11L}, + {8.5L, 7.409824494487057593465e+11L}, + {9.0L, -5.314481262655516003966e+11L}, + {9.5L, -3.340793205142771008909e+11L}}; + Table[{100, 50}] = {{0.0L, 2.012866090973193229952e+40L}, + {0.5L, 7.443534596401896641062e+39L}, + {1.0L, 2.669138788507098282403e+39L}, + {1.5L, 9.249176234411826049060e+38L}, + {2.0L, 3.084251014932004955122e+38L}, + {2.5L, 9.845130400028869442630e+37L}, + {3.0L, 2.987893170071304389611e+37L}, + {3.5L, 8.543446586522976381325e+36L}, + {4.0L, 2.272446170635370969488e+36L}, + {4.5L, 5.516531148637918790250e+35L}, + {5.0L, 1.184363719341424931765e+35L}, + {5.5L, 2.115939411877618808980e+34L}, + {6.0L, 2.679064308539579214357e+33L}, + {6.5L, 6.930372893313020587546e+31L}, + {7.0L, -7.006151656594741382119e+31L}, + {7.5L, -1.654854053550525565935e+31L}, + {8.0L, 5.124436639236123841092e+29L}, + {8.5L, 1.212547004925553863363e+30L}, + {9.0L, 2.130085958501174466698e+29L}, + {9.5L, -6.983885241918135725559e+28L}}; + Table[{100, 100}] = {{0.0L, 9.054851465610328117048e+58L}, + {0.5L, 5.505859466697853241594e+58L}, + {1.0L, 3.331262988031657892344e+58L}, + {1.5L, 2.005242577321074531325e+58L}, + {2.0L, 1.200686309789629297794e+58L}, + {2.5L, 7.150282668151546656718e+57L}, + {3.0L, 4.234184715733205334267e+57L}, + {3.5L, 2.492794079123697660987e+57L}, + {4.0L, 1.458763178926480348871e+57L}, + {4.5L, 8.483418739565801730640e+56L}, + {5.0L, 4.901672141399097794021e+56L}, + {5.5L, 2.813178471492599516284e+56L}, + {6.0L, 1.603299095004078549498e+56L}, + {6.5L, 9.071348729442847404988e+55L}, + {7.0L, 5.093722441315075925245e+55L}, + {7.5L, 2.837656817120185205382e+55L}, + {8.0L, 1.567791758256770324753e+55L}, + {8.5L, 8.587159575060838631949e+54L}, + {9.0L, 4.660753669629245418356e+54L}, + {9.5L, 2.505548374372461097140e+54L}}; + compareWithTable(Table, RelTolerance, AbsTolerance); +} + +int main(int, char**) { + compareWithBoost(1e-6f, 0.0f); + compareWithBoost(1e-12, 0.0); + + compareWithGCC(1e-6f, 0.0f); + compareWithGCC(1e-12, 0.0); +} diff --git a/libcxx/test/std/numerics/laguerre/assoc_laguerre_exceptions.pass.cpp b/libcxx/test/std/numerics/laguerre/assoc_laguerre_exceptions.pass.cpp new file mode 100644 --- /dev/null +++ b/libcxx/test/std/numerics/laguerre/assoc_laguerre_exceptions.pass.cpp @@ -0,0 +1,49 @@ +//===----------------------------------------------------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +// UNSUPPORTED: c++98, c++03, c++11, c++14, libcpp-no-exceptions + +// + +// floating_point assoc_laguerre (unsigned n, unsigned m, arithmetic x); +// float assoc_laguerref(unsigned n, unsigned m, float x); +// long double assoc_laguerrel(unsigned n, unsigned m, long double x); + +#include +#include +#include +#include + +/// \brief std::assoc_laguerre(l, m, x) throws if x < 0 and not NaN. +/// +/// \tparam T one of float, double, or long double +template +void testInvalidArguments() { + const T InvalidArguments[] = {T(-5), -std::numeric_limits::min()}; + for (unsigned l = 0; l < 128; ++l) { + for (unsigned m = 0; m < 128; ++m) { + for (T x : InvalidArguments) { + bool ThrowsDomainError = false; + try { + const T Value = std::assoc_laguerre(l, m, x); + (void)Value; + } catch (const std::domain_error&) { + ThrowsDomainError = true; + } catch (...) { + } + assert(ThrowsDomainError); + } + } + } +} + +int main(int, char**) { + testInvalidArguments(); + testInvalidArguments(); + testInvalidArguments(); +} diff --git a/libcxx/test/std/numerics/laguerre/laguerre_basic.pass.cpp b/libcxx/test/std/numerics/laguerre/laguerre_basic.pass.cpp new file mode 100644 --- /dev/null +++ b/libcxx/test/std/numerics/laguerre/laguerre_basic.pass.cpp @@ -0,0 +1,169 @@ +//===----------------------------------------------------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +// UNSUPPORTED: c++98, c++03, c++11, c++14 + +// + +// floating_point laguerre (unsigned n, arithmetic x); +// float laguerref(unsigned n, float x); +// long double laguerrel(unsigned n, long double x); + +#include +#include +#include +#include + +template +Real evaluatePolynomial(Iterator CoefficientsBegin, Iterator CoefficientsEnd, + Real x) { + Real Result = Real(0); + for (Iterator Iter = CoefficientsBegin; Iter != CoefficientsEnd; ++Iter) { + Result = Result * x + *Iter; + } + return Result; +} + +template +Real evaluatePolynomial(const Container& Coefficients, Real x) { + return evaluatePolynomial(std::begin(Coefficients), std::end(Coefficients), + x); +} + +/// \brief std::assoc_laguerre(n,m,x) must silently return NaN if x is NaN. +/// +/// \tparam T one of float, double, or long double +template +void testNaNPropagation() { + for (unsigned n = 0; n < 128; ++n) { + const T QuietNaN = std::numeric_limits::quiet_NaN(); + const T ExpectNaN = std::laguerre(n, QuietNaN); + assert(std::isnan(ExpectNaN)); + } +} + +/// \brief Generates values of the Laguerre polynomials from +/// analytic expressions for n = 0..6. +/// +/// See Table 13.1 on p. 819 of: +/// Arfken, G. B., & Weber, H. J. (2005). +/// Mathematical Methods for Physicists. 6th edn (Elsevier Academic Press) +/// \tparam T one of float, double, or long double +/// \tparam TestFunction a callable with arguments (unsigned,T,T) +/// \param [in] x argument of laguerre +/// \param [in] ValueSink will be invoked with the arguments of +/// laguerre and its expected outcome, i.e., (n,x,expected). +template +void compareWithAnalyticExpressions(T x, TestFunction&& ValueSink) { + std::vector > VectorOfCoeffients; + VectorOfCoeffients.push_back({1.0L}); + VectorOfCoeffients.push_back({-1.0L, 1.0L}); + VectorOfCoeffients.push_back({1.0L, -4.0L, 2.0L}); + VectorOfCoeffients.push_back({-1.0L, 9.0L, -18.0L, 6.0L}); + VectorOfCoeffients.push_back({1.0L, -16.0L, 72.0L, -96.0L, 24.0L}); + VectorOfCoeffients.push_back( + {-1.0L, 25.0L, -200.0L, 600.0L, -600.0L, 120.0L}); + VectorOfCoeffients.push_back( + {1.0L, -36.0L, 450.0L, -2400.0L, 5400.0L, -4320.0L, 720.0L}); + + long double Scale = 1.0L; + for (size_t i = 2; i < VectorOfCoeffients.size(); ++i) { + Scale /= static_cast(i); + for (long double& y : VectorOfCoeffients[i]) + y *= Scale; + } + + for (const auto& Coefficients : VectorOfCoeffients) { + assert(!Coefficients.empty()); + const unsigned int n = static_cast(Coefficients.size()) - 1; + const T ExpectedValue = static_cast( + evaluatePolynomial(Coefficients, static_cast(x))); + ValueSink(n, x, ExpectedValue); + } +} + +template +void compareWithAnalyticExpressions(T RelTolerance, T AbsTolerance) { + auto testFunction = [RelTolerance, AbsTolerance](unsigned n, T x, + T ExpectedValue) { + const T Result = std::laguerre(n, x); + const T Tolerance = std::abs(ExpectedValue) * RelTolerance + AbsTolerance; + const T Error = std::abs(Result - ExpectedValue); + assert(Error <= Tolerance); + }; + + const T Samples[] = {T(0), T(0.1), T(0.2), T(0.5), T(0.9), T(1.0)}; + for (T x : Samples) + compareWithAnalyticExpressions(x, testFunction); +} + +/// \brief Test if an analytic upper bound holds. +/// +/// For \f$ x \geq 0 \f$, there is a simple upper limit for all Legendre +/// polynomials, +/// \f[ +/// \left| L_{n}(x) \right| \leq \exp \left(\tfrac{1}{2}x\right). +/// \f] +/// See 22.14.12 on p. 786 of Abramowitz, M., & Stegun, I. A. (1965). +/// Handbook of mathematical functions with formulas, graphs, and mathematical table. +/// In US Department of Commerce. National Bureau of Standards Applied Mathematics series 55. +/// \param [in] x argument of laguerre, must be >= 0. +/// \param [in] RelTolerance relative tolerance for the upper limit of laguerre(x) +/// \param [in] AbsTolerance absolute tolerance for the upper limit of laguerre(x) +template +void testLaguerreUpperLimit(unsigned n, T x, T RelTolerance, T AbsTolerance) { + const T UpperLimit = std::exp(x / T(2)); + if (std::isinf(UpperLimit)) + return; + + const T UpperLimitWithTolerance = + UpperLimit + UpperLimit * RelTolerance + AbsTolerance; + const T Result = std::abs(std::laguerre(n, x)); + assert(Result <= UpperLimitWithTolerance); +} + +template +void testLaguerreUpperLimit(T RelTolerance, T AbsTolerance) { + const T LargeX = T(2) * std::log(std::numeric_limits::max()); + const T Samples[] = {T(0), T(0.5), T(1.0), T(10.0), T(100), T(0.99) * LargeX}; + for (unsigned n = 0u; n < 128u; ++n) { + for (T x : Samples) + testLaguerreUpperLimit(n, x, RelTolerance, AbsTolerance); + } +} + +/// \brief Tests the value of laguerre at x=0 +/// +/// We expect +/// \f[ +/// L_{n}(x) = 1 \quad \forall n \in \mathbb{N}, +/// \f] +/// which follows immediately from the recurrence relation. +template +void testValueAtZero(T RelTolerance, T AbsTolerance) { + for (unsigned n = 0u; n < 128u; ++n) { + const T Result = std::laguerre(n, T(0)); + const T Error = std::abs(Result - T(1)); + assert(Error <= RelTolerance + AbsTolerance); + } +} + +int main(int, char**) { + testNaNPropagation(); + testNaNPropagation(); + testNaNPropagation(); + + compareWithAnalyticExpressions(1e-6f, 0.0f); + compareWithAnalyticExpressions(1e-12, 0.0); + + testLaguerreUpperLimit(1e-6f, 0.0f); + testLaguerreUpperLimit(1e-14, 0.0); + + testValueAtZero(1e-6f, 0.0); + testValueAtZero(1e-14, 0.0); +} diff --git a/libcxx/test/std/numerics/laguerre/laguerre_comparisons.pass.cpp b/libcxx/test/std/numerics/laguerre/laguerre_comparisons.pass.cpp new file mode 100644 --- /dev/null +++ b/libcxx/test/std/numerics/laguerre/laguerre_comparisons.pass.cpp @@ -0,0 +1,221 @@ +//===----------------------------------------------------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +// UNSUPPORTED: c++98, c++03, c++11, c++14 + +// + +// floating_point laguerre (unsigned n, arithmetic x); +// float laguerref(unsigned n, float x); +// long double laguerrel(unsigned n, long double x); + +#include +#include +#include + +template +void compareWithTable( + const std::map >& Table, + T RelTolerance, T AbsTolerance) { + for (const auto& Entry : Table) { + const unsigned n = Entry.first; + for (const auto& xy : Entry.second) { + const T Value = std::laguerre(n, xy.first); + const T Expected = xy.second; + const T Error = std::abs(Expected - Value); + const T Tolerance = std::abs(Expected) * RelTolerance + AbsTolerance; + assert(Error <= Tolerance); + } + } +} + +template +void compareWithGCC(T RelTolerance, T AbsTolerance) { + // Generated with gcc 8.3.1 + std::map > Table; + + Table[10] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, -3.893744141378520447379e-01L}, + {1.0L, 4.189459325396825395394e-01L}, + {1.5L, 4.180177743094308035138e-01L}, + {2.0L, -3.090652557319223986312e-01L}, + {2.5L, -8.802526766660982968137e-01L}, + {3.0L, -7.000223214285714282861e-01L}, + {3.5L, 2.211135054223331404218e-01L}, + {4.0L, 1.379259259259259259400e+00L}, + {4.5L, 2.078780713762555804134e+00L}, + {5.0L, 1.756276179453262787115e+00L}, + {5.5L, 2.235215177990141372649e-01L}, + {6.0L, -2.222857142857142857166e+00L}, + {6.5L, -4.837352781211678103002e+00L}, + {7.0L, -6.589936342592592592466e+00L}, + {7.5L, -6.421207427978515625434e+00L}, + {8.0L, -3.538694885361552027944e+00L}, + {8.5L, 2.300054675712787284787e+00L}, + {9.0L, 1.057033482142857142971e+01L}, + {9.5L, 1.987481935988868565236e+01L}}; + Table[20] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, 3.120174870800154148060e-01L}, + {1.0L, -1.642588118277929596923e-01L}, + {1.5L, -3.415329100484649494903e-01L}, + {2.0L, 5.007215675600304566318e-01L}, + {2.5L, 4.695899854819089235349e-01L}, + {3.0L, -5.575093252171908237303e-01L}, + {3.5L, -1.070298248145694926204e+00L}, + {4.0L, -1.507914281021107644285e-01L}, + {4.5L, 1.437239150257817378586e+00L}, + {5.0L, 2.020225744476913576700e+00L}, + {5.5L, 5.842074183986535973554e-01L}, + {6.0L, -2.226540971456672472891e+00L}, + {6.5L, -4.371331381805792055777e+00L}, + {7.0L, -3.714934120411960052216e+00L}, + {7.5L, 4.510955524973794079566e-01L}, + {8.0L, 6.568596673285710330764e+00L}, + {8.5L, 1.118547863037195619872e+01L}, + {9.0L, 1.048071357226878702867e+01L}, + {9.5L, 2.334272934394408317502e+00L}}; + Table[50] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, -3.181388060269979059503e-01L}, + {1.0L, 2.329893629018530612886e-01L}, + {1.5L, -2.568364397519475615444e-01L}, + {2.0L, 4.416561042694050874145e-01L}, + {2.5L, -5.522468012177843677381e-01L}, + {3.0L, 1.463786536365346789817e-01L}, + {3.5L, 7.431381187021829679918e-01L}, + {4.0L, -8.261038021062498654763e-01L}, + {4.5L, -7.795068145562651417205e-01L}, + {5.0L, 1.473525881943055793988e+00L}, + {5.5L, 1.221520285433004828976e+00L}, + {6.0L, -2.083161364452417261411e+00L}, + {6.5L, -2.649222481840662724245e+00L}, + {7.0L, 1.961343956455813747775e+00L}, + {7.5L, 5.387138713354476469639e+00L}, + {8.0L, 8.689419213752169954268e-01L}, + {8.5L, -7.772031120599373795778e+00L}, + {9.0L, -8.576161635248512745203e+00L}, + {9.5L, 3.641234604888676069934e+00L}}; + Table[100] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, 1.868226036769227864408e-01L}, + {1.0L, 2.706659209933101103019e-01L}, + {1.5L, 6.389168619635640683941e-02L}, + {2.0L, -3.031283692749120643544e-01L}, + {2.5L, 4.257949400323675469694e-01L}, + {3.0L, -4.794265341305422037854e-01L}, + {3.5L, 3.853183234834782489313e-01L}, + {4.0L, 2.339473046150550509105e-02L}, + {4.5L, -7.907760977504203261110e-01L}, + {5.0L, 1.455527162532881262341e+00L}, + {5.5L, -9.898201026724628026772e-01L}, + {6.0L, -1.086042272127976560819e+00L}, + {6.5L, 2.881956231985904298958e+00L}, + {7.0L, -9.678893576691628101581e-01L}, + {7.5L, -3.937973983229931876360e+00L}, + {8.0L, 3.987406017379465362777e+00L}, + {8.5L, 4.500447855060778233616e+00L}, + {9.0L, -7.925065331096693160232e+00L}, + {9.5L, -5.601628252101942569488e+00L}}; + + compareWithTable(Table, RelTolerance, AbsTolerance); +} + +template +void compareWithBoost(T RelTolerance, T AbsTolerance) { + // Generated with gcc 8.3.1 and boost 1.66.0 + std::map > Table; + + Table[10] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, -3.893744141378520446024e-01L}, + {1.0L, 4.189459325396825396478e-01L}, + {1.5L, 4.180177743094308038933e-01L}, + {2.0L, -3.090652557319223986041e-01L}, + {2.5L, -8.802526766660982968679e-01L}, + {3.0L, -7.000223214285714285572e-01L}, + {3.5L, 2.211135054223331402863e-01L}, + {4.0L, 1.379259259259259258966e+00L}, + {4.5L, 2.078780713762555803267e+00L}, + {5.0L, 1.756276179453262786681e+00L}, + {5.5L, 2.235215177990141367499e-01L}, + {6.0L, -2.222857142857142857383e+00L}, + {6.5L, -4.837352781211678102568e+00L}, + {7.0L, -6.589936342592592592900e+00L}, + {7.5L, -6.421207427978515625000e+00L}, + {8.0L, -3.538694885361552027511e+00L}, + {8.5L, 2.300054675712787284570e+00L}, + {9.0L, 1.057033482142857142797e+01L}, + {9.5L, 1.987481935988868565236e+01L}}; + Table[20] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, 3.120174870800154149416e-01L}, + {1.0L, -1.642588118277929597872e-01L}, + {1.5L, -3.415329100484649494361e-01L}, + {2.0L, 5.007215675600304564149e-01L}, + {2.5L, 4.695899854819089239144e-01L}, + {3.0L, -5.575093252171908238387e-01L}, + {3.5L, -1.070298248145694926638e+00L}, + {4.0L, -1.507914281021107639542e-01L}, + {4.5L, 1.437239150257817378478e+00L}, + {5.0L, 2.020225744476913575616e+00L}, + {5.5L, 5.842074183986535974638e-01L}, + {6.0L, -2.226540971456672473108e+00L}, + {6.5L, -4.371331381805792056210e+00L}, + {7.0L, -3.714934120411960050481e+00L}, + {7.5L, 4.510955524973794070892e-01L}, + {8.0L, 6.568596673285710329897e+00L}, + {8.5L, 1.118547863037195619958e+01L}, + {9.0L, 1.048071357226878702780e+01L}, + {9.5L, 2.334272934394408312297e+00L}}; + Table[50] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, -3.181388060269979063839e-01L}, + {1.0L, 2.329893629018530610989e-01L}, + {1.5L, -2.568364397519475610836e-01L}, + {2.0L, 4.416561042694050873874e-01L}, + {2.5L, -5.522468012177843670876e-01L}, + {3.0L, 1.463786536365346786293e-01L}, + {3.5L, 7.431381187021829678833e-01L}, + {4.0L, -8.261038021062498686205e-01L}, + {4.5L, -7.795068145562651415578e-01L}, + {5.0L, 1.473525881943055794421e+00L}, + {5.5L, 1.221520285433004828868e+00L}, + {6.0L, -2.083161364452417261195e+00L}, + {6.5L, -2.649222481840662724029e+00L}, + {7.0L, 1.961343956455813756449e+00L}, + {7.5L, 5.387138713354476475276e+00L}, + {8.0L, 8.689419213752169962399e-01L}, + {8.5L, -7.772031120599373790140e+00L}, + {9.0L, -8.576161635248512751274e+00L}, + {9.5L, 3.641234604888676076005e+00L}}; + Table[100] = {{0.0L, 1.000000000000000000000e+00L}, + {0.5L, 1.868226036769227888667e-01L}, + {1.0L, 2.706659209933101110067e-01L}, + {1.5L, 6.389168619635640514534e-02L}, + {2.0L, -3.031283692749120640291e-01L}, + {2.5L, 4.257949400323675464002e-01L}, + {3.0L, -4.794265341305422027554e-01L}, + {3.5L, 3.853183234834782501781e-01L}, + {4.0L, 2.339473046150550478781e-02L}, + {4.5L, -7.907760977504203251352e-01L}, + {5.0L, 1.455527162532881259196e+00L}, + {5.5L, -9.898201026724628030025e-01L}, + {6.0L, -1.086042272127976562228e+00L}, + {6.5L, 2.881956231985904306981e+00L}, + {7.0L, -9.678893576691628229517e-01L}, + {7.5L, -3.937973983229931877661e+00L}, + {8.0L, 3.987406017379465369499e+00L}, + {8.5L, 4.500447855060778231881e+00L}, + {9.0L, -7.925065331096693160232e+00L}, + {9.5L, -5.601628252101942564284e+00L}}; + + compareWithTable(Table, RelTolerance, AbsTolerance); +} + +int main(int, char**) { + compareWithBoost(1e-6f, 0.0f); + compareWithBoost(1e-12, 0.0); + + compareWithGCC(1e-6f, 0.0f); + compareWithGCC(1e-12, 0.0); +} diff --git a/libcxx/test/std/numerics/laguerre/laguerre_exceptions.pass.cpp b/libcxx/test/std/numerics/laguerre/laguerre_exceptions.pass.cpp new file mode 100644 --- /dev/null +++ b/libcxx/test/std/numerics/laguerre/laguerre_exceptions.pass.cpp @@ -0,0 +1,47 @@ +//===----------------------------------------------------------------------===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// + +// UNSUPPORTED: c++98, c++03, c++11, c++14, libcpp-no-exceptions + +// + +// floating_point laguerre (unsigned n, arithmetic x); +// float laguerref(unsigned n, float x); +// long double laguerrel(unsigned n, long double x); + +#include +#include +#include +#include + +/// \brief std::laguerre(n,x) throws if x < 0 and not NaN. +/// +/// \tparam T one of float, double, or long double +template +void testInvalidArguments() { + T InvalidArguments[] = {T(-5), -std::numeric_limits::min()}; + for (unsigned n = 0; n < 128; ++n) { + for (T x : InvalidArguments) { + bool ThrowsDomainError = false; + try { + const T Value = std::laguerre(n, x); + (void)Value; + } catch (const std::domain_error&) { + ThrowsDomainError = true; + } catch (...) { + } + assert(ThrowsDomainError); + } + } +} + +int main(int, char**) { + testInvalidArguments(); + testInvalidArguments(); + testInvalidArguments(); +}