code: public static double PowerSeries(double a, double b, double x) {
          double s, t, u, v, n, t1, z, ai;
  
          ai = 1.0 / a;
          u = (1.0 - b) * x;
          v = u / (a + 1.0);
          t1 = v;
          t = u;
          n = 2.0;
          s = 0.0;
          z = Constants.DoubleEpsilon * ai;
          while (Math.abs(v) > z) {
              u = (n - b) * x / n;
              t *= u;
              v = t / (a + n);
              s += v;
              n += 1.0;
          }
          s += t1;
          s += ai;
  
          u = a * Math.log(x);
          if ((a + b) < Gamma.GammaMax && Math.abs(u) < Constants.LogMax) {
              t = Gamma.Function(a + b) / (Gamma.Function(a) * Gamma.Function(b));
              s = s * t * Math.pow(x, a);
          } else {
              t = Gamma.Log(a + b) - Gamma.Log(a) - Gamma.Log(b) + u + Math.log(s);
              if (t < Constants.LogMin) s = 0.0;
              else s = Math.exp(t);
          }
          return s;
      }
docs: Power series for incomplete beta integral. Use when b*x is small and x not too close to 1.
  
  @param a Value.
  @param b Value.
  @param x Value.
  @return Result.