计算方法数值积分

本文最后更新于：2024年3月18日 凌晨

计算方法数值积分

梯形公式

算法

Newton-Cotes公式

算法

复合求积

算法

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
public class 复合梯形公式 {
    public static void main(String[] args) {
        double Y = Trapezoid(2, 3, 800);
        System.out.println("计算结果为：" + Y);
        System.out.println("精确解为：" + fx());
        System.out.println("绝对误差为：" + Math.abs(Y - fx()));
    }

    static double Trapezoid(double a, double b, double n) {
        double Y = 0;
        double h = (b - a) / n;
        Y += f(a) + f(b);
        for (int i = 1; i < n; i++) {
            Y += 2 * f(a + h * i);
        }
        Y *= h / 2;
        return Y;
    }

    static double f(double x) {
        double y = 1.0 / (x * x - 1);
        return y;
    }

    static double fx() {
        return (Math.log(2) - Math.log(3)) / -2;
    }
}

变步长求积

算法

龙贝格积分

算法

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
public class 龙贝格积分 {
    public static void main(String[] args) {
        double Y = romberg(2, 3, 5e-8);
        System.out.println("计算结果为：" + Y);
        System.out.println("精确解为：" + fx());
        System.out.println("绝对误差为：" + Math.abs(fx() - romberg(2, 3, 5e-5)));
    }

    static double f(double x) {
        double y = 1.0 / (x * x - 1);
        return y;
    }

    static double fx() {
        return (Math.log(2) - Math.log(3)) / -2;
    }

    static double romberg(double a, double b, double eps) {
        int n = 1, k;
        double h = b - a, x, temp;
        double T1, T2, S1 = 0, S2, C1 = 0, C2, R1 = 0, R2;
        T1 = (b - a) / 2 * (f(a) + f(b));
        while (true) {
            temp = 0;
            for (k = 0; k <= n - 1; k++) {
                x = a + k * h + h / 2;
                temp += f(x);
            }

            T2 = (T1 + temp * h) / 2;
            if (Math.abs(T2 - T1) < eps) return T2;
            S2 = T2 + (T2 - T1) / 3;
            if (n == 1) {
                T1 = T2;
                S1 = S2;
                h /= 2;
                n *= 2;
                continue;
            }
            C2 = S2 + (S2 - S1) / 15;
            if (n == 2) {
                C1 = C2;
                T1 = T2;
                S1 = S2;
                h /= 2;
                n *= 2;
                continue;
            }
            R2 = C2 + (C2 - C1) / 63;
            if (n == 4) {
                R1 = R2;
                C1 = C2;
                T1 = T2;
                S1 = S2;
                h /= 2;
                n *= 2;
                continue;
            }
            if (Math.abs(R2 - R1) < eps) return R2;
            R1 = R2;
            C1 = C2;
            T1 = T2;
            S1 = S2;
            h /= 2;
            n *= 2;
        }
    }
}

例题