diff --git a/doc/equations/digamma1.png b/doc/equations/digamma1.png index 295415c27..408c17be2 100644 Binary files a/doc/equations/digamma1.png and b/doc/equations/digamma1.png differ diff --git a/doc/equations/digamma1.svg b/doc/equations/digamma1.svg index 034d066d3..5feb7f95b 100644 --- a/doc/equations/digamma1.svg +++ b/doc/equations/digamma1.svg @@ -1,2 +1,2 @@ -ψ(x)=ddxln(Γ(x))=Γ(x)Γ(x) \ No newline at end of file +ψ(x)=ddxln(Γ(x))=Γ(x)Γ(x) \ No newline at end of file diff --git a/doc/equations/digamma2.svg b/doc/equations/digamma2.svg index b4db806c0..9b7ea8fd5 100644 --- a/doc/equations/digamma2.svg +++ b/doc/equations/digamma2.svg @@ -1,2 +1,2 @@ -ψ(x)=ln(x)+12xn=1B2n2nx2n \ No newline at end of file +ψ(x)=ln(x)+12xn=1B2n2nx2n \ No newline at end of file diff --git a/doc/equations/digamma3.png b/doc/equations/digamma3.png index cbf8f6e91..789e1021b 100644 Binary files a/doc/equations/digamma3.png and b/doc/equations/digamma3.png differ diff --git a/doc/equations/digamma3.svg b/doc/equations/digamma3.svg index 7ab0de08d..0200e6495 100644 --- a/doc/equations/digamma3.svg +++ b/doc/equations/digamma3.svg @@ -1,2 +1,2 @@ -ψ(x)=z12x+g12+ln(x+g12)+P(x)P(x)Q(x)Q(x)1 \ No newline at end of file +ψ(x)=z12x+g12+ln(x+g12)+P(x)P(x)Q(x)Q(x)1 \ No newline at end of file diff --git a/doc/equations/digamma4.mml b/doc/equations/digamma4.mml new file mode 100644 index 000000000..16e227469 --- /dev/null +++ b/doc/equations/digamma4.mml @@ -0,0 +1,47 @@ + +]> + +digamma4 + + + + + + ψ + + + n + + + = + + + + k + = + 1 + + + n + + 1 + + + + 1 + k + + + γ + + ; + + n + + + + + + \ No newline at end of file diff --git a/doc/equations/digamma4.png b/doc/equations/digamma4.png new file mode 100644 index 000000000..5de4146f7 Binary files /dev/null and b/doc/equations/digamma4.png differ diff --git a/doc/equations/digamma4.svg b/doc/equations/digamma4.svg new file mode 100644 index 000000000..05fac0afa --- /dev/null +++ b/doc/equations/digamma4.svg @@ -0,0 +1,2 @@ + +ψ(n)=k=1n11kγ;n \ No newline at end of file diff --git a/doc/equations/digamma5.mml b/doc/equations/digamma5.mml new file mode 100644 index 000000000..ad6c9b313 --- /dev/null +++ b/doc/equations/digamma5.mml @@ -0,0 +1,78 @@ + +]> + +digamma5 + + + + + + ψ + + + + 1 + 2 + + + n + + + = + + + + k + = + 1 + + + n + + 1 + + + + 1 + k + + + + + + + k + = + n + + + 2 + n + + 1 + + + + 2 + k + + + log + + + 4 + + + + γ + + ; + + n + + + + + + \ No newline at end of file diff --git a/doc/equations/digamma5.png b/doc/equations/digamma5.png new file mode 100644 index 000000000..79c42400f Binary files /dev/null and b/doc/equations/digamma5.png differ diff --git a/doc/equations/digamma5.svg b/doc/equations/digamma5.svg new file mode 100644 index 000000000..8509f8828 --- /dev/null +++ b/doc/equations/digamma5.svg @@ -0,0 +1,2 @@ + +ψ(12n)=k=1n11k+k=n2n12klog(4)γ;n \ No newline at end of file diff --git a/doc/equations/polygamma1.mml b/doc/equations/polygamma1.mml new file mode 100644 index 000000000..1b6e92ed8 --- /dev/null +++ b/doc/equations/polygamma1.mml @@ -0,0 +1,94 @@ + +]> + +polygamma1 + + + + + + + ψ + + + n + + + + + + x + + + = + + + + + 1 + + + + n + + + 1 + + + n + ! + + + + k + = + 0 + + + + + 1 + + + + + k + + + x + + + + n + + + 1 + + + + + = + + + + + n + + ψ + + + x + + + + + + + n + + x + + + + + + \ No newline at end of file diff --git a/doc/equations/polygamma1.png b/doc/equations/polygamma1.png new file mode 100644 index 000000000..8eaaca1ac Binary files /dev/null and b/doc/equations/polygamma1.png differ diff --git a/doc/equations/polygamma1.svg b/doc/equations/polygamma1.svg new file mode 100644 index 000000000..737a855aa --- /dev/null +++ b/doc/equations/polygamma1.svg @@ -0,0 +1,2 @@ + +ψ(n)(x)=(1)n+1n!k=01(k+x)n+1=nψ(x)nx \ No newline at end of file diff --git a/doc/equations/polygamma2.mml b/doc/equations/polygamma2.mml new file mode 100644 index 000000000..60c9744fd --- /dev/null +++ b/doc/equations/polygamma2.mml @@ -0,0 +1,90 @@ + +]> + +polygamma2 + + + + + + + ψ + + + + n + + + + + + + 1 + + x + + + = + + + + + 1 + + + n + + + ψ + + + + n + + + + + + + x + + + + + + + + + 1 + + + n + + π + + + + + n + + cot + + + π + x + + + + + + + x + n + + + + + + + \ No newline at end of file diff --git a/doc/equations/polygamma2.png b/doc/equations/polygamma2.png new file mode 100644 index 000000000..c91e5f274 Binary files /dev/null and b/doc/equations/polygamma2.png differ diff --git a/doc/equations/polygamma2.svg b/doc/equations/polygamma2.svg new file mode 100644 index 000000000..4a1c0d2c7 --- /dev/null +++ b/doc/equations/polygamma2.svg @@ -0,0 +1,2 @@ + +ψ(n)(1x)=(1)nψ(n)(x)+(1)nπncot(πx)xn \ No newline at end of file diff --git a/doc/equations/polygamma3.mml b/doc/equations/polygamma3.mml new file mode 100644 index 000000000..215ecf2bc --- /dev/null +++ b/doc/equations/polygamma3.mml @@ -0,0 +1,211 @@ + +]> + +polygamma3 + + + + + + + + + + + + + n + + cot + + + π + x + + + + + + + x + n + + + + + + = + + + + + + π + n + + + + + sin + + n + + + 1 + + + + + π + x + + + + + + + + k + = + 0 + + + + + n + + 2 + + 2 + + + + + C + + n + , + k + + + cos + + + + + 1 + + + 2 + k + + + π + x + + + + + ; + + + n + 2 + + + + + + + + + + + = + + + + + + π + n + + + + + sin + + n + + + 1 + + + + + π + x + + + + + + + + C + + n + , + 0 + + + + + + + + k + = + 1 + + + + n + 2 + + + + + C + + n + , + k + + + cos + + + 2 + k + π + x + + + + + + + ; + + + n + 2 + + + + + + + + + + \ No newline at end of file diff --git a/doc/equations/polygamma3.png b/doc/equations/polygamma3.png new file mode 100644 index 000000000..156285b33 Binary files /dev/null and b/doc/equations/polygamma3.png differ diff --git a/doc/equations/polygamma3.svg b/doc/equations/polygamma3.svg new file mode 100644 index 000000000..c364b7cf6 --- /dev/null +++ b/doc/equations/polygamma3.svg @@ -0,0 +1,2 @@ + +ncot(πx)xn=πnsinn+1(πx)k=0n22Cn,kcos((1+2k)πx);n2=πnsinn+1(πx)(Cn,0+k=1n2Cn,kcos(2kπx));n2 \ No newline at end of file diff --git a/doc/equations/polygamma4.mml b/doc/equations/polygamma4.mml new file mode 100644 index 000000000..2f7570b69 --- /dev/null +++ b/doc/equations/polygamma4.mml @@ -0,0 +1,113 @@ + +]> + +polygamma4 + + + + + + + ψ + + + + n + + + + + + + x + + + = + + + + + + + 1 + + + + n + + 1 + + + n + ! + + + x + + n + + + 1 + + + + + + + + + k + = + 0 + + + + + + + + + + 1 + + + + k + + + n + + + 1 + + + + + k + + + n + + + ! + ζ + + + k + + + n + + + 1 + + + + x + k + + + + k + ! + + + + + + \ No newline at end of file diff --git a/doc/equations/polygamma4.png b/doc/equations/polygamma4.png new file mode 100644 index 000000000..9f013bf8a Binary files /dev/null and b/doc/equations/polygamma4.png differ diff --git a/doc/equations/polygamma4.svg b/doc/equations/polygamma4.svg new file mode 100644 index 000000000..5a2860a2e --- /dev/null +++ b/doc/equations/polygamma4.svg @@ -0,0 +1,2 @@ + +ψ(n)(x)=(1)n1n!xn+1+k=0(1)k+n+1(k+n)!ζ(k+n+1)xkk! \ No newline at end of file diff --git a/doc/equations/polygamma5.mml b/doc/equations/polygamma5.mml new file mode 100644 index 000000000..1104b9560 --- /dev/null +++ b/doc/equations/polygamma5.mml @@ -0,0 +1,134 @@ + +]> + +polygamma5 + + + + + + + ψ + + + + n + + + + + + + x + + + + + + + + + + 1 + + + + n + + 1 + + + + + n + + 1 + + + ! + + + n + + + 2 + x + + + + + 2 + + x + + n + + + 1 + + + + + + + + + + 1 + + + n + + + + + k + = + 1 + + + + + + + + 2 + k + + + n + + 1 + + + ! + + + + + 2 + k + + + ! + + x + + 2 + k + + + n + + + + + + B + + 2 + k + + + + + + \ No newline at end of file diff --git a/doc/equations/polygamma5.png b/doc/equations/polygamma5.png new file mode 100644 index 000000000..aaf99373d Binary files /dev/null and b/doc/equations/polygamma5.png differ diff --git a/doc/equations/polygamma5.svg b/doc/equations/polygamma5.svg new file mode 100644 index 000000000..7bd89991e --- /dev/null +++ b/doc/equations/polygamma5.svg @@ -0,0 +1,2 @@ + +ψ(n)(x)(1)n1(n1)!(n+2x)2xn+1(1)nk=1(2k+n1)!(2k)!x2k+nB2k \ No newline at end of file diff --git a/doc/equations/polygamma6.mml b/doc/equations/polygamma6.mml new file mode 100644 index 000000000..47938c868 --- /dev/null +++ b/doc/equations/polygamma6.mml @@ -0,0 +1,88 @@ + +]> + +polygamma6 + + + + + + + ψ + + + + n + + + + + + + x + + m + + + = + + ψ + + + + n + + + + + + + x + + + + + + + + 1 + + + n + + n + ! + + + + k + = + 1 + + m + + + 1 + + + + + x + + k + + + + n + + + 1 + + + + + + + + \ No newline at end of file diff --git a/doc/equations/polygamma6.png b/doc/equations/polygamma6.png new file mode 100644 index 000000000..d634ae162 Binary files /dev/null and b/doc/equations/polygamma6.png differ diff --git a/doc/equations/polygamma6.svg b/doc/equations/polygamma6.svg new file mode 100644 index 000000000..0ce063c0a --- /dev/null +++ b/doc/equations/polygamma6.svg @@ -0,0 +1,2 @@ + +ψ(n)(xm)=ψ(n)(x)(1)nn!k=1m1(xk)n+1 \ No newline at end of file diff --git a/doc/equations/trigamma1.mml b/doc/equations/trigamma1.mml new file mode 100644 index 000000000..89d09624e --- /dev/null +++ b/doc/equations/trigamma1.mml @@ -0,0 +1,69 @@ + +]> + +trigamma1 + + + + + + + ψ + + + 1 + + + + + + x + + + = + + + + k + = + 0 + + + + + 1 + + + + + k + + + x + + + 2 + + + + = + + + + ψ + + + x + + + + + + x + + + + + + \ No newline at end of file diff --git a/doc/equations/trigamma1.png b/doc/equations/trigamma1.png new file mode 100644 index 000000000..9e966e099 Binary files /dev/null and b/doc/equations/trigamma1.png differ diff --git a/doc/equations/trigamma1.svg b/doc/equations/trigamma1.svg new file mode 100644 index 000000000..9ff41b44c --- /dev/null +++ b/doc/equations/trigamma1.svg @@ -0,0 +1,2 @@ + +ψ(1)(x)=k=01(k+x)2=ψ(x)x \ No newline at end of file diff --git a/doc/equations/trigamma2.mml b/doc/equations/trigamma2.mml new file mode 100644 index 000000000..37ac715c7 --- /dev/null +++ b/doc/equations/trigamma2.mml @@ -0,0 +1,65 @@ + +]> + +trigamma2 + + + + + + + ψ + + + 1 + + + + + + 1 + + x + + + = + + + + π + 2 + + + + + sin + 2 + + + + π + x + + + + + + + ψ + + + 1 + + + + + + x + + + + + + \ No newline at end of file diff --git a/doc/equations/trigamma2.png b/doc/equations/trigamma2.png new file mode 100644 index 000000000..45676ec18 Binary files /dev/null and b/doc/equations/trigamma2.png differ diff --git a/doc/equations/trigamma2.svg b/doc/equations/trigamma2.svg new file mode 100644 index 000000000..68f4d9510 --- /dev/null +++ b/doc/equations/trigamma2.svg @@ -0,0 +1,2 @@ + +ψ(1)(1x)=π2sin2(πx)ψ(1)(x) \ No newline at end of file diff --git a/doc/equations/trigamma3.mml b/doc/equations/trigamma3.mml new file mode 100644 index 000000000..934572d82 --- /dev/null +++ b/doc/equations/trigamma3.mml @@ -0,0 +1,54 @@ + +]> + +trigamma3 + + + + + + + ψ + + + 1 + + + + + + x + + + = + + ψ + + + 1 + + + + + + 1 + + + x + + + + + + 1 + + + x + 2 + + + + + + + \ No newline at end of file diff --git a/doc/equations/trigamma3.png b/doc/equations/trigamma3.png new file mode 100644 index 000000000..38aceee22 Binary files /dev/null and b/doc/equations/trigamma3.png differ diff --git a/doc/equations/trigamma3.svg b/doc/equations/trigamma3.svg new file mode 100644 index 000000000..121ab5517 --- /dev/null +++ b/doc/equations/trigamma3.svg @@ -0,0 +1,2 @@ + +ψ(1)(x)=ψ(1)(1+x)+1x2 \ No newline at end of file diff --git a/doc/equations/trigamma4.mml b/doc/equations/trigamma4.mml new file mode 100644 index 000000000..557d16d0e --- /dev/null +++ b/doc/equations/trigamma4.mml @@ -0,0 +1,49 @@ + +]> + +trigamma4 + + + + + + + ψ + + + 1 + + + + + + x + + + = + + + + + C + + + R + + + x + + + + + + + x + 2 + + + + + + \ No newline at end of file diff --git a/doc/equations/trigamma4.png b/doc/equations/trigamma4.png new file mode 100644 index 000000000..2a527a420 Binary files /dev/null and b/doc/equations/trigamma4.png differ diff --git a/doc/equations/trigamma4.svg b/doc/equations/trigamma4.svg new file mode 100644 index 000000000..1ccac5dd3 --- /dev/null +++ b/doc/equations/trigamma4.svg @@ -0,0 +1,2 @@ + +ψ(1)(x)=(C+R(x))x2 \ No newline at end of file diff --git a/doc/equations/trigamma5.mml b/doc/equations/trigamma5.mml new file mode 100644 index 000000000..4a8736435 --- /dev/null +++ b/doc/equations/trigamma5.mml @@ -0,0 +1,49 @@ + +]> + +trigamma5 + + + + + + + ψ + + + 1 + + + + + + x + + + = + + + + + 1 + + + R + + + + 1 + x + + + + + + + x + + + + + \ No newline at end of file diff --git a/doc/equations/trigamma5.png b/doc/equations/trigamma5.png new file mode 100644 index 000000000..671a4b324 Binary files /dev/null and b/doc/equations/trigamma5.png differ diff --git a/doc/equations/trigamma5.svg b/doc/equations/trigamma5.svg new file mode 100644 index 000000000..3813ab577 --- /dev/null +++ b/doc/equations/trigamma5.svg @@ -0,0 +1,2 @@ + +ψ(1)(x)=(1+R(1x))x \ No newline at end of file diff --git a/doc/html/index.html b/doc/html/index.html index dcdaa714c..65323f0af 100644 --- a/doc/html/index.html +++ b/doc/html/index.html @@ -111,7 +111,7 @@ This manual is also available in -

