blob: 9f1ef8b41635d586d39aa9d557c91a26464eef0d [file] [log] [blame]
[template tr1_overview[]
Many of the special functions included in this library are also a part
of the either the [C99] or the [tr1]. Therefore this library
includes a thin wrapper header `boost/math/tr1.hpp` that provides
compatibility with these two standards.
There are various pros and cons to using the library in this way:
Pros:
* The header to include is lightweight (i.e. fast to compile).
* The functions have extern "C" linkage,
and so are usable from other languages (not just C and C++).
* C99 and C++ TR1 Standard compatibility.
Cons:
* You will need to compile and link to the external Boost.Math libraries.
* Limited to support for the types, `float`, `double` and `long double`.
* Error handling is handled via setting ::errno and returning NaN's and
infinities: this may be less flexible than an C++ exception based approach.
[note The separate libraries are required *only* if you choose to use
boost/math/tr1.hpp rather than some other Boost.Math header, the rest
of Boost.Math remains header-only.]
The separate libraries required in order to use tr1.hpp can be compiled
using bjam from within the libs/math/build directory, or from the Boost
root directory using the usual Boost-wide install procedure.
Alternatively the source files are located in libs/math/src and each have
the same name as the function they implement. The various libraries are
named as follows:
[table
[[Name][Type][Functions]]
[[boost_math_c99f-<suffix>][float][C99 Functions]]
[[boost_math_c99-<suffix>][double][C99 Functions]]
[[boost_math_c99l-<suffix>][long double][C99 Functions]]
[[boost_math_tr1f-<suffix>][float][TR1 Functions]]
[[boost_math_tr1-<suffix>][double][TR1 Functions]]
[[boost_math_tr1l-<suffix>][long double][TR1 Functions]]
]
Where `<suffix>` encodes the compiler and build options used to
build the libraries: for example "libboost_math_tr1-vc80-mt-gd.lib"
would be the statically linked TR1 library to use with Visual C++ 8.0,
in multithreading debug mode, with the DLL VC++ runtime, where as
"boost_math_tr1-vc80-mt.lib" would be import library for the TR1 DLL
to be used with Visual C++ 8.0 with the release multithreaded DLL
VC++ runtime.
Refer to the getting started guide for a
[@http://www.boost.org/doc/libs/1_35_0/more/getting_started/windows.html#library-naming
full explanation of the `<suffix>` meanings].
[note
Visual C++ users will typically have the correct library variant to link
against selected for them by boost/math/tr1.hpp based on your compiler
settings.
Users will need to define BOOST_MATH_TR1_DYN_LINK when building
their code if they want to link against the DLL versions of these libraries
rather than the static versions.
Users can disable auto-linking by defining BOOST_MATH_TR1_NO_LIB when
building: this is typically only used when linking against a customised
build of the libraries.]
[note Linux and Unix users will generally only have one variant of
these libraries installed, and can generally just link against
-lboost_math_tr1 etc.]
[h4 Usage Recomendations]
This library now presents the user with a choice:
* To include the header only versions of the functions and have
an easier time linking, but a longer compile time.
* To include the TR1 headers and link against an external library.
Which option you choose depends largely on how you prefer to work
and how your system is set up.
For example a casual user who just needs the acosh function, would
probably be better off including `<boost/math/special_functions/acosh.hpp>`
and using `boost::math::acosh(x)` in their code.
However, for large scale
software development where compile times are significant, and where the
Boost libraries are already built and installed on the system, then
including `<boost/math/tr1.hpp>` and using `boost::math::tr1::acosh(x)`
will speed up compile times, reduce object files sizes (since there are
no templates being instantiated any more), and also speed up debugging
runtimes - since the externally compiled libraries can be
compiler optimised, rather than built using full settings - the difference
in performance between
[link math_toolkit.perf.getting_best
release and debug builds can be as much as 20 times],
so for complex applications this can be a big win.
[h4 Supported C99 Functions]
See also the [link math_toolkit.extern_c.c99
quick reference guide for these functions].
namespace boost{ namespace math{ namespace tr1{ extern "C"{
typedef unspecified float_t;
typedef unspecified double_t;
double acosh(double x);
float acoshf(float x);
long double acoshl(long double x);
double asinh(double x);
float asinhf(float x);
long double asinhl(long double x);
double atanh(double x);
float atanhf(float x);
long double atanhl(long double x);
double cbrt(double x);
float cbrtf(float x);
long double cbrtl(long double x);
double copysign(double x, double y);
float copysignf(float x, float y);
long double copysignl(long double x, long double y);
double erf(double x);
float erff(float x);
long double erfl(long double x);
double erfc(double x);
float erfcf(float x);
long double erfcl(long double x);
double expm1(double x);
float expm1f(float x);
long double expm1l(long double x);
double fmax(double x, double y);
float fmaxf(float x, float y);
long double fmaxl(long double x, long double y);
double fmin(double x, double y);
float fminf(float x, float y);
long double fminl(long double x, long double y);
double hypot(double x, double y);
float hypotf(float x, float y);
long double hypotl(long double x, long double y);
double lgamma(double x);
float lgammaf(float x);
long double lgammal(long double x);
long long llround(double x);
long long llroundf(float x);
long long llroundl(long double x);
double log1p(double x);
float log1pf(float x);
long double log1pl(long double x);
long lround(double x);
long lroundf(float x);
long lroundl(long double x);
double nextafter(double x, double y);
float nextafterf(float x, float y);
long double nextafterl(long double x, long double y);
double nexttoward(double x, long double y);
float nexttowardf(float x, long double y);
long double nexttowardl(long double x, long double y);
double round(double x);
float roundf(float x);
long double roundl(long double x);
double tgamma(double x);
float tgammaf(float x);
long double tgammal(long double x);
double trunc(double x);
float truncf(float x);
long double truncl(long double x);
}}}} // namespaces
[h4 Supported TR1 Functions]
See also the [link math_toolkit.extern_c.tr1
quick reference guide for these functions].
namespace boost{ namespace math{ namespace tr1{ extern "C"{
// [5.2.1.1] associated Laguerre polynomials:
double assoc_laguerre(unsigned n, unsigned m, double x);
float assoc_laguerref(unsigned n, unsigned m, float x);
long double assoc_laguerrel(unsigned n, unsigned m, long double x);
// [5.2.1.2] associated Legendre functions:
double assoc_legendre(unsigned l, unsigned m, double x);
float assoc_legendref(unsigned l, unsigned m, float x);
long double assoc_legendrel(unsigned l, unsigned m, long double x);
// [5.2.1.3] beta function:
double beta(double x, double y);
float betaf(float x, float y);
long double betal(long double x, long double y);
// [5.2.1.4] (complete) elliptic integral of the first kind:
double comp_ellint_1(double k);
float comp_ellint_1f(float k);
long double comp_ellint_1l(long double k);
// [5.2.1.5] (complete) elliptic integral of the second kind:
double comp_ellint_2(double k);
float comp_ellint_2f(float k);
long double comp_ellint_2l(long double k);
// [5.2.1.6] (complete) elliptic integral of the third kind:
double comp_ellint_3(double k, double nu);
float comp_ellint_3f(float k, float nu);
long double comp_ellint_3l(long double k, long double nu);
// [5.2.1.8] regular modified cylindrical Bessel functions:
double cyl_bessel_i(double nu, double x);
float cyl_bessel_if(float nu, float x);
long double cyl_bessel_il(long double nu, long double x);
// [5.2.1.9] cylindrical Bessel functions (of the first kind):
double cyl_bessel_j(double nu, double x);
float cyl_bessel_jf(float nu, float x);
long double cyl_bessel_jl(long double nu, long double x);
// [5.2.1.10] irregular modified cylindrical Bessel functions:
double cyl_bessel_k(double nu, double x);
float cyl_bessel_kf(float nu, float x);
long double cyl_bessel_kl(long double nu, long double x);
// [5.2.1.11] cylindrical Neumann functions;
// cylindrical Bessel functions (of the second kind):
double cyl_neumann(double nu, double x);
float cyl_neumannf(float nu, float x);
long double cyl_neumannl(long double nu, long double x);
// [5.2.1.12] (incomplete) elliptic integral of the first kind:
double ellint_1(double k, double phi);
float ellint_1f(float k, float phi);
long double ellint_1l(long double k, long double phi);
// [5.2.1.13] (incomplete) elliptic integral of the second kind:
double ellint_2(double k, double phi);
float ellint_2f(float k, float phi);
long double ellint_2l(long double k, long double phi);
// [5.2.1.14] (incomplete) elliptic integral of the third kind:
double ellint_3(double k, double nu, double phi);
float ellint_3f(float k, float nu, float phi);
long double ellint_3l(long double k, long double nu, long double phi);
// [5.2.1.15] exponential integral:
double expint(double x);
float expintf(float x);
long double expintl(long double x);
// [5.2.1.16] Hermite polynomials:
double hermite(unsigned n, double x);
float hermitef(unsigned n, float x);
long double hermitel(unsigned n, long double x);
// [5.2.1.18] Laguerre polynomials:
double laguerre(unsigned n, double x);
float laguerref(unsigned n, float x);
long double laguerrel(unsigned n, long double x);
// [5.2.1.19] Legendre polynomials:
double legendre(unsigned l, double x);
float legendref(unsigned l, float x);
long double legendrel(unsigned l, long double x);
// [5.2.1.20] Riemann zeta function:
double riemann_zeta(double);
float riemann_zetaf(float);
long double riemann_zetal(long double);
// [5.2.1.21] spherical Bessel functions (of the first kind):
double sph_bessel(unsigned n, double x);
float sph_besself(unsigned n, float x);
long double sph_bessell(unsigned n, long double x);
// [5.2.1.22] spherical associated Legendre functions:
double sph_legendre(unsigned l, unsigned m, double theta);
float sph_legendref(unsigned l, unsigned m, float theta);
long double sph_legendrel(unsigned l, unsigned m, long double theta);
// [5.2.1.23] spherical Neumann functions;
// spherical Bessel functions (of the second kind):
double sph_neumann(unsigned n, double x);
float sph_neumannf(unsigned n, float x);
long double sph_neumannl(unsigned n, long double x);
}}}} // namespaces
In addition sufficient additional overloads of the `double` versions of the
above functions are provided, so that calling the function with any mixture
of `float`, `double`, `long double`, or /integer/ arguments is supported, with the
return type determined by the __arg_pomotion_rules.
[h4 Currently Unsupported C99 Functions]
double exp2(double x);
float exp2f(float x);
long double exp2l(long double x);
double fdim(double x, double y);
float fdimf(float x, float y);
long double fdiml(long double x, long double y);
double fma(double x, double y, double z);
float fmaf(float x, float y, float z);
long double fmal(long double x, long double y, long double z);
int ilogb(double x);
int ilogbf(float x);
int ilogbl(long double x);
long long llrint(double x);
long long llrintf(float x);
long long llrintl(long double x);
double log2(double x);
float log2f(float x);
long double log2l(long double x);
double logb(double x);
float logbf(float x);
long double logbl(long double x);
long lrint(double x);
long lrintf(float x);
long lrintl(long double x);
double nan(const char *str);
float nanf(const char *str);
long double nanl(const char *str);
double nearbyint(double x);
float nearbyintf(float x);
long double nearbyintl(long double x);
double remainder(double x, double y);
float remainderf(float x, float y);
long double remainderl(long double x, long double y);
double remquo(double x, double y, int *pquo);
float remquof(float x, float y, int *pquo);
long double remquol(long double x, long double y, int *pquo);
double rint(double x);
float rintf(float x);
long double rintl(long double x);
double scalbln(double x, long ex);
float scalblnf(float x, long ex);
long double scalblnl(long double x, long ex);
double scalbn(double x, int ex);
float scalbnf(float x, int ex);
long double scalbnl(long double x, int ex);
[h4 Currently Unsupported TR1 Functions]
// [5.2.1.7] confluent hypergeometric functions:
double conf_hyperg(double a, double c, double x);
float conf_hypergf(float a, float c, float x);
long double conf_hypergl(long double a, long double c, long double x);
// [5.2.1.17] hypergeometric functions:
double hyperg(double a, double b, double c, double x);
float hypergf(float a, float b, float c, float x);
long double hypergl(long double a, long double b, long double c,
long double x);
]
[/
Copyright 2008 John Maddock and Paul A. Bristow.
Distributed under the Boost Software License, Version 1.0.
(See accompanying file LICENSE_1_0.txt or copy at
http://www.boost.org/LICENSE_1_0.txt).
]