function atinfinityplus( n, x ) {
	var partTerm; // Value of current term excluding the Bernoulli number part
	var xsquared;
	var term; // Value of current term to be added to sum
	var sum; // Current value of accumulated sum
	var nlx;
	var k2;
	var k;

	if ( n+x === x ) {
		// If `x` is very large, just concentrate on the first part of the expression and use logs:
		if ( n === 1 ) {
			return 1.0 / x;
		}
		nlx = n * ln( x );
		if ( nlx < MAX_LN && n < MAX_FACTORIAL ) {
			return ( (n & 1) ? 1.0 : -1.0 ) * factorial( n-1 ) * pow( x, -n );
		}
		return ( (n & 1) ? 1.0 : -1.0 ) * exp( gammaln( n ) - ( n*ln(x) ) );
	}
	xsquared = x * x;

	// Start by setting `partTerm` to `(n-1)! / x^(n+1)`, which is common to both the first term of the series (with k = 1) and to the leading part. We can then get to the leading term by: `partTerm * (n + 2 * x) / 2` and to the first term in the series (excluding the Bernoulli number) by: `partTerm n * (n + 1) / (2x)`. If either the factorial would over- or the power term underflow, set `partTerm` to 0 and then we know that we have to use logs for the initial terms:
	if ( n > MAX_FACTORIAL && n*n > MAX_LN ) {
		partTerm = 0.0;
	} else {
		partTerm = factorial( n-1 ) * pow( x, -n-1 );
	}
	if ( partTerm === 0.0 ) {
		// Either `n` is very large, or the power term underflows. Set the initial values of `partTerm`, `term`, and `sum` via logs:
		partTerm = gammaln(n) - ( (n+1) * ln(x) );
		sum = exp( partTerm + ln( n + (2.0*x) ) - LN_TWO );
		partTerm += ln( n*(n+1) ) - LN_TWO - ln(x);
		partTerm = exp( partTerm );
	} else {
		sum = partTerm * ( n+(2.0*x) ) / 2.0;
		partTerm *= ( n*(n+1) ) / 2.0;
		partTerm /= x;
	}
	// If the leading term is 0, so is the result:
	if ( sum === 0.0 ) {
		return sum;
	}
	for ( k = 1; ; ) {
		term = partTerm * bernoulli( k*2 );
		sum += term;

		// Normal termination condition:
		if ( abs( term/sum ) < EPS ) {
			break;
		}

		// Increment our counter, and move `partTerm` on to the next value:
		k += 1;
		k2 = 2 * k;
		partTerm *= ( n+k2-2 ) * ( n-1+k2 );
		partTerm /= ( k2-1 ) * k2;
		partTerm /= xsquared;
		if ( k > MAX_SERIES_ITERATIONS ) {
			debug( 'Series did not converge, closest value was: %d.', sum );
			return NaN;
		}
	}
	if ( ( n-1 ) & 1 ) {
		sum = -sum;
	}
	return sum;
}
Evaluates the polygamma function for large values of `x` such as for `x > 400`.

@private
@param {PositiveInteger} n - derivative to evaluate
@param {number} x - input
@returns {number} (n+1)'th derivative
@see {@link http://functions.wolfram.com/GammaBetaErf/PolyGamma2/06/02/0001/}