octave-antonio: libinterp/corefcn/psi.cc annotate

annotate libinterp/corefcn/psi.cc @ 20161:65e22ba879f0

psi: add support to compute the polygamma function (kth-derivative). * libinterp/corefcn/psi.cc: previously, only the digamma function, k == 0, was being computed. Add support for polygamma function, add tests, and improve documentation. * liboctave/cruft/slatec-fn/dpsifn.f, liboctave/cruft/slatec-fn/psifn.f: the two functions that actually compute the the polygamma functions, copied verbatim from SLATEC, and under public domain. * liboctave/cruft/slatec-fn/module.mk: add dpsifn.f and psifn.f to the build system. * liboctave/numeric/lo-specfun.cc: add new signature for function psi to compute polygamma function that wraps the Fortran DPSIFN and PSIFN functions. * liboctave/numeric/lo-specfun.h: declare new function and document all psi() with doxygen.

author	Carnë Draug <carandraug@octave.org>
date	Sun, 03 May 2015 22:52:07 +0100
parents	bd565f3e0ecb
children

rev	line source
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	1 /*
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	2
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	3 Copyright (C) 2015 Carnë Draug
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	4
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	5 This file is part of Octave.
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	6
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	7 Octave is free software; you can redistribute it and/or modify it
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	8 under the terms of the GNU General Public License as published by the
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	9 Free Software Foundation; either version 3 of the License, or (at your
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	10 option) any later version.
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	11
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	12 Octave is distributed in the hope that it will be useful, but WITHOUT
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	13 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	15 for more details.
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	16
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	17 You should have received a copy of the GNU General Public License
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	18 along with Octave; see the file COPYING. If not, see
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	19 <http://www.gnu.org/licenses/>.
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	20
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	21 */
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	22
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	23 #ifdef HAVE_CONFIG_H
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	24 #include <config.h>
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	25 #endif
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	26
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	27 #include "ov.h"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	28 #include "defun.h"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	29 #include "error.h"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	30 #include "dNDArray.h"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	31 #include "fNDArray.h"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	32
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	33 #include "lo-specfun.h"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	34
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	35 DEFUN (psi, args, ,
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	36 "-- texinfo --\n\
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	37 @deftypefn {Function File} {} psi (@var{z})\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	38 @deftypefnx {Function File} {} psi (@var{k}, @var{z})\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	39 Compute the psi (polygamma) function.\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	40 \n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	41 The polygamma functions are the @var{k}th derivative of the logarithm\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	42 of the gamma function. If unspecified, @var{k} defaults to zero. A value\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	43 of zero computes the digamma function, a value of 1, the trigamma function,\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	44 and so on.\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	45 \n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	46 The digamma function is defined:\n\
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	47 \n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	48 @tex\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	49 $$\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	50 \\Psi (z) = {d (log (\\Gamma (z))) \\over dx}\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	51 $$\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	52 @end tex\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	53 @ifnottex\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	54 @example\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	55 @group\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	56 psi (z) = d (log (gamma (z))) / dx\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	57 @end group\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	58 @end example\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	59 @end ifnottex\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	60 \n\
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	61 When computing the digamma function (when @var{k} equals zero), @var{z}\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	62 can have any value real or complex value. However, for polygamma functions\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	63 (@var{k} higher than 0), @var{z} must be real and non-negative.\n\
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	64 \n\
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	65 @seealso{gamma, gammainc, gammaln}\n\
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	66 @end deftypefn")
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	67 {
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	68 octave_value retval;
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	69
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	70 const octave_idx_type nargin = args.length ();
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	71 if (nargin < 1 \|\| nargin > 2)
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	72 {
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	73 print_usage ();
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	74 return retval;
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	75 }
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	76
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	77 const octave_value oct_z = (nargin == 1) ? args(0) : args(1);
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	78 const octave_idx_type k = (nargin == 1) ? 0 : args(0).idx_type_value ();
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	79 if (error_state \|\| k < 0)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	80 {
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	81 error ("psi: K must be a non-negative integer");
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	82 return retval;
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	83 }
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	84 else if (k == 0)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	85 {
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	86 #define FLOAT_BRANCH(T, A, M, E) \
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	87 if (oct_z.is_ ## T ##_type ()) \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	88 { \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	89 const A ## NDArray z = oct_z.M ## array_value (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	90 A ## NDArray psi_z (z.dims ()); \
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	91 \
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	92 const E* zv = z.data (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	93 E* psi_zv = psi_z.fortran_vec (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	94 const octave_idx_type n = z.numel (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	95 for (octave_idx_type i = 0; i < n; i++) \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	96 psi_zv++ = psi (zv++); \
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	97 \
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	98 retval = psi_z; \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	99 }
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	100
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	101 if (oct_z.is_complex_type ())
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	102 {
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	103 FLOAT_BRANCH(double, Complex, complex_, Complex)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	104 else FLOAT_BRANCH(single, FloatComplex, float_complex_, FloatComplex)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	105 else
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	106 {
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	107 error ("psi: Z must be a floating point");
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	108 }
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	109 }
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	110 else
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	111 {
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	112 FLOAT_BRANCH(double, , , double)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	113 else FLOAT_BRANCH(single, Float, float_, float)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	114 else
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	115 {
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	116 error ("psi: Z must be a floating point");
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	117 }
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	118 }
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	119
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	120 #undef FLOAT_BRANCH
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	121 }
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	122 else
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	123 {
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	124 if (! oct_z.is_real_type ())
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	125 {
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	126 error ("psi: Z must be real value for polygamma (K > 0)");
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	127 return retval;
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	128 }
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	129
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	130 #define FLOAT_BRANCH(T, A, M, E) \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	131 if (oct_z.is_ ## T ##_type ()) \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	132 { \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	133 const A ## NDArray z = oct_z.M ## array_value (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	134 A ## NDArray psi_z (z.dims ()); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	135 \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	136 const E* zv = z.data (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	137 E* psi_zv = psi_z.fortran_vec (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	138 const octave_idx_type n = z.numel (); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	139 for (octave_idx_type i = 0; i < n; i++) \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	140 { \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	141 if (*zv < 0) \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	142 { \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	143 error ("psi: Z must be non-negative for polygamma (K > 0)"); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	144 return retval; \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	145 } \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	146 psi_zv++ = psi (k, zv++); \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	147 } \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	148 retval = psi_z; \
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	149 }
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	150
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	151 FLOAT_BRANCH(double, , , double)
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	152 else FLOAT_BRANCH(single, Float, float_, float)
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	153 else
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	154 {
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	155 error ("psi: Z must be a floating point for polygamma (K > 0)");
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	156 }
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	157
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	158 #undef FLOAT_BRANCH
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	159 }
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	160
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	161 return retval;
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	162 }
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	163
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	164 /*
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	165 %!shared em
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	166 %! em = 0.577215664901532860606512090082402431042; # Euler-Mascheroni Constant
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	167
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	168 %!assert (psi (ones (7, 3, 5)), repmat (-em, [7 3 5]))
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	169 %!assert (psi ([0 1]), [-Inf -em])
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	170 %!assert (psi ([-20:1]), [repmat(-Inf, [1 21]) -em])
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	171 %!assert (psi (single ([0 1])), single ([-Inf -em]))
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	172
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	173 ## Abramowitz and Stegun, page 258, eq 6.3.5
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	174 %!test
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	175 %! z = [-10:.1:-.1 .1:.1:20]; # drop the 0
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	176 %! assert (psi (z + 1), psi (z) + 1 ./ z, eps*1000)
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	177
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	178 ## Abramowitz and Stegun, page 258, eq 6.3.2
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	179 %!assert (psi (1), -em)
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	180
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	181 ## Abramowitz and Stegun, page 258, eq 6.3.3
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	182 %!assert (psi (1/2), -em - 2 * log (2))
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	183
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	184 ## The following tests are from Pascal Sebah and Xavier Gourdon (2002)
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	185 ## "Introduction to the Gamma Function"
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	186
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	187 ## Interesting identities of the digamma function, in section of 5.1.3
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	188 %!assert (psi (1/3), - em - (3/2) * log(3) - ((sqrt (3) / 6) * pi), eps*10)
20156 bd565f3e0ecb psi: refactor to reduce code duplication. Carnë Draug <carandraug@octave.org> parents: 20155 diff changeset	189 %!assert (psi (1/4), - em -3 * log (2) - pi /2)
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	190 %!assert (psi (1/6), - em -2 * log (2) - (3/2) * log (3) - ((sqrt (3) / 2) * pi), eps*10)
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	191
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	192 ## First 6 zeros of the digamma function, in section of 5.1.5 (and also on
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	193 ## Abramowitz and Stegun, page 258, eq 6.3.19)
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	194 %!assert (psi ( 1.46163214496836234126265954232572132846819620400644), 0, eps)
20156 bd565f3e0ecb psi: refactor to reduce code duplication. Carnë Draug <carandraug@octave.org> parents: 20155 diff changeset	195 %!assert (psi (-0.504083008264455409258269304533302498955385182368579), 0, eps)
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	196 %!assert (psi (-1.573498473162390458778286043690434612655040859116846), 0, eps)
45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	197 %!assert (psi (-2.610720868444144650001537715718724207951074010873480), 0, eps*10)
20156 bd565f3e0ecb psi: refactor to reduce code duplication. Carnë Draug <carandraug@octave.org> parents: 20155 diff changeset	198 %!assert (psi (-3.635293366436901097839181566946017713948423861193530), 0, eps*10)
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	199 %!assert (psi (-4.653237761743142441714598151148207363719069416133868), 0, eps*100)
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	200
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	201 ## Tests for complex values
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	202 %!shared z
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	203 %! z = [-10:.1:-.1 .1:.1:20]; # drop the 0
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	204
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	205 ## Abramowitz and Stegun, page 259 eq 6.3.10
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	206 %!assert (real (psi (iz)), real (psi (1 - iz)))
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	207
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	208 ## Abramowitz and Stegun, page 259 eq 6.3.11
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	209 %!assert (imag (psi (iz)), 1/2 . 1./z + 1/2 * pi * coth (pi * z), eps *10)
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	210
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	211 ## Abramowitz and Stegun, page 259 eq 6.3.12
20156 bd565f3e0ecb psi: refactor to reduce code duplication. Carnë Draug <carandraug@octave.org> parents: 20155 diff changeset	212 %!assert (imag (psi (1/2 + iz)), 1/2 pi * tanh (pi * z), eps)
20155 1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	213
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	214 ## Abramowitz and Stegun, page 259 eq 6.3.13
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	215 %!assert (imag (psi (1 + iz)), - 1./(2z) + 1/2 * pi * coth (pi * z), eps*10)
1fae49e34a1a psi: add support for complex numbers. Carnë Draug <carandraug@octave.org> parents: 20154 diff changeset	216
20161 65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	217 ## Abramowitz and Stegun, page 260 eq 6.4.5
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	218 %!test
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	219 %! for z = 0:20
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	220 %! assert (psi (1, z + 0.5), 0.5 * (pi^2) - 4 * sum ((2(1:z) -1) .^(-2)), eps10)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	221 %! endfor
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	222
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	223 ## Abramowitz and Stegun, page 260 eq 6.4.6
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	224 %!test
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	225 %! z = 0.1:0.1:20;
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	226 %! for n = 0:8
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	227 %! ## our precision goes down really quick when computing n is too high,
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	228 %! assert (psi (n, z+1), psi (n, z) + ((-1)^n) * factorial (n) * (z.^(-n-1)), 0.1)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	229 %! endfor
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	230
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	231 ## Test input validation
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	232 %!error psi ()
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	233 %!error psi (1, 2, 3)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	234 %!error <Z must be> psi ("non numeric")
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	235 %!error <K must be a non-negative integer> psi (-5, 1)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	236 %!error <Z must be non-negative for polygamma> psi (5, -1)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	237 %!error <Z must be a floating point> psi (5, uint8 (-1))
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	238 %!error <Z must be real value for polygamma> psi (5, 5i)
65e22ba879f0 psi: add support to compute the polygamma function (kth-derivative). Carnë Draug <carandraug@octave.org> parents: 20156 diff changeset	239
20154 45565ecec019 New function psi to compute the digamma function. Carnë Draug <carandraug@octave.org> parents: diff changeset	240 */

Mercurial > octave-antonio

annotate libinterp/corefcn/psi.cc @ 20161:65e22ba879f0