ADMB Documentation  11.2.2828
 All Classes Files Functions Variables Typedefs Friends Defines
cbetacf.cpp
Go to the documentation of this file.
00001 
00007 #include <fvar.hpp>
00008 #include <math.h>
00009 #define EPS 1.0e-9
00010 #define FPMIN 1.0e-30
00011 
00022 double betacf(const double a, const double b, const double x, int MAXIT)
00023 {
00024   int m,m2;
00025   double aa,c,d,del,h,qab,qam,qap;
00026 
00027   qab=a+b;
00028   qap=a+1.0;
00029   qam=a-1.0;
00030   c=1.0;
00031   d=1.0-qab*x/qap;
00032   if (fabs(d) < FPMIN) d=FPMIN;
00033   d=1.0/d;
00034   h=d;
00035   for (m=1;m<=MAXIT;m++) {
00036     m2=2*m;
00037     aa=m*(b-m)*x/((qam+m2)*(a+m2));
00038     d=1.0+aa*d;
00039     if (fabs(d) < FPMIN) d=FPMIN;
00040     c=1.0+aa/c;
00041     if (fabs(c) < FPMIN) c=FPMIN;
00042     d=1.0/d;
00043     h *= d*c;
00044     aa = -(a+m)*(qab+m)*x/((a+m2)*(qap+m2));
00045     d=1.0+aa*d;
00046     if (fabs(d) < FPMIN) d=FPMIN;
00047     c=1.0+aa/c;
00048     if (fabs(c) < FPMIN) c=FPMIN;
00049     d=1.0/d;
00050     del=d*c;
00051     h *= del;
00052     if (fabs(del-1.0) < EPS) break;
00053   }
00054   if (m > MAXIT) cerr << "a or b too big, or MAXIT too small in betacf"
00055          << endl;
00056   return h;
00057 }
00058 #undef MAXIT
00059 #undef EPS
00060 #undef FPMIN