使用 Chudnovsky 算法计算 pi 的错误结果

Wrong result computing pi with the Chudnovsky algorithm

我正在尝试使用 Chudnovsky algorithm 计算圆周率。我做了一些研究并使用了 BigDecimals。但是当我 运行 我的代码时,它给了我错误的值。我检查了我的代码很多次,但没有发现任何改进。

public class PI {
    private static void Chudnovsky(int numSigDigits) {
        long time1 = System.nanoTime();
        int numDigits = numSigDigits + 10;
        MathContext mc = new MathContext(numDigits,RoundingMode.HALF_EVEN), mcx=new MathContext(numDigits*20, RoundingMode.HALF_EVEN);

        BigDecimal ratC = new BigDecimal(426880L, mc),
                irratC = new BigDecimal(10005L, mcx),
                incrL = new BigDecimal(545140134L, mc),
                l = new BigDecimal(13591409L, mcx),
                mulX = new BigDecimal(-262537412640768000L, mc),
                x = new BigDecimal(1L, mcx),
                incrK = new BigDecimal(12L, mc),
                k = new BigDecimal(6L, mc),
                m = new BigDecimal(1L, mc),
                b16 = new BigDecimal(16L, mc), j;

        BigDecimal C = ratC.add(irratC.sqrt(mcx), mcx);
        BigDecimal sum = new BigDecimal(0, mcx);

        for (int i = 0; i < numDigits;i++) {
            j = new BigDecimal(i + 1, mc);

            sum = sum.add(m.multiply(l, mc).divide(x, mc), mc);

            l = l.add(incrL, mcx);
            x = x.multiply(mulX, mcx);
            m = m.multiply(k.pow(3, mc).subtract(b16.multiply(k, mc), mc).divide(j.pow(3, mc), mc), mc);
            k = k.add(incrK, mc);
        }
        
        long time2 = System.nanoTime();
        System.out.println((double) (time2 - time1) / 1000000000);

        System.out.println(C.divide(sum, mcx).setScale(numSigDigits, RoundingMode.HALF_EVEN).toEngineeringString());
    }

    public static void main(String[] args) {
        Chudnovsky(15);
        System.out.println(Math.PI);
        //System.out.println(str);
    }
}

输出:

0.270327757  //time
0.031415434926348  //my value of pi
3.141592653589793  //Math.PI

我从 https://en.wikipedia.org/wiki/Chudnovsky_algorithm

中获取了这个算法

你应该改变

BigDecimal C = ratC.add(irratC.sqrt(mcx), mcx);

BigDecimal C = ratC.multiply(irratC.sqrt(mcx), mcx);