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 |
class doubleIntegral: def __init__(self, xMin, xMax, yMin, yMax, func): self.xMin = xMin self.xMax = xMax self.yMin = yMin self.yMax = yMax self.func = func def trapezoid(self, nX, nY): result = 0 betweenX = (self.xMax - self.xMin) / nX betweenY = (self.yMax - self.yMin) / nY sums = [] for i in range(nY + 1): temp = 0 for j in range(nX): temp += betweenX / 2 * (self.func(self.xMin + j * betweenX, self.yMin + i * betweenY) + self.func(self.xMin + (j + 1) * betweenX, self.yMin + i * betweenY)) sums.append(temp) ans = 0 for i in range(nY): ans += betweenY / 2 * (sums[i] + sums[i + 1]) return ans def simpsons(self, nX, nY): result = 0 betweenX = (self.xMax - self.xMin) / nX betweenY = (self.yMax - self.yMin) / nY / 2 sums = [] for i in range(2 * nY + 1): temp = 0 for j in range(nX): a = self.xMin + j * betweenX b = self.xMin + (j + 1) * betweenX c = (a + b) / 2 temp += betweenX / 6 * (self.func(a, self.yMin + i * betweenY) + self.func( b, self.yMin + i * betweenY) + 4 * self.func(c, self.yMin + i * betweenY)) sums.append(temp) ans = 0 for i in range(nY): ans += betweenY / 3 * \ (sums[i * 2] + sums[(i + 1) * 2] + 4 * sums[i * 2 + 1]) return ans if __name__ == "__main__": # 个人不知道为什么梯形法求出来不对,有兴趣的可以自己去看看 # 输入格式:xMin, xMax, yMin, yMax, func # 输入格式': nx, ny print(doubleIntegral(0, 1, 0, 2, lambda x, y: x**2+2*y).trapezoid(4, 8)) print(doubleIntegral(0, 1, 0, 2, lambda x, y: x**2+2*y).simpsons(4, 8)) |