nd4j-common/src/main/java/org/nd4j/linalg/util/BigDecimalMath.java
static public BigDecimal log(int n, final MathContext mc) {
        /* the value is undefined if x is negative.
         */
        if (n <= 0) {
            throw new ArithmeticException("Cannot take log of negative " + n);
        } else if (n == 1) {
            return BigDecimal.ZERO;
        } else if (n == 2) {
            if (mc.getPrecision() < LOG2.precision()) {
                return LOG2.round(mc);
            } else {
                /* Broadhurst \protect\vrule width0pt\protect\href{http://arxiv.org/abs/math/9803067}{arXiv:math/9803067}
                 * Error propagation: the error in log(2) is twice the error in S(2,-5,...).
                 */
                int[] a = {2, -5, -2, -7, -2, -5, 2, -3};
                BigDecimal S = broadhurstBBP(2, 1, a, new MathContext(1 + mc.getPrecision()));
                S = S.multiply(new BigDecimal(8));
                S = sqrt(divideRound(S, 3));
                return S.round(mc);
            }
        } else if (n == 3) {
            /* summation of a series roughly proportional to (7/500)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.013^k <= 10^(-precision), so k*log10(0.013) <= -precision
             * so k>= precision/1.87.
             */
            int kmax = (int) (mc.getPrecision() / 1.87);
            MathContext mcloc = new MathContext(mc.getPrecision() + 1 + (int) (Math.log10(kmax * 0.693 / 1.098)));
            BigDecimal log3 = multiplyRound(log(2, mcloc), 19);
            /* log3 is roughly 1, so absolute and relative error are the same. The
             * result will be divided by 12, so a conservative error is the one
             * already found in mc
             */
            double eps = prec2err(1.098, mc.getPrecision()) / kmax;
            Rational r = new Rational(7153, 524288);
            Rational pk = new Rational(7153, 524288);
            for (int k = 1;; k++) {
                Rational tmp = pk.divide(k);
                if (tmp.doubleValue() < eps) {
                    break;
                }
                /* how many digits of tmp do we need in the sum?
                 */
                mcloc = new MathContext(err2prec(tmp.doubleValue(), eps));
                BigDecimal c = pk.divide(k).BigDecimalValue(mcloc);
                if (k % 2 != 0) {
                    log3 = log3.add(c);
                } else {
                    log3 = log3.subtract(c);
                }
                pk = pk.multiply(r);
            }
            log3 = divideRound(log3, 12);
            return log3.round(mc);
        } else if (n == 5) {
            /* summation of a series roughly proportional to (7/160)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.046^k <= 10^(-precision), so k*log10(0.046) <= -precision
             * so k>= precision/1.33.
             */
            int kmax = (int) (mc.getPrecision() / 1.33);
            MathContext mcloc = new MathContext(mc.getPrecision() + 1 + (int) (Math.log10(kmax * 0.693 / 1.609)));
            BigDecimal log5 = multiplyRound(log(2, mcloc), 14);
            /* log5 is roughly 1.6, so absolute and relative error are the same. The
             * result will be divided by 6, so a conservative error is the one
             * already found in mc
             */
            double eps = prec2err(1.6, mc.getPrecision()) / kmax;
            Rational r = new Rational(759, 16384);
            Rational pk = new Rational(759, 16384);
            for (int k = 1;; k++) {
                Rational tmp = pk.divide(k);
                if (tmp.doubleValue() < eps) {
                    break;
                }
                /* how many digits of tmp do we need in the sum?
                 */
                mcloc = new MathContext(err2prec(tmp.doubleValue(), eps));
                BigDecimal c = pk.divide(k).BigDecimalValue(mcloc);
                log5 = log5.subtract(c);
                pk = pk.multiply(r);
            }
            log5 = divideRound(log5, 6);
            return log5.round(mc);
        } else if (n == 7) {
            /* summation of a series roughly proportional to (1/8)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.125^k <= 10^(-precision), so k*log10(0.125) <= -precision
             * so k>= precision/0.903.
             */
            int kmax = (int) (mc.getPrecision() / 0.903);
            MathContext mcloc = new MathContext(mc.getPrecision() + 1 + (int) (Math.log10(kmax * 3 * 0.693 / 1.098)));
            BigDecimal log7 = multiplyRound(log(2, mcloc), 3);
            /* log7 is roughly 1.9, so absolute and relative error are the same.
             */
            double eps = prec2err(1.9, mc.getPrecision()) / kmax;
            Rational r = new Rational(1, 8);
            Rational pk = new Rational(1, 8);
            for (int k = 1;; k++) {
                Rational tmp = pk.divide(k);
                if (tmp.doubleValue() < eps) {
                    break;
                }
                /* how many digits of tmp do we need in the sum?
                 */
                mcloc = new MathContext(err2prec(tmp.doubleValue(), eps));
                BigDecimal c = pk.divide(k).BigDecimalValue(mcloc);
                log7 = log7.subtract(c);
                pk = pk.multiply(r);
            }
            return log7.round(mc);
        } else {
            /* At this point one could either forward to the log(BigDecimal) signature (implemented)
             * or decompose n into Ifactors and use an implemenation of all the prime bases.
             * Estimate of the result; convert the mc argument to an absolute error eps
             * log(n+errn) = log(n)+errn/n = log(n)+eps
             */
            double res = Math.log((double) n);
            double eps = prec2err(res, mc.getPrecision());
            /* errn = eps*n, convert absolute error in result to requirement on absolute error in input
             */
            eps *= n;

            /* Convert this absolute requirement of error in n to a relative error in n
             */
            final MathContext mcloc = new MathContext(1 + err2prec((double) n, eps));
            /* Padd n with a number of zeros to trigger the required accuracy in
             * the standard signature method
             */
            BigDecimal nb = scalePrec(new BigDecimal(n), mcloc);
            return log(nb);
        }
    }
