[comment {-*- tcl -*- doctools manpage}]
[manpage_begin math::special n 0.4]
[keywords {Bessel functions}]
[keywords {error function}]
[keywords math]
[keywords {special functions}]
[copyright {2004 Arjen Markus <arjenmarkus@users.sourceforge.net>}]
[moddesc   {Tcl Math Library}]
[titledesc {Special mathematical functions}]
[category  Mathematics]
[require Tcl [opt 8.5]]
[require math::special [opt 0.5]]

[description]
[para]
This package implements several so-called special functions, like
the Gamma function, the Bessel functions and such.

[para]
Each function is implemented by a procedure that bears its name (well,
in close approximation):

[list_begin itemized]
[item]
J0 for the zeroth-order Bessel function of the first kind

[item]
J1 for the first-order Bessel function of the first kind

[item]
Jn for the nth-order Bessel function of the first kind

[item]
J1/2 for the half-order Bessel function of the first kind

[item]
J-1/2 for the minus-half-order Bessel function of the first kind

[item]
I_n for the modified Bessel function of the first kind of order n

[item]
Gamma for the Gamma function, erf and erfc for the error function and
the complementary error function

[item]
fresnel_C and fresnel_S for the Fresnel integrals

[item]
elliptic_K and elliptic_E (complete elliptic integrals)

[item]
exponent_Ei and other functions related to the so-called exponential
integrals

[item]
legendre, hermite: some of the classical orthogonal polynomials.

[list_end]

[section OVERVIEW]

In the following table several characteristics of the functions in this
package are summarized: the domain for the argument, the values for the
parameters and error bounds.

[example {
Family       | Function    | Domain x    | Parameter   | Error bound
-------------+-------------+-------------+-------------+--------------
Bessel       | J0, J1,     | all of R    | n = integer |   < 1.0e-8
             | Jn          |             |             |  (|x|<20, n<20)
Bessel       | J1/2, J-1/2,|  x > 0      | n = integer |   exact
Bessel       | I_n         | all of R    | n = integer |   < 1.0e-6
             |             |             |             |
Elliptic     | cn          | 0 <= x <= 1 |     --      |   < 1.0e-10
functions    | dn          | 0 <= x <= 1 |     --      |   < 1.0e-10
             | sn          | 0 <= x <= 1 |     --      |   < 1.0e-10
Elliptic     | K           | 0 <= x < 1  |     --      |   < 1.0e-6
integrals    | E           | 0 <= x < 1  |     --      |   < 1.0e-6
             |             |             |             |
Error        | erf         |             |     --      |
functions    | erfc        |             |             |
             |             |             |             |
Inverse      | invnorm     | 0 < x < 1   |     --      |   < 1.2e-9
normal       |             |             |             |
distribution |             |             |             |
             |             |             |             |
Exponential  | Ei          |  x != 0     |     --      |   < 1.0e-10 (relative)
integrals    | En          |  x >  0     |     --      |   as Ei
             | li          |  x > 0      |     --      |   as Ei
             | Chi         |  x > 0      |     --      |   < 1.0e-8
             | Shi         |  x > 0      |     --      |   < 1.0e-8
             | Ci          |  x > 0      |     --      |   < 2.0e-4
             | Si          |  x > 0      |     --      |   < 2.0e-4
             |             |             |             |
Fresnel      | C           |  all of R   |     --      |   < 2.0e-3
integrals    | S           |  all of R   |     --      |   < 2.0e-3
             |             |             |             |
general      | Beta        | (see Gamma) |     --      |   < 1.0e-9
             | Gamma       |  x != 0,-1, |     --      |   < 1.0e-9
             |             |  -2, ...    |             |
             | incBeta     |             |  a, b > 0   |   < 1.0e-9
             | regIncBeta  |             |  a, b > 0   |   < 1.0e-9
             | digamma     |  x != 0,-1  |             |   < 1.0e-9
             |             |  -2, ...    |             |
             |             |             |             |
             | sinc        |  all of R   |     --      |   exact
             |             |             |             |
orthogonal   | Legendre    |  all of R   | n = 0,1,... |   exact
polynomials  | Chebyshev   |  all of R   | n = 0,1,... |   exact
             | Laguerre    |  all of R   | n = 0,1,... |   exact
             |             |             | alpha el. R |
             | Hermite     |  all of R   | n = 0,1,... |   exact
}]

[emph Note:] Some of the error bounds are estimated, as no
"formal" bounds were available with the implemented approximation
method, others hold for the auxiliary functions used for estimating
the primary functions.

[para]
The following well-known functions are currently missing from the package:
[list_begin itemized]
[item]
Bessel functions of the second kind (Y_n, K_n)
[item]
Bessel functions of arbitrary order (and hence the Airy functions)
[item]
Chebyshev polynomials of the second kind (U_n)
[item]
The incomplete gamma function
[list_end]

[section "PROCEDURES"]

The package defines the following public procedures:

[list_begin definitions]

[call [cmd ::math::special::eulerNumber] [arg index]]

