code: protected double quantile(double p, double xmin, double xmax, double eps) {
          if (eps <= 0.0) {
              throw new IllegalArgumentException("Invalid epsilon: " + eps);            
          }
          
          while (Math.abs(xmax - xmin) > eps) {
              double xmed = (xmax + xmin) / 2;
              if (cdf(xmed) > p) {
                  xmax = xmed;
              } else {
                  xmin = xmed;
              }
          }
  
          return xmin;
      }
docs: Inversion of CDF by bisection numeric root finding of "cdf(x) = p"
  for continuous distribution.