How to rewrite the digamma function to nonrecursive

Question

I ran a java program and got a stackoverflow error and this error was caused in a digamma function, the code is

public static double digamma(double x) {
    if (x >= 0 && x < GAMMA_MINX) {
        x = GAMMA_MINX;
    }
    if (x < DIGAMMA_MINNEGX) {
        return digamma(DIGAMMA_MINNEGX + GAMMA_MINX);
    }
    if (x > 0 && x <= S_LIMIT) {
        return -GAMMA - 1 / x;
    }

    if (x >= C_LIMIT) {
        double inv = 1 / (x * x);
        return Math.log(x) - 0.5 / x - inv
                * ((1.0 / 12) + inv * (1.0 / 120 - inv / 252));
    }
    return digamma(x + 1) - 1 / x;
}

the constant is defined as followed(the private static final is omitted for space consideration)

double GAMMA = 0.577215664901532860606512090082;
double GAMMA_MINX = 1.e-12;
double DIGAMMA_MINNEGX = -1250;
double C_LIMIT = 49;
double S_LIMIT = 1e-5;

As the stackoverflow error is caused by the level of the recursion is too deep so as to exceed the stack size of the thread(right?)so I increased the Xss parameter to 20M,but the error persists, so I think I have to rewrite it to nonrecursive, but I have no experience to rewite a recursive function to nonrecursive.

Answer 1

Untested but here it is:

public static double digamma(double x) {

    double value = 0;

    while (true){

        if (x >= 0 && x < GAMMA_MINX) {
            x = GAMMA_MINX;
        }
        if (x < DIGAMMA_MINNEGX) {
            x = DIGAMMA_MINNEGX + GAMMA_MINX;
            continue;
        }
        if (x > 0 && x <= S_LIMIT) {
            return value + -GAMMA - 1 / x;
        }

        if (x >= C_LIMIT) {
            double inv = 1 / (x * x);
            return value + Math.log(x) - 0.5 / x - inv
                    * ((1.0 / 12) + inv * (1.0 / 120 - inv / 252));
        }

        value -= 1 / x;
        x = x + 1;
    }

}

Since the code is almost tail-recursive, the trick is to throw a loop over the whole body.

The catch is the - 1 / x you have at the end. But since it's additive, you can just subtract 1 / x from the result prior to starting the next iteration.