Index: libcxx/include/cmath =================================================================== --- libcxx/include/cmath +++ 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,128 @@ return isfinite(__lcpp_x); } +_LIBCPP_FUNC_VIS +void __libcpp_cmath_domain_error() _NOEXCEPT; + +#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); + for (unsigned __i = 2; __i <= __n; ++__i) { + const _Real __li = ((_Real(2) * __i + __alpham1 - __x) * __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) _NOEXCEPT { + if (_VSTD::isnan(__x)) + return _VSTD::numeric_limits<_Real>::quiet_NaN(); + + if (__x < _Real(0)) { + __libcpp_cmath_domain_error(); + return _VSTD::numeric_limits<_Real>::quiet_NaN(); + } + + return __libcpp_generalized_laguerre_recurrence(__n, _Real(__m), __x); +} + +template +inline _LIBCPP_INLINE_VISIBILITY +_Real __libcpp_laguerre(unsigned __n, _Real __x) _NOEXCEPT { + if (_VSTD::isnan(__x)) + return _VSTD::numeric_limits<_Real>::quiet_NaN(); + + if (__x < _Real(0)) { + __libcpp_cmath_domain_error(); + return _VSTD::numeric_limits<_Real>::quiet_NaN(); + } + + 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) _NOEXCEPT { + 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) _NOEXCEPT { + 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) _NOEXCEPT { + 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) _NOEXCEPT { + 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) _NOEXCEPT { + return laguerre(__lcpp_n, __lcpp_x); +} + +inline _LIBCPP_INLINE_VISIBILITY +long double laguerrel(unsigned __lcpp_n, long double __lcpp_x) _NOEXCEPT { + return laguerre(__lcpp_n, __lcpp_x); +} + +#endif // _LIBCPP_STD_VER > 14 + + + #if _LIBCPP_STD_VER > 17 template constexpr Index: libcxx/src/CMakeLists.txt =================================================================== --- libcxx/src/CMakeLists.txt +++ libcxx/src/CMakeLists.txt @@ -7,6 +7,7 @@ bind.cpp charconv.cpp chrono.cpp + cmath.cpp condition_variable.cpp debug.cpp exception.cpp Index: libcxx/src/cmath.cpp =================================================================== --- /dev/null +++ libcxx/src/cmath.cpp @@ -0,0 +1,35 @@ +//===------------------------------ cmath.cpp -----------------------------===// +// +// 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 +// +//===----------------------------------------------------------------------===// + +#include "cmath" +#include "cerrno" +#include "cfenv" + +_LIBCPP_BEGIN_NAMESPACE_STD + +/// It is unspecified whether math_errhandling is a macro or an +/// identifier with external linkage (see J.1 of the C11 standard). +inline _LIBCPP_INLINE_VISIBILITY bool __libcpp_cmath_use_errno() _NOEXCEPT { + return (math_errhandling & MATH_ERRNO) != 0; +} + +inline _LIBCPP_INLINE_VISIBILITY bool +__libcpp_cmath_use_floating_point_exceptions() _NOEXCEPT { + return (math_errhandling & MATH_ERREXCEPT) != 0; +} + +/// \brief Signal errors from mathematics functions according to +/// section 7.12.1 of the C11 standard. +void __libcpp_cmath_domain_error() _NOEXCEPT { + if (__libcpp_cmath_use_floating_point_exceptions()) + _VSTD::feraiseexcept(FE_INVALID); + if (__libcpp_cmath_use_errno()) + errno = EDOM; +} + +_LIBCPP_END_NAMESPACE_STD Index: libcxx/test/std/numerics/c.math/c.math.laguerre/assoc_laguerre_basic.pass.cpp =================================================================== --- /dev/null +++ libcxx/test/std/numerics/c.math/c.math.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); +} Index: libcxx/test/std/numerics/c.math/c.math.laguerre/assoc_laguerre_comparisons.pass.cpp =================================================================== --- /dev/null +++ libcxx/test/std/numerics/c.math/c.math.laguerre/assoc_laguerre_comparisons.pass.cpp @@ -0,0 +1,263 @@ +//===----------------------------------------------------------------------===// +// +// 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 n = 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(n, 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 + // Confirmed with msvc-192027508 and boost-1.69.0 + // using boost::multiprecision::cpp_bin_float_oct + std::map, std::map > Table; + Table[{10u, 10u}] = { + {0.000000000e+00f, 1.847560000000000000000000000000000000e+05L}, + {5.000000000e-01f, 1.150042533944161228401951058201058201e+05L}, + {1.000000000e+00f, 6.835330763365299823633156966490299824e+04L}, + {1.500000000e+00f, 3.819497840595790318080357142857142857e+04L}, + {2.000000000e+00f, 1.952449269841269841269841269841269841e+04L}, + {2.500000000e+00f, 8.617790029165816265018738977072310406e+03L}, + {3.000000000e+00f, 2.761132790178571428571428571428571429e+03L}, + {3.500000000e+00f, 2.647443095454463252314814814814814815e+01L}, + {4.000000000e+00f, -9.133717107583774250440917107583774250e+02L}, + {4.500000000e+00f, -9.351607177516392299107142857142857143e+02L}, + {5.000000000e+00f, -5.861136119378306878306878306878306878e+02L}, + {5.500000000e+00f, -1.796631517410951106426366843033509700e+02L}, + {6.000000000e+00f, 1.300342857142857142857142857142857143e+02L}, + {6.500000000e+00f, 2.913334329296293712797619047619047619e+02L}, + {7.000000000e+00f, 3.129948823302469135802469135802469136e+02L}, + {7.500000000e+00f, 2.340251775469098772321428571428571429e+02L}, + {8.000000000e+00f, 1.030916402116402116402116402116402116e+02L}, + {8.500000000e+00f, -3.454145314884353987544091710758377425e+01L}, + {9.000000000e+00f, -1.443173883928571428571428571428571429e+02L}, + {9.500000000e+00f, -2.054278306854086578207671957671957672e+02L}}; + Table[{20u, 10u}] = { + {0.000000000e+00f, 3.004501500000000000000000000000000000e+07L}, + {5.000000000e-01f, 1.138930927314853253029395339203192700e+07L}, + {1.000000000e+00f, 3.651201251651593485217380578844588588e+06L}, + {1.500000000e+00f, 8.473698098486793286200561209213255411e+05L}, + {2.000000000e+00f, 4.566251493573249583609543057909799018e+04L}, + {2.500000000e+00f, -7.338348878522165050610581250007527262e+04L}, + {3.000000000e+00f, -3.208538228843214609843437052774735216e+04L}, + {3.500000000e+00f, 5.614549693109722522609450388942484022e+03L}, + {4.000000000e+00f, 1.414790961935669154114686973704524845e+04L}, + {4.500000000e+00f, 6.838363861624681720329730407685479615e+03L}, + {5.000000000e+00f, -1.795486420939170528844540670782701058e+03L}, + {5.500000000e+00f, -5.316437926431718982436197192994482959e+03L}, + {6.000000000e+00f, -3.887263891795776094758845753980650044e+03L}, + {6.500000000e+00f, -4.294692178238768076149262089473078895e+02L}, + {7.000000000e+00f, 2.300650817841901842105515027636280296e+03L}, + {7.500000000e+00f, 3.004218471787209729032542514375576839e+03L}, + {8.000000000e+00f, 1.857655658194353204937035610305568061e+03L}, + {8.500000000e+00f, -1.126184178317543084796255368866824306e+02L}, + {9.000000000e+00f, -1.763110242143649015988952169016299932e+03L}, + {9.500000000e+00f, -2.355059370043133450966538634289651068e+03L}}; + Table[{20u, 20u}] = { + {0.000000000e+00f, 1.378465288200000000000000000000000000e+11L}, + {5.000000000e-01f, 8.468295361673511391094370546313916410e+10L}, + {1.000000000e+00f, 5.077912460175697549246274616499530651e+10L}, + {1.500000000e+00f, 2.960576939629002842690542206765752763e+10L}, + {2.000000000e+00f, 1.669670827385554016835834451769714764e+10L}, + {2.500000000e+00f, 9.043858989194316966879219754251941826e+09L}, + {3.000000000e+00f, 4.656130531018169693851517930219106413e+09L}, + {3.500000000e+00f, 2.241401912717946666170865381396122973e+09L}, + {4.000000000e+00f, 9.800663873313004746557925074295885128e+08L}, + {4.500000000e+00f, 3.659874821411927665563127613828386379e+08L}, + {5.000000000e+00f, 9.651825546133988803537110943821953031e+07L}, + {5.500000000e+00f, -2.223016973336849317163394181832929881e+06L}, + {6.000000000e+00f, -2.510182459882302288600386787029865314e+07L}, + {6.500000000e+00f, -2.030857764686113736095824731440825658e+07L}, + {7.000000000e+00f, -9.790386453177882433133401337967622734e+06L}, + {7.500000000e+00f, -1.710558937907162579124982278576880798e+06L}, + {8.000000000e+00f, 2.291830176433364632498633077422097904e+06L}, + {8.500000000e+00f, 3.127649252979586820518042226411356033e+06L}, + {9.000000000e+00f, 2.276440826891129583391996555247052815e+06L}, + {9.500000000e+00f, 9.499740595601278431761433792649449840e+05L}}; + Table[{50u, 10u}] = { + {0.000000000e+00f, 7.539402756600000000000000000000000000e+10L}, + {5.000000000e-01f, 5.365943119336201378134035359882466723e+09L}, + {1.000000000e+00f, -6.058846804287141383677660607814191704e+07L}, + {1.500000000e+00f, -5.428935628098645254618824494613585085e+06L}, + {2.000000000e+00f, 6.491875483166804221363183930819364769e+06L}, + {2.500000000e+00f, -3.125723663735501903981495567951833431e+06L}, + {3.000000000e+00f, 3.215290936254617409894258836545016817e+05L}, + {3.500000000e+00f, 7.737700630158116543339424386027985510e+05L}, + {4.000000000e+00f, -3.911037482683599918682066744552496798e+05L}, + {4.500000000e+00f, -2.527733547889387020232815706465346664e+05L}, + {5.000000000e+00f, 2.294630650721308950777906418429615007e+05L}, + {5.500000000e+00f, 1.446393186292154792785216573589521607e+05L}, + {6.000000000e+00f, -1.282699408320068905328231784843734305e+05L}, + {6.500000000e+00f, -1.258756273224852413426599162361074252e+05L}, + {7.000000000e+00f, 5.353708357092831456614142683643980840e+04L}, + {7.500000000e+00f, 1.185268439074285104446392014084015298e+05L}, + {8.000000000e+00f, 1.536260133508163855344880117997141767e+04L}, + {8.500000000e+00f, -9.111765379659439550361616786814448742e+04L}, + {9.000000000e+00f, -7.288173794915159594043582343607864974e+04L}, + {9.500000000e+00f, 2.992603738587986627890005773037749811e+04L}}; + Table[{50u, 20u}] = { + {0.000000000e+00f, 1.618846036626578760000000000000000000e+17L}, + {5.000000000e-01f, 4.682199741184386787306309877587393591e+16L}, + {1.000000000e+00f, 1.197682287075144468099789990123105742e+16L}, + {1.500000000e+00f, 2.562563860642318863254850569301798172e+15L}, + {2.000000000e+00f, 4.025512140617272670581356478051263732e+14L}, + {2.500000000e+00f, 2.558241994957465181484554568176374268e+13L}, + {3.000000000e+00f, -7.107135583922841746185410103050560329e+12L}, + {3.500000000e+00f, -1.527944352509347694013650035435234584e+12L}, + {4.000000000e+00f, 4.962109792862451580365778910400865798e+11L}, + {4.500000000e+00f, 2.207360501475066801653910395634842569e+11L}, + {5.000000000e+00f, -5.495990181724055273037968440911833256e+10L}, + {5.500000000e+00f, -5.395164965408203375217225726208229013e+10L}, + {6.000000000e+00f, 2.160824469618739428466565631898423015e+09L}, + {6.500000000e+00f, 1.625278916791839990306321835579021240e+10L}, + {7.000000000e+00f, 4.561836324614418239176370604718990970e+09L}, + {7.500000000e+00f, -4.311154213993686890627225570074389169e+09L}, + {8.000000000e+00f, -3.772710079809093438010731705317810654e+09L}, + {8.500000000e+00f, 6.289714539106999181133286266989056814e+07L}, + {9.000000000e+00f, 1.862043401563646102238160053582182304e+09L}, + {9.500000000e+00f, 1.093508422577525383092594327891707390e+09L}}; + Table[{50u, 50u}] = { + {0.000000000e+00f, 1.008913445455641933348124972560000000e+29L}, + {5.000000000e-01f, 6.150335228601113622655274327476517634e+28L}, + {1.000000000e+00f, 3.712351462622653280180388413124362151e+28L}, + {1.500000000e+00f, 2.217375806127585064244642685788011723e+28L}, + {2.000000000e+00f, 1.309702308652025846947762457749189316e+28L}, + {2.500000000e+00f, 7.643899665374721928130027177404479787e+27L}, + {3.000000000e+00f, 4.404413646826269803469256712809515693e+27L}, + {3.500000000e+00f, 2.502994762342186408737089901513642102e+27L}, + {4.000000000e+00f, 1.401305429217811544713965889548878920e+27L}, + {4.500000000e+00f, 7.718403057066373185582148275575174840e+26L}, + {5.000000000e+00f, 4.176012189047808124673913500291806295e+26L}, + {5.500000000e+00f, 2.215242482188571358317665570141917798e+26L}, + {6.000000000e+00f, 1.149522800633857380266071746750882174e+26L}, + {6.500000000e+00f, 5.818720365519067242354593399696793284e+25L}, + {7.000000000e+00f, 2.862862686765511838081738288136282216e+25L}, + {7.500000000e+00f, 1.362737172361160252482872914114289566e+25L}, + {8.000000000e+00f, 6.236142696693560631999275859574346499e+24L}, + {8.500000000e+00f, 2.718917011811881112775517614181728596e+24L}, + {9.000000000e+00f, 1.113959548616688992412790077416342014e+24L}, + {9.500000000e+00f, 4.190270485617718536487729672007575863e+23L}}; + Table[{100u, 10u}] = { + {0.000000000e+00f, 4.689763662398100000000000000000000000e+13L}, + {5.000000000e-01f, -4.294878427688371024531851061319102037e+09L}, + {1.000000000e+00f, 3.944092408530536255021300956559379230e+09L}, + {1.500000000e+00f, -2.013797229225323655763608919100640149e+08L}, + {2.000000000e+00f, -2.157299751648581342915805462513841785e+07L}, + {2.500000000e+00f, 1.459430547397130525591946163074198253e+07L}, + {3.000000000e+00f, -1.016058418748330157007669790028306271e+06L}, + {3.500000000e+00f, -6.908288622716346507642979053656636551e+06L}, + {4.000000000e+00f, 9.787558305163868596268774054090610805e+06L}, + {4.500000000e+00f, -8.176816207612928381922046863047783064e+06L}, + {5.000000000e+00f, 3.377298430833765118518934385161482380e+06L}, + {5.500000000e+00f, 1.735890291779104318734228859236942215e+06L}, + {6.000000000e+00f, -3.807192070466917535765572118312904221e+06L}, + {6.500000000e+00f, 1.734983554763113700932087132127203438e+06L}, + {7.000000000e+00f, 1.706476023943294065526440991182059084e+06L}, + {7.500000000e+00f, -2.355218427716766162855871644138864339e+06L}, + {8.000000000e+00f, -3.058164836881507326950573697268896555e+05L}, + {8.500000000e+00f, 2.149540815590451995873247829981956679e+06L}, + {9.000000000e+00f, -3.345374513274238428141122383104290356e+05L}, + {9.500000000e+00f, -1.920373498643499778407713232713532093e+06L}}; + Table[{100u, 20u}] = { + {0.000000000e+00f, 2.946222729117663571812600000000000000e+22L}, + {5.000000000e-01f, 2.261308869389457759796235383014297726e+21L}, + {1.000000000e+00f, 9.475799170461032113561945315589383841e+19L}, + {1.500000000e+00f, -3.637696286243510429249175972473207315e+17L}, + {2.000000000e+00f, -2.432420166605434024442549600364329977e+16L}, + {2.500000000e+00f, 9.128749654970749467372896617177986571e+15L}, + {3.000000000e+00f, -2.745178778987279971099338286175840260e+15L}, + {3.500000000e+00f, 7.321257828016063552325086442582636844e+14L}, + {4.000000000e+00f, -1.137164349493819840359152415741650520e+14L}, + {4.500000000e+00f, -3.827471412315585203853000359476516937e+13L}, + {5.000000000e+00f, 4.041283789325254422512662801481139629e+13L}, + {5.500000000e+00f, -1.225276122801313933983424613922360422e+13L}, + {6.000000000e+00f, -4.723168970277397798041090588121513695e+12L}, + {6.500000000e+00f, 5.788209026781585843863266023986448368e+12L}, + {7.000000000e+00f, -6.214739115929799348058182370662955557e+11L}, + {7.500000000e+00f, -2.016398224166447073157533140373654060e+12L}, + {8.000000000e+00f, 8.280834859221748541813477062071812611e+11L}, + {8.500000000e+00f, 7.409824494487057609609125215091391500e+11L}, + {9.000000000e+00f, -5.314481262655515997135741587208250356e+11L}, + {9.500000000e+00f, -3.340793205142771007300366607840246511e+11L}}; + Table[{100u, 50u}] = { + {0.000000000e+00f, 2.012866090973193229424023438092931575e+40L}, + {5.000000000e-01f, 7.443534596401896623181000999994559399e+39L}, + {1.000000000e+00f, 2.669138788507098280725287624431078965e+39L}, + {1.500000000e+00f, 9.249176234411826068277724127925644582e+38L}, + {2.000000000e+00f, 3.084251014932004953935189222010128877e+38L}, + {2.500000000e+00f, 9.845130400028869444876995434105101614e+37L}, + {3.000000000e+00f, 2.987893170071304389357855943231028346e+37L}, + {3.500000000e+00f, 8.543446586522976373254811628719854682e+36L}, + {4.000000000e+00f, 2.272446170635370969567731464878334522e+36L}, + {4.500000000e+00f, 5.516531148637918792441174751622387006e+35L}, + {5.000000000e+00f, 1.184363719341424928994025359792846372e+35L}, + {5.500000000e+00f, 2.115939411877618805051987208371818521e+34L}, + {6.000000000e+00f, 2.679064308539579231700268292049484381e+33L}, + {6.500000000e+00f, 6.930372893313020732819287711452508203e+31L}, + {7.000000000e+00f, -7.006151656594741380171312277322640990e+31L}, + {7.500000000e+00f, -1.654854053550525555488608020539259653e+31L}, + {8.000000000e+00f, 5.124436639236123928610297788419312856e+29L}, + {8.500000000e+00f, 1.212547004925553863405248897314245135e+30L}, + {9.000000000e+00f, 2.130085958501174467534602948195510460e+29L}, + {9.500000000e+00f, -6.983885241918135688715753929980327426e+28L}}; + Table[{100u, 100u}] = { + {0.000000000e+00f, 9.054851465610328116540417707748416387e+58L}, + {5.000000000e-01f, 5.505859466697853236589965073011607970e+58L}, + {1.000000000e+00f, 3.331262988031657888739721478546839075e+58L}, + {1.500000000e+00f, 2.005242577321074531921026075680558361e+58L}, + {2.000000000e+00f, 1.200686309789629297925004812989445025e+58L}, + {2.500000000e+00f, 7.150282668151546649652638611573067949e+57L}, + {3.000000000e+00f, 4.234184715733205326717133474264394384e+57L}, + {3.500000000e+00f, 2.492794079123697662955583674749617412e+57L}, + {4.000000000e+00f, 1.458763178926480348998289096865283510e+57L}, + {4.500000000e+00f, 8.483418739565801741761908158999327307e+56L}, + {5.000000000e+00f, 4.901672141399097797235587554127969820e+56L}, + {5.500000000e+00f, 2.813178471492599512166658564973355635e+56L}, + {6.000000000e+00f, 1.603299095004078549094596123165367278e+56L}, + {6.500000000e+00f, 9.071348729442847395435232033972470741e+55L}, + {7.000000000e+00f, 5.093722441315075919170677202982054817e+55L}, + {7.500000000e+00f, 2.837656817120185205018156958262388528e+55L}, + {8.000000000e+00f, 1.567791758256770324224113366837280217e+55L}, + {8.500000000e+00f, 8.587159575060838641628619363348101396e+54L}, + {9.000000000e+00f, 4.660753669629245412882516558571331123e+54L}, + {9.500000000e+00f, 2.505548374372461099929474613080159001e+54L}}; + compareWithTable(Table, RelTolerance, AbsTolerance); +} + +int main(int, char**) { + compareWithBoost(1e-6f, 0.0f); + compareWithBoost(1e-12, 0.0); +} Index: libcxx/test/std/numerics/c.math/c.math.laguerre/assoc_laguerre_exceptions.pass.cpp =================================================================== --- /dev/null +++ libcxx/test/std/numerics/c.math/c.math.laguerre/assoc_laguerre_exceptions.pass.cpp @@ -0,0 +1,53 @@ +//===----------------------------------------------------------------------===// +// +// 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 + +/// \brief std::assoc_laguerre(n, m, x) throws if x < 0 and not NaN. +/// +/// \tparam T one of float, double, or long double +template +void testErrnoOnInvalidArgument() { + const T InvalidArguments[] = {T(-5), -std::numeric_limits::min()}; + for (unsigned n = 0; n < 128; ++n) { + for (unsigned m = 0; m < 128; ++m) { + for (T x : InvalidArguments) { + errno = 0; + const T Value = std::assoc_laguerre(n, m, x); + assert(errno == EDOM); + // The values returned by the mathematics functions on domain errors + // are implementation-defined according to J.3.12 of the C11 standard. + // We are going for NaN here. + assert(std::isnan(Value)); + } + } + } +} + +bool doMathematicsFunctionsSetErrno() { + return (math_errhandling & MATH_ERRNO) != 0; +} + +int main(int, char**) { + if (doMathematicsFunctionsSetErrno()) { + testErrnoOnInvalidArgument(); + testErrnoOnInvalidArgument(); + testErrnoOnInvalidArgument(); + } +} Index: libcxx/test/std/numerics/c.math/c.math.laguerre/laguerre_basic.pass.cpp =================================================================== --- /dev/null +++ libcxx/test/std/numerics/c.math/c.math.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); +} Index: libcxx/test/std/numerics/c.math/c.math.laguerre/laguerre_comparisons.pass.cpp =================================================================== --- /dev/null +++ libcxx/test/std/numerics/c.math/c.math.laguerre/laguerre_comparisons.pass.cpp @@ -0,0 +1,136 @@ +//===----------------------------------------------------------------------===// +// +// 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 ExpectedValue = static_cast(xy.second); + if (std::isinf(ExpectedValue)) + continue; + const T Value = std::laguerre(n, static_cast(xy.first)); + const T Error = std::abs(ExpectedValue - Value); + const T Tolerance = std::abs(ExpectedValue) * RelTolerance + AbsTolerance; + assert(Error <= Tolerance); + } + } +} + +template +void compareWithBoost(T RelTolerance, T AbsTolerance) { + // Generated with gcc-8.3.1 and boost-1.66.0 + // Confirmed with msvc-192027508 and boost-1.69.0 + // using boost::multiprecision::cpp_bin_float_oct + std::map > Table; + + Table[10u] = { + {0.000000000e+00f, 1.000000000000000000000000000000000000e+00L}, + {5.000000000e-01f, -3.893744141378520447530864197530864198e-01L}, + {1.000000000e+00f, 4.189459325396825396825396825396825397e-01L}, + {1.500000000e+00f, 4.180177743094308035714285714285714286e-01L}, + {2.000000000e+00f, -3.090652557319223985890652557319223986e-01L}, + {2.500000000e+00f, -8.802526766660982969576719576719576720e-01L}, + {3.000000000e+00f, -7.000223214285714285714285714285714286e-01L}, + {3.500000000e+00f, 2.211135054223331404320987654320987654e-01L}, + {4.000000000e+00f, 1.379259259259259259259259259259259259e+00L}, + {4.500000000e+00f, 2.078780713762555803571428571428571429e+00L}, + {5.000000000e+00f, 1.756276179453262786596119929453262787e+00L}, + {5.500000000e+00f, 2.235215177990141369047619047619047619e-01L}, + {6.000000000e+00f, -2.222857142857142857142857142857142857e+00L}, + {6.500000000e+00f, -4.837352781211678102954144620811287478e+00L}, + {7.000000000e+00f, -6.589936342592592592592592592592592593e+00L}, + {7.500000000e+00f, -6.421207427978515625000000000000000000e+00L}, + {8.000000000e+00f, -3.538694885361552028218694885361552028e+00L}, + {8.500000000e+00f, 2.300054675712787285052910052910052910e+00L}, + {9.000000000e+00f, 1.057033482142857142857142857142857143e+01L}, + {9.500000000e+00f, 1.987481935988868565365961199294532628e+01L}}; + Table[20u] = { + {0.000000000e+00f, 1.000000000000000000000000000000000000e+00L}, + {5.000000000e-01f, 3.120174870800154148915399248893113635e-01L}, + {1.000000000e+00f, -1.642588118277929599533376074856864704e-01L}, + {1.500000000e+00f, -3.415329100484649492676608230779311708e-01L}, + {2.000000000e+00f, 5.007215675600304565600885558368226549e-01L}, + {2.500000000e+00f, 4.695899854819089238467687303409922944e-01L}, + {3.000000000e+00f, -5.575093252171908239855566812665441591e-01L}, + {3.500000000e+00f, -1.070298248145694926198919867695008861e+00L}, + {4.000000000e+00f, -1.507914281021107637394178001324043627e-01L}, + {4.500000000e+00f, 1.437239150257817378525582974722170738e+00L}, + {5.000000000e+00f, 2.020225744476913576329950559141910796e+00L}, + {5.500000000e+00f, 5.842074183986535967964184734247564266e-01L}, + {6.000000000e+00f, -2.226540971456672473921478786582722894e+00L}, + {6.500000000e+00f, -4.371331381805792056672744105257936708e+00L}, + {7.000000000e+00f, -3.714934120411960051278165321215939774e+00L}, + {7.500000000e+00f, 4.510955524973794070410344366311559571e-01L}, + {8.000000000e+00f, 6.568596673285710330027732102272119709e+00L}, + {8.500000000e+00f, 1.118547863037195619875283558017243641e+01L}, + {9.000000000e+00f, 1.048071357226878702777347597117610386e+01L}, + {9.500000000e+00f, 2.334272934394408312494397341015663529e+00L}}; + Table[50u] = { + {0.000000000e+00f, 1.000000000000000000000000000000000000e+00L}, + {5.000000000e-01f, -3.181388060269979064951118308575628227e-01L}, + {1.000000000e+00f, 2.329893629018530615494699091219080976e-01L}, + {1.500000000e+00f, -2.568364397519475611466100723987632743e-01L}, + {2.000000000e+00f, 4.416561042694050876596607737004876946e-01L}, + {2.500000000e+00f, -5.522468012177843674236457431971580103e-01L}, + {3.000000000e+00f, 1.463786536365346789185198903782973029e-01L}, + {3.500000000e+00f, 7.431381187021829667512955593487581021e-01L}, + {4.000000000e+00f, -8.261038021062498676073778350077467278e-01L}, + {4.500000000e+00f, -7.795068145562651416494321484050019245e-01L}, + {5.000000000e+00f, 1.473525881943055793840385397812127879e+00L}, + {5.500000000e+00f, 1.221520285433004828956607124502786323e+00L}, + {6.000000000e+00f, -2.083161364452417261418198051651629275e+00L}, + {6.500000000e+00f, -2.649222481840662726160614188991140748e+00L}, + {7.000000000e+00f, 1.961343956455813751282062845447903523e+00L}, + {7.500000000e+00f, 5.387138713354476470625032695471140719e+00L}, + {8.000000000e+00f, 8.689419213752169947355208916454008896e-01L}, + {8.500000000e+00f, -7.772031120599373791734250499154946734e+00L}, + {9.000000000e+00f, -8.576161635248512747973711940416216024e+00L}, + {9.500000000e+00f, 3.641234604888676080221365466131811541e+00L}}; + Table[100u] = { + {0.000000000e+00f, 1.000000000000000000000000000000000000e+00L}, + {5.000000000e-01f, 1.868226036769227880076698172652800357e-01L}, + {1.000000000e+00f, 2.706659209933101112753498185332200218e-01L}, + {1.500000000e+00f, 6.389168619635640581860885965780472302e-02L}, + {2.000000000e+00f, -3.031283692749120637896753135855788028e-01L}, + {2.500000000e+00f, 4.257949400323675465623240248736101500e-01L}, + {3.000000000e+00f, -4.794265341305422028896879920086323876e-01L}, + {3.500000000e+00f, 3.853183234834782507208095016909928949e-01L}, + {4.000000000e+00f, 2.339473046150550347899313144984138234e-02L}, + {4.500000000e+00f, -7.907760977504203253534083320149925756e-01L}, + {5.000000000e+00f, 1.455527162532881260683435926668633937e+00L}, + {5.500000000e+00f, -9.898201026724628019853876600903030499e-01L}, + {6.000000000e+00f, -1.086042272127976558789672262834712399e+00L}, + {6.500000000e+00f, 2.881956231985904305362319452251036810e+00L}, + {7.000000000e+00f, -9.678893576691628160513659890082145166e-01L}, + {7.500000000e+00f, -3.937973983229931874952783679158484738e+00L}, + {8.000000000e+00f, 3.987406017379465374348023436284327674e+00L}, + {8.500000000e+00f, 4.500447855060778235662666322390440610e+00L}, + {9.000000000e+00f, -7.925065331096693164026460157953977998e+00L}, + {9.500000000e+00f, -5.601628252101942555995709594792649050e+00L}}; + compareWithTable(Table, RelTolerance, AbsTolerance); +} + +int main(int, char**) { + compareWithBoost(1e-6f, 0.0f); + compareWithBoost(1e-12, 0.0); +} Index: libcxx/test/std/numerics/c.math/c.math.laguerre/laguerre_exceptions.pass.cpp =================================================================== --- /dev/null +++ libcxx/test/std/numerics/c.math/c.math.laguerre/laguerre_exceptions.pass.cpp @@ -0,0 +1,52 @@ +//===----------------------------------------------------------------------===// +// +// 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 + +/// \brief Tests if std::laguerre(n, x) raises a domain error via errno +/// if x < 0 and not NaN. +/// +/// \tparam T one of float, double, or long double +template +void testErrnoOnInvalidArgument() { + T InvalidArguments[] = {T(-5), -std::numeric_limits::min()}; + for (unsigned n = 0; n < 128; ++n) { + for (T x : InvalidArguments) { + errno = 0; + const T Value = std::laguerre(n, x); + assert(errno == EDOM); + // The values returned by the mathematics functions on domain errors + // are implementation-defined according to J.3.12 of the C11 standard. + // We are going for NaN here. + assert(std::isnan(Value)); + } + } +} + +bool doMathematicsFunctionsSetErrno() { + return (math_errhandling & MATH_ERRNO) != 0; +} + +int main(int, char**) { + if (doMathematicsFunctionsSetErrno()) { + testErrnoOnInvalidArgument(); + testErrnoOnInvalidArgument(); + testErrnoOnInvalidArgument(); + } +} \ No newline at end of file Index: libcxx/test/std/numerics/c.math/cmath.pass.cpp =================================================================== --- libcxx/test/std/numerics/c.math/cmath.pass.cpp +++ 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 + template ()))> std::true_type has_abs_imp(int); template @@ -1553,6 +1558,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(); @@ -1626,5 +1671,11 @@ test_tgamma(); test_trunc(); + +#if _LIBCPP_STD_VER > 14 + test_assoc_laguerre(); + test_laguerre(); +#endif + return 0; }