Dear Per,

        thanks for the information.
I can reproduce your bug using python / pygsl and pure C (see
attachment). As pygsl is a pure wrapper (able to produce bugs of its
own, but not in this case), this bug should be reported to gsl.

If you want, I can do it for you, but given the fact that I lack the
required mathematical understanding, I would prefer if you could do it
yourself. The attached file reproduces the bug. The author of the code
was aware of the problem and his opinion it is outside the domain (the
error  code corresponds to sqrt(-1)). Info on the GSL bug mailling list
is to be found here:
        http://lists.gnu.org/mailman/listinfo/bug-gsl

Please drop me a line how you want to proceed.

Sincerely yours
        Pierre

> hi all,
> 
> i have been using the function for the pdf of the dirichlet distribution
> in PyGSL.
> 
> the function works correctly on most values, but appears to return nan
> when it shouldn't on more extreme cases. for example, the pdf evaluated
> on the set of values [0, 0, 1] with parameter settings [1/3, 1/3, 1/3]
> returns NaN, even though [0, 0, 1] is a perfectly fine value for the
> dirichlet pdf that should have non-zero probability. this is true for
> both the dirichlet_pdf and dirichlet_lnpdf (log of pdf) functions. 
> 
> in python notation, this is:
> 
> dirichlet_pdf([0, 0, 1], [1/3., 1/3., 1/3.])  (evaluates to NaN)
> dirichlet_lnpdf([0, 0, 1], [1/3., 1/3., 1/3.])  (evaluates to NaN as well)
> 
> any idea why this is or how i can fix it?
> 
> thank you.
> 
> 
> ------------------------------------------------------------------------
> 
> ------------------------------------------------------------------------------
> Enter the BlackBerry Developer Challenge  
> This is your chance to win up to $100,000 in prizes! For a limited time, 
> vendors submitting new applications to BlackBerry App World(TM) will have
> the opportunity to enter the BlackBerry Developer Challenge. See full prize  
> details at: http://p.sf.net/sfu/Challenge
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> pygsl-discuss mailing list
> [email protected]
> https://lists.sourceforge.net/lists/listinfo/pygsl-discuss

#include <stdio.h>

int 
main(int argc, char * argv[])
{
	const size_t K = 3;
	const double alpha[] = {  0.0,   0.0,   1.0};
	const double theta[] = {1./3., 1./3., 1./3.};
	double result;

	result = gsl_ran_dirichlet_pdf(K, alpha, theta);
	printf("result = %e\n", result);
	return 0;
}
------------------------------------------------------------------------------
Enter the BlackBerry Developer Challenge  
This is your chance to win up to $100,000 in prizes! For a limited time, 
vendors submitting new applications to BlackBerry App World(TM) will have
the opportunity to enter the BlackBerry Developer Challenge. See full prize  
details at: http://p.sf.net/sfu/Challenge
_______________________________________________
pygsl-discuss mailing list
[email protected]
https://lists.sourceforge.net/lists/listinfo/pygsl-discuss

Reply via email to