import math
from scipy import integrate
import matplotlib.pyplot as plt


def erf_inner(x):
    return math.e ** (- x ** 2)


def solve0(a, b, N):
    """
    标准函数求积分
    """
    return integrate.quad(erf_inner, a, b)[0]


def solve1(a, b, N):
    """
    复化梯形积分公式
    """
    h = (b - a) / N
    s = erf_inner(a)
    for i in range(1, N):
        s += 2 * erf_inner(a + h * i)
    s += erf_inner(b)
    return s * h / 2


def solve2(a, b, N):
    """
    复化Simpson 1/3
    """
    h = (b - a) / N
    s = 0
    s = erf_inner(a)
    for i in range(1, N):
        s += 2 * erf_inner(a + h * i)
    for i in range(0, N):
        s += 4 * erf_inner(a + h * (i + 0.5))
    s += erf_inner(b)
    return s * h / 6


def solve3(a, b, N):
    """
    复化Simpson 3/8
    """
    h = (b - a) / N
    s = 0
    s = 7 * erf_inner(a)
    for i in range(1, N):
        s += 32 * erf_inner(a + h * (i - 0.75)) + 12 * erf_inner(a + h * (i - 0.5)) + 32 * erf_inner(
            a + h * (i - 0.25)) + 14 * erf_inner(a + h * i)
    s += 7 * erf_inner(b)
    return s * h / 90


def show():
    x = [i for i in range(60)]
    y = [[0 for j in range(60)] for _ in range(4)]
    solves = [solve0, solve1, solve2, solve3]
    for s in range(4):
        for i in x:
            y[s][i] = solves[s](0, i, 10)
        plt.plot(x, y[s], label=str(s))
    plt.legend(loc="lower right")
    plt.show()


def showN():
    x = [i for i in range(1, 30)]
    y = [[0 for j in range(30)] for _ in range(4)]
    solves = [solve0, solve1, solve2, solve3]
    for s in range(4):
        for i in x:
            y[s][i] = solves[s](0, 10, i)
        plt.plot(x, y[s][1:], label=str(s))
    plt.legend(loc="lower right")
    plt.show()


if __name__ == '__main__':
    # a = 0
    # b = 20
    # n = 90
    # print(solve1(a, b, n))
    # print(solve2(a, b, n))
    # print(solve3(a, b, n))
    # print(solve0(a, b, n))
    showN()