static public BigDecimal log(int n, final MathContext mc) {
        /* the value is undefined if x is negative.
         */
        if (n <= 0) {
            throw new ArithmeticException("Cannot take log of negative " + n);
        } else if (n == 1) {
            return BigDecimal.ZERO;
        } else if (n == 2) {
            if (mc.getPrecision() < LOG2.precision()) {
                return LOG2.round(mc);
            } else {
                /* Broadhurst \protect\vrule width0pt\protect\href{http://arxiv.org/abs/math/9803067}{arXiv:math/9803067}
                 * Error propagation: the error in log(2) is twice the error in S(2,-5,...).
                 */
                int[] a = {2, -5, -2, -7, -2, -5, 2, -3};
                BigDecimal S = broadhurstBBP(2, 1, a, new MathContext(1 + mc.getPrecision()));
                S = S.multiply(new BigDecimal(8));
                S = sqrt(divideRound(S, 3));
                return S.round(mc);
            }
        } else if (n == 3) {
            /* summation of a series roughly proportional to (7/500)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.013^k <= 10^(-precision), so k*log10(0.013) <= -precision
             * so k>= precision/1.87.
             */
            int kmax = (int) (mc.getPrecision() / 1.87);
            MathContext mcloc = new MathContext(mc.getPrecision() + 1 + (int) (Math.log10(kmax * 0.693 / 1.098)));
            BigDecimal log3 = multiplyRound(log(2, mcloc), 19);
            /* log3 is roughly 1, so absolute and relative error are the same. The
             * result will be divided by 12, so a conservative error is the one
             * already found in mc
             */
            double eps = prec2err(1.098, mc.getPrecision()) / kmax;
            Rational r = new Rational(7153, 524288);
            Rational pk = new Rational(7153, 524288);
            for (int k = 1;; k++) {
                Rational tmp = pk.divide(k);
                if (tmp.doubleValue() < eps) {
                    break;
                }
                /* how many digits of tmp do we need in the sum?
                 */
                mcloc = new MathContext(err2prec(tmp.doubleValue(), eps));
                BigDecimal c = pk.divide(k).BigDecimalValue(mcloc);
                if (k % 2 != 0) {
                    log3 = log3.add(c);
                } else {
                    log3 = log3.subtract(c);
                }
                pk = pk.multiply(r);
            }
            log3 = divideRound(log3, 12);
            return log3.round(mc);
        } else if (n == 5) {
            /* summation of a series roughly proportional to (7/160)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.046^k <= 10^(-precision), so k*log10(0.046) <= -precision
             * so k>= precision/1.33.
             */
            int kmax = (int) (mc.getPrecision() / 1.33);
            MathContext mcloc = new MathContext(mc.getPrecision() + 1 + (int) (Math.log10(kmax * 0.693 / 1.609)));
            BigDecimal log5 = multiplyRound(log(2, mcloc), 14);
            /* log5 is roughly 1.6, so absolute and relative error are the same. The
             * result will be divided by 6, so a conservative error is the one
             * already found in mc
             */
            double eps = prec2err(1.6, mc.getPrecision()) / kmax;
            Rational r = new Rational(759, 16384);
            Rational pk = new Rational(759, 16384);
            for (int k = 1;; k++) {
                Rational tmp = pk.divide(k);
                if (tmp.doubleValue() < eps) {
                    break;
                }
                /* how many digits of tmp do we need in the sum?
                 */
                mcloc = new MathContext(err2prec(tmp.doubleValue(), eps));
                BigDecimal c = pk.divide(k).BigDecimalValue(mcloc);
                log5 = log5.subtract(c);
                pk = pk.multiply(r);
            }
            log5 = divideRound(log5, 6);
            return log5.round(mc);
        } else if (n == 7) {
            /* summation of a series roughly proportional to (1/8)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.125^k <= 10^(-precision), so k*log10(0.125) <= -precision
             * so k>= precision/0.903.
             */
            int kmax = (int) (mc.getPrecision() / 0.903);
            MathContext mcloc = new MathContext(mc.getPrecision() + 1 + (int) (Math.log10(kmax * 3 * 0.693 / 1.098)));
            BigDecimal log7 = multiplyRound(log(2, mcloc), 3);
            /* log7 is roughly 1.9, so absolute and relative error are the same.
             */
            double eps = prec2err(1.9, mc.getPrecision()) / kmax;
            Rational r = new Rational(1, 8);
            Rational pk = new Rational(1, 8);
            for (int k = 1;; k++) {
                Rational tmp = pk.divide(k);
                if (tmp.doubleValue() < eps) {
                    break;
                }
                /* how many digits of tmp do we need in the sum?
                 */
                mcloc = new MathContext(err2prec(tmp.doubleValue(), eps));
                BigDecimal c = pk.divide(k).BigDecimalValue(mcloc);
                log7 = log7.subtract(c);
                pk = pk.multiply(r);
            }
            return log7.round(mc);
        } else {
            /* At this point one could either forward to the log(BigDecimal) signature (implemented)
             * or decompose n into Ifactors and use an implemenation of all the prime bases.
             * Estimate of the result; convert the mc argument to an absolute error eps
             * log(n+errn) = log(n)+errn/n = log(n)+eps
             */
            double res = Math.log((double) n);
            double eps = prec2err(res, mc.getPrecision());
            /* errn = eps*n, convert absolute error in result to requirement on absolute error in input
             */
            eps *= n;

            /* Convert this absolute requirement of error in n to a relative error in n
             */
            final MathContext mcloc = new MathContext(1 + err2prec((double) n, eps));
            /* Padd n with a number of zeros to trigger the required accuracy in
             * the standard signature method
             */
            BigDecimal nb = scalePrec(new BigDecimal(n), mcloc);
            return log(nb);
        }
    }
