35d34
< 
45c44
< #define BETAINV_INIT_ERR .001
---
> #define BETAINV_INIT_ERR .01
394c393
<   while ( relerr > BETAINV_INIT_ERR)
---
>   while ( relerr > BETAINV_INIT_ERR && n_iter < 100)
435,447d433
<       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;
< 	}
598c584
<   while ( relerr > BETAINV_INIT_ERR)
---
>   while ( relerr > BETAINV_INIT_ERR && n_iter < 100)
639,651d624
<       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;
< 	}