Return the index'th Euler number (note: these are integer values). As the
size of these numbers grows very fast, only a limited number are available.

[list_begin arguments]
[arg_def int index] Index of the number to be returned (should be between 0 and 54)
[list_end]

[para]

[call [cmd ::math::special::bernoulliNumber] [arg index]]

Return the index'th Bernoulli number. As the size of the numbers grows very fast,
only a limited number are available.

[list_begin arguments]
[arg_def int index] Index of the number to be returned (should be between 0 and 52)
[list_end]

[para]

[call [cmd ::math::special::Beta] [arg x] [arg y]]

Compute the Beta function for arguments "x" and "y"

[list_begin arguments]
[arg_def float x] First argument for the Beta function

[arg_def float y] Second argument for the Beta function
[list_end]

[para]

[call [cmd ::math::special::incBeta] [arg a] [arg b] [arg x]]

Compute the incomplete Beta function for argument "x" with parameters "a" and "b"

[list_begin arguments]
[arg_def float a] First parameter for the incomplete Beta function, a > 0
[arg_def float b] Second parameter for the incomplete Beta function, b > 0
[arg_def float x] Argument for the incomplete Beta function
[list_end]

[para]

[call [cmd ::math::special::regIncBeta] [arg a] [arg b] [arg x]]

Compute the regularized incomplete Beta function for argument "x" with parameters "a" and "b"

[list_begin arguments]
[arg_def float a] First parameter for the incomplete Beta function, a > 0
[arg_def float b] Second parameter for the incomplete Beta function, b > 0
[arg_def float x] Argument for the regularized incomplete Beta function
[list_end]

[para]

[call [cmd ::math::special::Gamma] [arg x]]

Compute the Gamma function for argument "x"

[list_begin arguments]
[arg_def float x] Argument for the Gamma function
[list_end]

[para]

[call [cmd ::math::special::digamma] [arg x]]

Compute the digamma function (psi) for argument "x"

[list_begin arguments]
[arg_def float x] Argument for the digamma function
[list_end]

[para]

[call [cmd ::math::special::erf] [arg x]]

Compute the error function for argument "x"

[list_begin arguments]
[arg_def float x] Argument for the error function
[list_end]

[para]

[call [cmd ::math::special::erfc] [arg x]]

Compute the complementary error function for argument "x"

[list_begin arguments]
[arg_def float x] Argument for the complementary error function
[list_end]

[para]

[call [cmd ::math::special::invnorm] [arg p]]

Compute the inverse of the normal distribution function for argument "p"

[list_begin arguments]
[arg_def float p] Argument for the inverse normal distribution function
(p must be greater than 0 and lower than 1)
[list_end]

[para]

[call [cmd ::math::special::J0] [arg x]]

Compute the zeroth-order Bessel function of the first kind for the
argument "x"

[list_begin arguments]
[arg_def float x] Argument for the Bessel function
[list_end]

[call [cmd ::math::special::J1] [arg x]]

Compute the first-order Bessel function of the first kind for the
argument "x"

[list_begin arguments]
[arg_def float x] Argument for the Bessel function
[list_end]

[call [cmd ::math::special::Jn] [arg n] [arg x]]

Compute the nth-order Bessel function of the first kind for the
argument "x"

[list_begin arguments]
[arg_def integer n] Order of the Bessel function
[arg_def float x] Argument for the Bessel function
[list_end]

[call [cmd ::math::special::J1/2] [arg x]]

Compute the half-order Bessel function of the first kind for the
argument "x"

[list_begin arguments]
[arg_def float x] Argument for the Bessel function
[list_end]

[call [cmd ::math::special::J-1/2] [arg x]]

Compute the minus-half-order Bessel function of the first kind for the
argument "x"

[list_begin arguments]
[arg_def float x] Argument for the Bessel function
[list_end]

[call [cmd ::math::special::I_n] [arg x]]

Compute the modified Bessel function of the first kind of order n for
the argument "x"

[list_begin arguments]
[arg_def int x] Positive integer order of the function
[arg_def float x] Argument for the function
[list_end]

[call [cmd ::math::special::cn] [arg u] [arg k]]

Compute the elliptic function [emph cn] for the argument "u" and
parameter "k".

[list_begin arguments]
[arg_def float u] Argument for the function
[arg_def float k] Parameter
[list_end]

[call [cmd ::math::special::dn] [arg u] [arg k]]

Compute the elliptic function [emph dn] for the argument "u" and
parameter "k".

[list_begin arguments]
[arg_def float u] Argument for the function
[arg_def float k] Parameter
[list_end]

[call [cmd ::math::special::sn] [arg u] [arg k]]

Compute the elliptic function [emph sn] for the argument "u" and
parameter "k".

[list_begin arguments]
[arg_def float u] Argument for the function
[arg_def float k] Parameter
[list_end]

[call [cmd ::math::special::elliptic_K] [arg k]]

Compute the complete elliptic integral of the first kind
for the argument "k"