static public BigDecimal log ( int n , final MathContext mc ) { /* the value is undefined if x is negative.
         */ if ( n <= 0 ) { throw new ArithmeticException ( "Cannot take log of negative " + n ) ; } else if ( n == 1 ) { return BigDecimal . ZERO ; } else if ( n == 2 ) { if ( mc . getPrecision ( ) < LOG2 . precision ( ) ) { return LOG2 . round ( mc ) ; } else { /* Broadhurst \protect\vrule width0pt\protect\href{http://arxiv.org/abs/math/9803067}{arXiv:math/9803067}
                 * Error propagation: the error in log(2) is twice the error in S(2,-5,...).
                 */ int [ ] a = { 2 , - 5 , - 2 , - 7 , - 2 , - 5 , 2 , - 3 } ; BigDecimal S = broadhurstBBP ( 2 , 1 , a , new MathContext ( 1 + mc . getPrecision ( ) ) ) ; S = S . multiply ( new BigDecimal ( 8 ) ) ; S = sqrt ( divideRound ( S , 3 ) ) ; return S . round ( mc ) ; } } else if ( n == 3 ) { /* summation of a series roughly proportional to (7/500)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.013^k <= 10^(-precision), so k*log10(0.013) <= -precision
             * so k>= precision/1.87.
             */ int kmax = ( int ) ( mc . getPrecision ( ) / 1.87 ) ; MathContext mcloc = new MathContext ( mc . getPrecision ( ) + 1 + ( int ) ( Math . log10 ( kmax * 0.693 / 1.098 ) ) ) ; BigDecimal log3 = multiplyRound ( log ( 2 , mcloc ) , 19 ) ; /* log3 is roughly 1, so absolute and relative error are the same. The
             * result will be divided by 12, so a conservative error is the one
             * already found in mc
             */ double eps = prec2err ( 1.098 , mc . getPrecision ( ) ) / kmax ; Rational r = new Rational ( 7153 , 524288 ) ; Rational pk = new Rational ( 7153 , 524288 ) ; for ( int k = 1 ; ; k ++ ) { Rational tmp = pk . divide ( k ) ; if ( tmp . doubleValue ( ) < eps ) { break ; } /* how many digits of tmp do we need in the sum?
                 */ mcloc = new MathContext ( err2prec ( tmp . doubleValue ( ) , eps ) ) ; BigDecimal c = pk . divide ( k ) . BigDecimalValue ( mcloc ) ; if ( k % 2 != 0 ) { log3 = log3 . add ( c ) ; } else { log3 = log3 . subtract ( c ) ; } pk = pk . multiply ( r ) ; } log3 = divideRound ( log3 , 12 ) ; return log3 . round ( mc ) ; } else if ( n == 5 ) { /* summation of a series roughly proportional to (7/160)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.046^k <= 10^(-precision), so k*log10(0.046) <= -precision
             * so k>= precision/1.33.
             */ int kmax = ( int ) ( mc . getPrecision ( ) / 1.33 ) ; MathContext mcloc = new MathContext ( mc . getPrecision ( ) + 1 + ( int ) ( Math . log10 ( kmax * 0.693 / 1.609 ) ) ) ; BigDecimal log5 = multiplyRound ( log ( 2 , mcloc ) , 14 ) ; /* log5 is roughly 1.6, so absolute and relative error are the same. The
             * result will be divided by 6, so a conservative error is the one
             * already found in mc
             */ double eps = prec2err ( 1.6 , mc . getPrecision ( ) ) / kmax ; Rational r = new Rational ( 759 , 16384 ) ; Rational pk = new Rational ( 759 , 16384 ) ; for ( int k = 1 ; ; k ++ ) { Rational tmp = pk . divide ( k ) ; if ( tmp . doubleValue ( ) < eps ) { break ; } /* how many digits of tmp do we need in the sum?
                 */ mcloc = new MathContext ( err2prec ( tmp . doubleValue ( ) , eps ) ) ; BigDecimal c = pk . divide ( k ) . BigDecimalValue ( mcloc ) ; log5 = log5 . subtract ( c ) ; pk = pk . multiply ( r ) ; } log5 = divideRound ( log5 , 6 ) ; return log5 . round ( mc ) ; } else if ( n == 7 ) { /* summation of a series roughly proportional to (1/8)^k. Estimate count
             * of terms to estimate the precision (drop the favorable additional
             * 1/k here): 0.125^k <= 10^(-precision), so k*log10(0.125) <= -precision
             * so k>= precision/0.903.
             */ int kmax = ( int ) ( mc . getPrecision ( ) / 0.903 ) ; MathContext mcloc = new MathContext ( mc . getPrecision ( ) + 1 + ( int ) ( Math . log10 ( kmax * 3 * 0.693 / 1.098 ) ) ) ; BigDecimal log7 = multiplyRound ( log ( 2 , mcloc ) , 3 ) ; /* log7 is roughly 1.9, so absolute and relative error are the same.
             */ double eps = prec2err ( 1.9 , mc . getPrecision ( ) ) / kmax ; Rational r = new Rational ( 1 , 8 ) ; Rational pk = new Rational ( 1 , 8 ) ; for ( int k = 1 ; ; k ++ ) { Rational tmp = pk . divide ( k ) ; if ( tmp . doubleValue ( ) < eps ) { break ; } /* how many digits of tmp do we need in the sum?
                 */ mcloc = new MathContext ( err2prec ( tmp . doubleValue ( ) , eps ) ) ; BigDecimal c = pk . divide ( k ) . BigDecimalValue ( mcloc ) ; log7 = log7 . subtract ( c ) ; pk = pk . multiply ( r ) ; } return log7 . round ( mc ) ; } else { /* At this point one could either forward to the log(BigDecimal) signature (implemented)
             * or decompose n into Ifactors and use an implemenation of all the prime bases.
             * Estimate of the result; convert the mc argument to an absolute error eps
             * log(n+errn) = log(n)+errn/n = log(n)+eps
             */ double res = Math . log ( ( double ) n ) ; double eps = prec2err ( res , mc . getPrecision ( ) ) ; /* errn = eps*n, convert absolute error in result to requirement on absolute error in input
             */ eps *= n ; /* Convert this absolute requirement of error in n to a relative error in n
             */ final MathContext mcloc = new MathContext ( 1 + err2prec ( ( double ) n , eps ) ) ; /* Padd n with a number of zeros to trigger the required accuracy in
             * the standard signature method
             */ BigDecimal nb = scalePrec ( new BigDecimal ( n ) , mcloc ) ; return log ( nb ) ; } }
The natural logarithm.

@param n  The main argument, a strictly positive integer.
@param mc The requirements on the precision.
@return ln(n).
The natural logarithm .