projects
/
pspp
/ blobdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
|
commitdiff
|
tree
raw
|
inline
| side by side
Fixed bug which crept in with the recent lexer changes.
[pspp]
/
lib
/
gsl-extras
/
betadistinv.c
diff --git
a/lib/gsl-extras/betadistinv.c
b/lib/gsl-extras/betadistinv.c
index a0838d27eb744b5ce54bddc051c294b264e55741..490914056d691e7fc7dd06dca94d299bde8db80f 100644
(file)
--- a/
lib/gsl-extras/betadistinv.c
+++ b/
lib/gsl-extras/betadistinv.c
@@
-1,6
+1,7
@@
/* cdf/betadistinv.c
*
/* cdf/betadistinv.c
*
- * Copyright (C) 2004 Jason H. Stover.
+ * Copyright (C) 2004 Free Software Foundation, Inc.
+ * Written by Jason H. Stover.
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
@@
-14,7
+15,7
@@
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
- * Foundation, Inc., 5
9 Temple Place, Suite 330, Boston, MA 02111-1307
, USA.
+ * Foundation, Inc., 5
1 Franklin Street, Fifth Floor, Boston, MA 02110-1301
, USA.
*/
/*
*/
/*
@@
-31,7
+32,7
@@
* Cornish-Fisher type," Annals of Mathematical Statistics, volume 39, number 8,
* August 1968, pages 1264-1273.
*/
* Cornish-Fisher type," Annals of Mathematical Statistics, volume 39, number 8,
* August 1968, pages 1264-1273.
*/
-#include <config.h>
+
#include <math.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_errno.h>
#include <math.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_errno.h>
@@
-40,7
+41,7
@@
#include <gsl/gsl_randist.h>
#include "gsl-extras.h"
#include <gsl/gsl_randist.h>
#include "gsl-extras.h"
-#define BETAINV_INIT_ERR .01
+#define BETAINV_INIT_ERR .0
0
1
#define BETADISTINV_N_TERMS 3
#define BETADISTINV_MAXITER 20
#define BETADISTINV_N_TERMS 3
#define BETADISTINV_MAXITER 20
@@
-389,7
+390,7
@@
gslextras_cdf_beta_Pinv ( const double p, const double a, const double b)
err = beta_result - p;
abserr = fabs(err);
relerr = abserr / p;
err = beta_result - p;
abserr = fabs(err);
relerr = abserr / p;
- while ( relerr > BETAINV_INIT_ERR
&& n_iter < 100
)
+ while ( relerr > BETAINV_INIT_ERR)
{
tmp = new_guess_P ( state, lower, upper,
p, a, b);
{
tmp = new_guess_P ( state, lower, upper,
p, a, b);
@@
-430,6
+431,19
@@
gslextras_cdf_beta_Pinv ( const double p, const double a, const double b)
result = state;
min_err = relerr;
}
result = state;
min_err = relerr;
}
+ else
+ {
+ /*
+ * Lagrange polynomial failed to reduce the error.
+ * This will happen with a very skewed beta density.
+ * Undo previous steps.
+ */
+ state = result;
+ beta_result = gsl_cdf_beta_P(state,a,b);
+ err = p - beta_result;
+ abserr = fabs(err);
+ relerr = abserr / p;
+ }
}
n_iter = 0;
}
n_iter = 0;
@@
-580,7
+594,7
@@
gslextras_cdf_beta_Qinv (double q, double a, double b)
err = beta_result - q;
abserr = fabs(err);
relerr = abserr / q;
err = beta_result - q;
abserr = fabs(err);
relerr = abserr / q;
- while ( relerr > BETAINV_INIT_ERR
&& n_iter < 100
)
+ while ( relerr > BETAINV_INIT_ERR)
{
n_iter++;
tmp = new_guess_Q ( state, lower, upper,
{
n_iter++;
tmp = new_guess_Q ( state, lower, upper,
@@
-621,6
+635,19
@@
gslextras_cdf_beta_Qinv (double q, double a, double b)
result = state;
min_err = relerr;
}
result = state;
min_err = relerr;
}
+ else
+ {
+ /*
+ * Lagrange polynomial failed to reduce the error.
+ * This will happen with a very skewed beta density.
+ * Undo previous steps.
+ */
+ state = result;
+ beta_result = gsl_cdf_beta_P(state,a,b);
+ err = q - beta_result;
+ abserr = fabs(err);
+ relerr = abserr / q;
+ }
}
/*
}
/*