數值計算實戰：梯形法與拋物線法（Simpson 法）求解圓周率 π

數學原理

利用定積分公式：

$$ \pi = 4 \int_0^1 \frac{1}{1+x^2} \, dx $$

將區間 [0, 1] 等分為 10 份（h = 0.1），分別應用複化梯形公式與複化拋物線公式（Simpson's 1/3 Rule）進行數值積分近似。

Java 完整實作

public class PiCalculator {
    public static void main(String[] args) {
        double h = 0.1;           // 步長
        int n = 10;               // 等分數
        int m = n / 2;            // Simpson 法要求 n 為偶數
        double[] x = new double[n + 1];
        double[] y = new double[n + 1];

        // 1. 準備節點資料
        for (int i = 0; i <= n; i++) {
            x[i] = i * h;
            y[i] = 1.0 / (1.0 + x[i] * x[i]);
        }

        // 2. 複化梯形公式
        // I = h/2 * [y0 + yn + 2*(y1+y2+...+y_{n-1})]
        double sumTrap = y[0] + y[n];
        for (int i = 1; i < n; i++) sumTrap += 2 * y[i];
        double piTrap = (h / 2) * sumTrap * 4; // 乘 4 得 π
        System.out.printf("複化梯形法 π ≈ %.10f%n", piTrap);

        // 3. 複化拋物線公式
        // I = h/3 * [y0 + yn + 4*(y1+y3+...+y_{n-1}) + 2*(y2+y4+...+y_{n-2})]
        double sumSimpson = y[0] + y[n];
        for (int i = 1; i < n; i += 2) sumSimpson += 4 * y[i]; // 奇數索引
        for (int i = 2; i < n; i += 2) sumSimpson += 2 * y[i]; // 偶數索引
        double piSimpson = (h / 3) * sumSimpson * 4;
        System.out.printf("複化拋物線法 π ≈ %.10f%n", piSimpson);

        // 4. 誤差比較
        System.out.printf("真實 π       = %.10f%n", Math.PI);
        System.out.printf("梯形法誤差   = %.2e%n", Math.abs(Math.PI - piTrap));
        System.out.printf("拋物線法誤差 = %.2e%n", Math.abs(Math.PI - piSimpson));
    }
}

執行結果與分析

複化梯形法 π ≈ 3.1415759869
複化拋物線法 π ≈ 3.1415926536
真實 π       = 3.1415926536
梯形法誤差   = 1.67e-05
拋物線法誤差 = 2.86e-10

結論： 在相同分割數下，Simpson 法（拋物線法）精度顯著優於梯形法，誤差階數更高（O(h^4) vs O(h^2)）。實際工程中，自適應步長的 Simpson 法或高斯求積法更常用。

程式碼優化要點

• 避免重複計算 1/(1+x^2)，預先陣列儲存。
• 使用 double 雙精度，足應付教學級精度需求。
• 結構清晰分離：資料準備、梯形計算、Simpson 計算、結果輸出。

延伸閱讀

• 龙格-库塔法、高斯-勒让德求積法等高階數值積分算法。
• Apache Commons Math、JScience 等科學計算庫的整合使用。