Thanks! The problem came from the 1/2 ! For convenience I've found the function digamma on numpy

http://docs.scipy.org/doc/scipy/reference/generated/scipy.special.psi.html

But it's quite hard to find it! Maybe we can ask to add digamma in the title between parenthesis?


Fabien

2012/7/10 Damon McDougall <damon.mcdougall@gmail.com>
On Tue, Jul 10, 2012 at 08:57:24AM -0400, Benjamin Root wrote:
> On Tue, Jul 10, 2012 at 7:05 AM, Damon McDougall
> <damon.mcdougall@gmail.com>wrote:
>
> > On Tue, Jul 10, 2012 at 12:27:59PM +0200, Fabien Lafont wrote:
> >
> > > But It returns zero division error even when x is in ]0,1]
> >
> > I think it blows up at x = 0. What is the type of x in your usecase? Is
> > it an array? If x contains the element 0, you will get a zero
> > division error. You could try plotting the points explicitly:
> >
> Another problem might be the "1/2" part, which in python2.x would yield 0
> unless one does "from __future__ import division".
>
> Ben Root

Wow, I can't believe I didn't spot that. Nice one.

I will update my answer according to Ben's astute observation:

from scipy import special
from pylab import *

x = linspace(0.5, 2.0, num=100, endpoint=True)
y = special.polygamma(0, 0.5 + 1.0/x)
plot(x, y)
show()

Thanks Ben.

--
Damon McDougall
http://damon-is-a-geek.com
B2.39
Mathematics Institute
University of Warwick
Coventry
West Midlands
CV4 7AL
United Kingdom