[list_begin arguments]
[arg_def float k] Argument for the function
[list_end]

[call [cmd ::math::special::elliptic_E] [arg k]]

Compute the complete elliptic integral of the second kind
for the argument "k"

[list_begin arguments]
[arg_def float k] Argument for the function
[list_end]

[call [cmd ::math::special::exponential_Ei] [arg x]]

Compute the exponential integral of the second kind
for the argument "x"

[list_begin arguments]
[arg_def float x] Argument for the function (x != 0)
[list_end]

[call [cmd ::math::special::exponential_En] [arg n] [arg x]]

Compute the exponential integral of the first kind
for the argument "x" and order n

[list_begin arguments]
[arg_def int n] Order of the integral (n >= 0)
[arg_def float x] Argument for the function (x >= 0)
[list_end]

[call [cmd ::math::special::exponential_li] [arg x]]

Compute the logarithmic integral for the argument "x"

[list_begin arguments]
[arg_def float x] Argument for the function (x > 0)
[list_end]

[call [cmd ::math::special::exponential_Ci] [arg x]]

Compute the cosine integral for the argument "x"

[list_begin arguments]
[arg_def float x] Argument for the function (x > 0)
[list_end]

[call [cmd ::math::special::exponential_Si] [arg x]]

Compute the sine integral for the argument "x"

[list_begin arguments]
[arg_def float x] Argument for the function (x > 0)
[list_end]

[call [cmd ::math::special::exponential_Chi] [arg x]]

Compute the hyperbolic cosine integral for the argument "x"

[list_begin arguments]
[arg_def float x] Argument for the function (x > 0)
[list_end]

[call [cmd ::math::special::exponential_Shi] [arg x]]

Compute the hyperbolic sine integral for the argument "x"

[list_begin arguments]
[arg_def float x] Argument for the function (x > 0)
[list_end]

[call [cmd ::math::special::fresnel_C] [arg x]]

Compute the Fresnel cosine integral for real argument x

[list_begin arguments]
[arg_def float x] Argument for the function
[list_end]

[call [cmd ::math::special::fresnel_S] [arg x]]

Compute the Fresnel sine integral for real argument x

[list_begin arguments]
[arg_def float x] Argument for the function
[list_end]

[call [cmd ::math::special::sinc] [arg x]]

Compute the sinc function for real argument x

[list_begin arguments]
[arg_def float x] Argument for the function
[list_end]

[call [cmd ::math::special::legendre] [arg n]]

Return the Legendre polynomial of degree n
(see [sectref "THE ORTHOGONAL POLYNOMIALS"])

[list_begin arguments]
[arg_def int n] Degree of the polynomial
[list_end]

[para]

[call [cmd ::math::special::chebyshev] [arg n]]

Return the Chebyshev polynomial of degree n (of the first kind)

[list_begin arguments]
[arg_def int n] Degree of the polynomial
[list_end]

[para]

[call [cmd ::math::special::laguerre] [arg alpha] [arg n]]

Return the Laguerre polynomial of degree n with parameter alpha

[list_begin arguments]
[arg_def float alpha] Parameter of the Laguerre polynomial
[arg_def int n] Degree of the polynomial
[list_end]

[para]

[call [cmd ::math::special::hermite] [arg n]]

Return the Hermite polynomial of degree n

[list_begin arguments]
[arg_def int n] Degree of the polynomial
[list_end]

[para]

[list_end]

[section "THE ORTHOGONAL POLYNOMIALS"]

For dealing with the classical families of orthogonal polynomials, the
package relies on the [emph math::polynomials] package. To evaluate the
polynomial at some coordinate, use the [emph evalPolyn] command:
[example {
   set leg2 [::math::special::legendre 2]
   puts "Value at x=$x: [::math::polynomials::evalPolyn $leg2 $x]"
}]

[para]
The return value from the [emph legendre] and other commands is actually
the definition of the corresponding polynomial as used in that package.

[section "REMARKS ON THE IMPLEMENTATION"]

It should be noted, that the actual implementation of J0 and J1 depends
on straightforward Gaussian quadrature formulas. The (absolute) accuracy
of the results is of the order 1.0e-4 or better. The main reason to
implement them like that was that it was fast to do (the formulas are
simple) and the computations are fast too.

[para]
The implementation of J1/2 does not suffer from this: this function can
be expressed exactly in terms of elementary functions.

[para]
The functions J0 and J1 are the ones you will encounter most frequently
in practice.

[para]
The computation of I_n is based on Miller's algorithm for computing the
minimal function from recurrence relations.

[para]
The computation of the Gamma and Beta functions relies on the
combinatorics package, whereas that of the error functions relies on the
statistics package.

[para]
The computation of the complete elliptic integrals uses the AGM
algorithm.

[para]
Much information about these functions can be found in:
[para]
Abramowitz and Stegun: [emph "Handbook of Mathematical Functions"]
(Dover, ISBN 486-61272-4)

[vset CATEGORY {math :: special}]
[include ../common-text/feedback.inc]
[manpage_end]