DEFAULT_ MAX_ ITERATIONS
The maximum number of terms permitted when summing the incomplete gamma series.
This is a runaway guard, not the thing that decides accuracy. The series' own convergence test terminates the loop, normally after a few dozen terms; this bound exists only so that an input which cannot converge fails instead of running forever.
The term count the series genuinely needs grows like the square root of the shape, so no constant is correct for every shape. It does not have to be. A million terms covers shapes up to roughly 1.0E10, while a gamma fitted to real data reaches a shape near 1.0E6 only when the sample's standard deviation is about a thousandth of its mean. Beyond the bound the routine throws, which is the intended behaviour for an input that far outside the range the library is meant to serve.
The previous value of 5000 was low enough to be reached by ordinary fitted distributions: a shape above about 1.1E6 evaluated near its mean exhausted it and threw.