Last revised: October 05, 2014 at 07:41:09 GMT

+

Last revised: November 08, 2014 at 11:40:24 GMT


diff --git a/doc/html/indexes/s01.html b/doc/html/indexes/s01.html index ee4cf538c..ac6bdf5c0 100644 --- a/doc/html/indexes/s01.html +++ b/doc/html/indexes/s01.html @@ -24,7 +24,7 @@

-Function Index

+Function Index

A B C D E F G H I J K L M N O P Q R S T U V Z

@@ -1863,6 +1863,7 @@
  • Non-Member Properties

  • +
  • polygamma

  • powm1

  • prime

    @@ -2360,6 +2361,7 @@
  • Ratios of Gamma Functions

  • +
  • trigamma

  • trunc

    -Class Index

    +Class Index

    B C D E F G H I L M N O P Q R S T U W

    diff --git a/doc/html/indexes/s03.html b/doc/html/indexes/s03.html index 653c0b95f..0a428f995 100644 --- a/doc/html/indexes/s03.html +++ b/doc/html/indexes/s03.html @@ -24,7 +24,7 @@

    -Typedef Index

    +Typedef Index

    A B C D E F G H I L N O P R S T U V W

    diff --git a/doc/html/indexes/s04.html b/doc/html/indexes/s04.html index bf61ac720..2a9511664 100644 --- a/doc/html/indexes/s04.html +++ b/doc/html/indexes/s04.html @@ -24,7 +24,7 @@

    -Macro Index

    +Macro Index

    B F

    diff --git a/doc/html/indexes/s05.html b/doc/html/indexes/s05.html index aa562383a..b7d3c63b1 100644 --- a/doc/html/indexes/s05.html +++ b/doc/html/indexes/s05.html @@ -23,7 +23,7 @@

    -Index

    +Index

    A B C D E F G H I J K L M N O P Q R S T U V W Z

    @@ -4256,6 +4256,8 @@
  • Weibull Distribution

  • +
  • Polygamma

  • +
  • polygamma

  • Polynomial and Rational Function Evaluation

  • +
  • Trigamma

  • +
  • trigamma

  • trunc

    - +

    This documentation aims to use of the following naming and formatting conventions. diff --git a/doc/html/math_toolkit/navigation.html b/doc/html/math_toolkit/navigation.html index e05e1a7d3..cd43fa309 100644 --- a/doc/html/math_toolkit/navigation.html +++ b/doc/html/math_toolkit/navigation.html @@ -27,7 +27,7 @@ Navigation

    - +

    Boost.Math documentation is provided in both HTML and PDF formats. diff --git a/doc/html/math_toolkit/sf_gamma.html b/doc/html/math_toolkit/sf_gamma.html index 4392aa74d..050b2c6bf 100644 --- a/doc/html/math_toolkit/sf_gamma.html +++ b/doc/html/math_toolkit/sf_gamma.html @@ -30,6 +30,8 @@

    Gamma
    Log Gamma
    Digamma
    +
    Trigamma
    +
    Polygamma
    Ratios of Gamma Functions
    Incomplete Gamma Functions
    Incomplete Gamma Function diff --git a/doc/html/math_toolkit/sf_gamma/digamma.html b/doc/html/math_toolkit/sf_gamma/digamma.html index f4babd503..5e6a74a82 100644 --- a/doc/html/math_toolkit/sf_gamma/digamma.html +++ b/doc/html/math_toolkit/sf_gamma/digamma.html @@ -7,7 +7,7 @@ - + @@ -20,7 +20,7 @@

    -PrevUpHomeNext +PrevUpHomeNext

    @@ -62,11 +62,6 @@ what level of precision to use etc. Refer to the policy documentation for more details.

    -

    - There is no fully generic version of this function: all the implementations - are tuned to specific accuracy levels, the most precise of which delivers - 34-digits of precision. -

    The return type of this function is computed using the result type calculation rules: the result is of type double when T is an integer type, and type @@ -319,6 +314,17 @@ and BIG=20 for 128-bit reals allows the series to truncated after a suitably small number of terms and evaluated as a polynomial in 1/(x*x).

    +

    + The arbitrary precision version of this function uses recurrence relations + until x > BIG, and then evaluation via the asymtotic expansion above. + As special cases integer and half integer arguments are handled via: +

    +

    + +

    +

    + +

    The rational approximation devised by JM in the range [1,2] is derived as follows. @@ -379,7 +385,7 @@


    -PrevUpHomeNext +PrevUpHomeNext
    diff --git a/doc/html/math_toolkit/sf_gamma/gamma_ratios.html b/doc/html/math_toolkit/sf_gamma/gamma_ratios.html index df1b8f5f8..4acc20a6a 100644 --- a/doc/html/math_toolkit/sf_gamma/gamma_ratios.html +++ b/doc/html/math_toolkit/sf_gamma/gamma_ratios.html @@ -6,7 +6,7 @@ - + @@ -20,7 +20,7 @@
    -PrevUpHomeNext +PrevUpHomeNext

    @@ -349,7 +349,7 @@
    -PrevUpHomeNext +PrevUpHomeNext
    diff --git a/doc/html/math_toolkit/sf_gamma/polygamma.html b/doc/html/math_toolkit/sf_gamma/polygamma.html new file mode 100644 index 000000000..0123bddc8 --- /dev/null +++ b/doc/html/math_toolkit/sf_gamma/polygamma.html @@ -0,0 +1,269 @@ + + + +Polygamma + + + + + + + + + + + + + + + +
    Boost C++ LibrariesHomeLibrariesPeopleFAQMore
    +
    +
    +PrevUpHomeNext +
    +
    + +
    + + Synopsis +
    +
    #include <boost/math/special_functions/polygamma.hpp>
    +
    +
    namespace boost{ namespace math{
    +
    +template <class T>
    +calculated-result-type polygamma(int n, T z);
    +
    +template <class T, class Policy>
    +calculated-result-type polygamma(int n, T z, const Policy&);
    +
    +}} // namespaces
    +
    +
    + + Description +
    +

    + Returns the polygamma function of x. Polygamma is defined + as the n'th derivative of the digamma function: +

    +

    + +

    +

    + +

    +

    + The final Policy argument is optional and can + be used to control the behaviour of the function: how it handles errors, + what level of precision to use etc. Refer to the policy + documentation for more details. +

    +

    + The return type of this function is computed using the result + type calculation rules: the result is of type double when T is an integer type, and type + T otherwise. +

    +
    + + Accuracy +
    +

    + The following table shows the peak errors (in units of epsilon) found on + various platforms with various floating point types. Unless otherwise specified + any floating point type that is narrower than the one shown will have effectively zero error. +

    +
    ++++++ + + + + + + + + + + + + + + + + + + + + + + + + + + +
    +

    + Significand Size +

    +
    +

    + Platform and Compiler +

    +
    +

    + Small-medium positive arguments +

    +
    +

    + Small-medium negative x +

    +
    +

    + 53 +

    +
    +

    + Win32 Visual C++ 12 +

    +
    +

    + Peak=5.0 Mean=1 +

    +
    +

    + Peak=1200 Mean=65 +

    +
    +

    + 64 +

    +
    +

    + Win64 Mingw GCC +

    +
    +

    + Peak=16 Mean=3 +

    +
    +

    + Peak=33 Mean=3 +

    +
    +

    + 113 +

    +
    +

    + Win64 Mingw GCC __float128 +

    +
    +

    + Peak=6.5 Mean=1 +

    +
    +

    + Peak=30 Mean=4 +

    +
    +

    + As shown above, error rates are generally very acceptible for moderately + sized arguments. Error rates generally increase with increasing n + - this is particularly true for negative x and even + n - in this situation there is a root between each negative + integer and cancellation errors incured in computing the result increase + with increasing n. By the time n=170 + the errors are so bad we can no longer even tell the sign of the result at + double precision. It should + also be noted that for large n the function becomes + increasingly badly behaved when x is negative and is + very sensitive to slight changes in input. +

    +

    + For these reasons results should be treated with extreme + caution when n is large and x negative. +

    +
    + + Testing +
    +

    + Testing is against Mathematica generated spot values to 35 digit precision. +

    +
    + + Implementation +
    +

    + For x < 0 the following reflection formula is used: +

    +

    + +

    +

    + The n'th derivative of cot(x) is tabulated for small + n, and for larger n has the general form: +

    +

    + +

    +

    + The coefficients of the cosine terms can be calculated iteratively starting + from C1,0 = -1: see polygamma.hpp for the full details. +

    +

    + Once x is positive then we have two methods available to us, for small x + we use the series expansion: +

    +

    + +

    +

    + Note that as this is an alternating series we can accelerate the series using + Algorithm 1 from "Convergence Acceleration of Alternating + Series", Henri Cohen, Fernando Rodriguez Villegas, and Don Zagier, Experimental + Mathematics, 1999. Also note that the evaluation of zeta functions + at integer values is essentially a table lookup as zeta + is optimized for those cases. +

    +

    + For large x we use the asymptotic expansion: +

    +

    + +

    +

    + For x in-between the two extremes we use the relation: +

    +

    + +

    +

    + to make x large enough for the asymptotic expansion to be used. +

    +
    + + + +
    +
    +
    +PrevUpHomeNext +
    + + diff --git a/doc/html/math_toolkit/sf_gamma/trigamma.html b/doc/html/math_toolkit/sf_gamma/trigamma.html new file mode 100644 index 000000000..deb7b9446 --- /dev/null +++ b/doc/html/math_toolkit/sf_gamma/trigamma.html @@ -0,0 +1,217 @@ + + + +Trigamma + + + + + + + + + + + + + + + +
    Boost C++ LibrariesHomeLibrariesPeopleFAQMore
    +
    +
    +PrevUpHomeNext +
    +
    + +
    + + Synopsis +
    +
    #include <boost/math/special_functions/trigamma.hpp>
    +
    +
    namespace boost{ namespace math{
    +
    +template <class T>
    +calculated-result-type trigamma(T z);
    +
    +template <class T, class Policy>
    +calculated-result-type trigamma(T z, const Policy&);
    +
    +}} // namespaces
    +
    +
    + + Description +
    +

    + Returns the trigamma function of x. Trigamma is defined + as the derivative of the digamma function: +

    +

    + +

    +

    + +

    +

    + The final Policy argument is optional and can + be used to control the behaviour of the function: how it handles errors, + what level of precision to use etc. Refer to the policy + documentation for more details. +

    +

    + The return type of this function is computed using the result + type calculation rules: the result is of type double when T is an integer type, and type + T otherwise. +

    +
    + + Accuracy +
    +

    + The following table shows the peak errors (in units of epsilon) found on + various platforms with various floating point types. Unless otherwise specified + any floating point type that is narrower than the one shown will have effectively zero error. +

    +
    +++++ + + + + + + + + + + + + + + + + + + + + + + +
    +

    + Significand Size +

    +
    +

    + Platform and Compiler +

    +
    +

    + Random Values +

    +
    +

    + 53 +

    +
    +

    + Win32 Visual C++ 12 +

    +
    +

    + Peak=1.0 Mean=0.4 +

    +
    +

    + 64 +

    +
    +

    + Win64 Mingw GCC +

    +
    +

    + Peak=1.4 Mean=0.4 +

    +
    +

    + 113 +

    +
    +

    + Win64 Mingw GCC __float128 +

    +
    +

    + Peak=1.0 Mean=0.5 +

    +
    +

    + As shown above, error rates are generally very low for built in types. For + multiprecision types, error rates are typically in the order of a few epsilon. +

    +
    + + Testing +
    +

    + Testing is against Mathematica generated spot values to 35 digit precision. +

    +
    + + Implementation +
    +

    + The arbitary precision version of this function simply calls polygamma. +

    +

    + For built in fixed precision types, negative arguments are first made positive + via: +

    +

    + +

    +

    + Then arguments in the range [0, 1) are shifted to >= 1 via: +

    +

    + +

    +

    + Then evaluation is via one of a number of rational approximations, for small + x these are of the form: +

    +

    + +

    +

    + and for large x of the form: +

    +

    + +

    +
    + + + +
    +
    +
    +PrevUpHomeNext +
    + + diff --git a/doc/html/special.html b/doc/html/special.html index abd86044b..97b72391d 100644 --- a/doc/html/special.html +++ b/doc/html/special.html @@ -40,6 +40,8 @@
    Gamma
    Log Gamma
    Digamma
    +
    Trigamma
    +
    Polygamma
    Ratios of Gamma Functions
    Incomplete Gamma Functions
    Incomplete Gamma Function diff --git a/doc/math.qbk b/doc/math.qbk index 1f684e520..32ad0807d 100644 --- a/doc/math.qbk +++ b/doc/math.qbk @@ -99,6 +99,8 @@ and use the function's name as the link text.] [/gammas] [def __lgamma [link math_toolkit.sf_gamma.lgamma lgamma]] [def __digamma [link math_toolkit.sf_gamma.digamma digamma]] +[def __trigamma [link math_toolkit.sf_gamma.trigamma trigamma]] +[def __polygamma [link math_toolkit.sf_gamma.polygamma polygamma]] [def __tgamma_ratio [link math_toolkit.sf_gamma.gamma_ratios tgamma_ratio]] [def __tgamma_delta_ratio [link math_toolkit.sf_gamma.gamma_ratios tgamma_delta_ratio]] [def __tgamma [link math_toolkit.sf_gamma.tgamma tgamma]] @@ -464,6 +466,8 @@ and as a CD ISBN 0-9504833-2-X 978-0-9504833-2-0, Classification 519.2-dc22. [include sf/tgamma.qbk] [include sf/lgamma.qbk] [include sf/digamma.qbk] +[include sf/trigamma.qbk] +[include sf/polygamma.qbk] [include sf/gamma_ratios.qbk] [include sf/igamma.qbk] [include sf/igamma_inv.qbk] diff --git a/doc/sf/digamma.qbk b/doc/sf/digamma.qbk index cf3b9c232..b0c142bfd 100644 --- a/doc/sf/digamma.qbk +++ b/doc/sf/digamma.qbk @@ -27,10 +27,6 @@ derivative of the gamma function: [optional_policy] -There is no fully generic version of this function: all the implementations -are tuned to specific accuracy levels, the most precise of which delivers -34-digits of precision. - The return type of this function is computed using the __arg_pomotion_rules: the result is of type `double` when T is an integer type, and type T otherwise. @@ -97,6 +93,14 @@ Choosing BIG=10 for up to 80-bit reals, and BIG=20 for 128-bit reals allows the series to truncated after a suitably small number of terms and evaluated as a polynomial in `1/(x*x)`. +The arbitrary precision version of this function uses recurrence relations until +x > BIG, and then evaluation via the asymtotic expansion above. As special cases +integer and half integer arguments are handled via: + +[equation digamma4] + +[equation digamma5] + The rational approximation [jm_rationals] in the range [1,2] is derived as follows. First a high precision approximation to digamma was constructed using a 60-term diff --git a/doc/sf/polygamma.qbk b/doc/sf/polygamma.qbk new file mode 100644 index 000000000..c61f8582c --- /dev/null +++ b/doc/sf/polygamma.qbk @@ -0,0 +1,104 @@ +[section:polygamma Polygamma] + +[h4 Synopsis] + +`` +#include +`` + + namespace boost{ namespace math{ + + template + ``__sf_result`` polygamma(int n, T z); + + template + ``__sf_result`` polygamma(int n, T z, const ``__Policy``&); + + }} // namespaces + +[h4 Description] + +Returns the polygamma function of /x/. Polygamma is defined as the n'th +derivative of the digamma function: + +[equation polygamma1] + +[graph polygamma] + +[optional_policy] + +The return type of this function is computed using the __arg_pomotion_rules: +the result is of type `double` when T is an integer type, and type T otherwise. + +[h4 Accuracy] + +The following table shows the peak errors (in units of epsilon) +found on various platforms with various floating point types. +Unless otherwise specified any floating point type that is narrower +than the one shown will have __zero_error. + +[table +[[Significand Size] [Platform and Compiler] [Small-medium positive arguments] [Small-medium negative x] ] +[[53] [Win32 Visual C++ 12] [Peak=5.0 Mean=1] [Peak=1200 Mean=65]] +[[64] [Win64 Mingw GCC] [Peak=16 Mean=3] [Peak=33 Mean=3] ] +[[113] [Win64 Mingw GCC __float128] [Peak=6.5 Mean=1][Peak=30 Mean=4] ] +] + +As shown above, error rates are generally very acceptible for moderately sized +arguments. Error rates generally increase with increasing /n/ - this is particularly true +for negative /x/ and even /n/ - in this situation there is a root between each negative integer +and cancellation errors incured in computing the result increase with increasing /n/. By the time +/n=170/ the errors are so bad we can no longer even tell the sign of the result at `double` precision. +It should also be noted that for large /n/ the function becomes increasingly badly behaved +when /x/ is negative and is very sensitive to slight changes in input. + +[*For these reasons results should be treated with extreme caution when /n/ is large and x negative]. + +[h4 Testing] + +Testing is against Mathematica generated spot values to 35 digit precision. + +[h4 Implementation] + +For x < 0 the following reflection formula is used: + +[equation polygamma2] + +The n'th derivative of ['cot(x)] is tabulated for small /n/, and for larger n +has the general form: + +[equation polygamma3] + +The coefficients of the cosine terms can be calculated iteratively starting +from ['C[sub 1,0] = -1]: see polygamma.hpp for the full details. + +Once x is positive then we have two methods available to us, for small x +we use the series expansion: + +[equation polygamma4] + +Note that as this is an alternating series we can accelerate the series using +Algorithm 1 from ['"Convergence Acceleration of Alternating Series", +Henri Cohen, Fernando Rodriguez Villegas, and Don Zagier, Experimental Mathematics, 1999.] +Also note that the evaluation of zeta functions at integer values is essentially a table lookup +as __zeta is optimized for those cases. + +For large x we use the asymptotic expansion: + +[equation polygamma5] + +For x in-between the two extremes we use the relation: + +[equation polygamma6] + +to make x large enough for the asymptotic expansion to be used. + +[endsect][/section:polygamma The Polygamma Function] + +[/ + Copyright 2006 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). +] + diff --git a/doc/sf/trigamma.qbk b/doc/sf/trigamma.qbk new file mode 100644 index 000000000..e2514963c --- /dev/null +++ b/doc/sf/trigamma.qbk @@ -0,0 +1,84 @@ +[section:trigamma Trigamma] + +[h4 Synopsis] + +`` +#include +`` + + namespace boost{ namespace math{ + + template + ``__sf_result`` trigamma(T z); + + template + ``__sf_result`` trigamma(T z, const ``__Policy``&); + + }} // namespaces + +[h4 Description] + +Returns the trigamma function of /x/. Trigamma is defined as the +derivative of the digamma function: + +[equation trigamma1] + +[graph digamma] + +[optional_policy] + +The return type of this function is computed using the __arg_pomotion_rules: +the result is of type `double` when T is an integer type, and type T otherwise. + +[h4 Accuracy] + +The following table shows the peak errors (in units of epsilon) +found on various platforms with various floating point types. +Unless otherwise specified any floating point type that is narrower +than the one shown will have __zero_error. + +[table +[[Significand Size] [Platform and Compiler] [Random Values] ] +[[53] [Win32 Visual C++ 12] [Peak=1.0 Mean=0.4] ] +[[64] [Win64 Mingw GCC] [Peak=1.4 Mean=0.4] ] +[[113] [Win64 Mingw GCC __float128] [Peak=1.0 Mean=0.5] ] +] + +As shown above, error rates are generally very low for built in types. +For multiprecision types, error rates are typically in the order of a +few epsilon. + +[h4 Testing] + +Testing is against Mathematica generated spot values to 35 digit precision. + +[h4 Implementation] + +The arbitary precision version of this function simply calls __polygamma. + +For built in fixed precision types, negative arguments are first made positive via: + +[equation trigamma2] + +Then arguments in the range \[0, 1) are shifted to >= 1 via: + +[equation trigamma3] + +Then evaluation is via one of a number of rational approximations, for small x these are +of the form: + +[equation trigamma4] + +and for large x of the form: + +[equation trigamma5] + +[endsect][/section:digamma The Trigamma Function] + +[/ + Copyright 2006 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). +] +