The R Project SVN R

Rev

Rev 42549 | Go to most recent revision | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 42549 Rev 50406
Line 37... Line 37...
37 
 
 37 
 
38 
#include "nmath.h"
 38 
#include "nmath.h"
39 
 
 39 
 
40 
double beta(double a, double b)
 40 
double beta(double a, double b)
41 
{
 41 
{
42 
    double val;
 - 
 
43 
 
 - 
 
44 
#ifdef NOMORE_FOR_THREADS
 42 
#ifdef NOMORE_FOR_THREADS
45 
    static double xmin, xmax = 0;/*-> typically = 171.61447887 for IEEE */
 43 
    static double xmin, xmax = 0;/*-> typically = 171.61447887 for IEEE */
46 
    static double lnsml = 0;/*-> typically = -708.3964185 */
 44 
    static double lnsml = 0;/*-> typically = -708.3964185 */
47 
 
 45 
 
48 
    if (xmax == 0) {
 46 
    if (xmax == 0) {
Line 74... Line 72...
74 
	return 0;
 72 
	return 0;
75 
    }
 73 
    }
76 
 
 74 
 
77 
    if (a + b < xmax) /* ~= 171.61 for IEEE */
 75 
    if (a + b < xmax) /* ~= 171.61 for IEEE */
78 
	return gammafn(a) * gammafn(b) / gammafn(a+b);
 76 
	return gammafn(a) * gammafn(b) / gammafn(a+b);
79 
 
 77 
    else {
80 
    val = lbeta(a, b);
 78 
	double val = lbeta(a, b);
81 
    if (val < lnsml) {
 79 
	if (val < lnsml) {
82 
	/* a and/or b so big that beta underflows */
 80 
	    /* a and/or b so big that beta underflows */
83 
	ML_ERROR(ME_UNDERFLOW, "beta");
 81 
	    ML_ERROR(ME_UNDERFLOW, "beta");
84 
	/* return ML_UNDERFLOW; pointless giving incorrect value */
 82 
	    /* return ML_UNDERFLOW; pointless giving incorrect value */
- 
 
 83 
	}
- 
 
 84 
	return exp(val);
85 
    }
 85 
    }
86 
    return exp(val);
 - 
 
87 
}
 86 
}