|
1 | 1 | """
|
2 | 2 | Numerical integration or quadrature for a smooth function f with known values at x_i
|
3 | 3 |
|
4 |
| -This method is the classical approach of suming 'Equally Spaced Abscissas' |
| 4 | +This method is the classical approach of summing 'Equally Spaced Abscissas' |
5 | 5 |
|
6 | 6 | method 2:
|
7 | 7 | "Simpson Rule"
|
8 | 8 |
|
9 | 9 | """
|
10 | 10 |
|
11 | 11 |
|
12 |
| -def method_2(boundary, steps): |
| 12 | +def method_2(boundary: list[int], steps: int) -> float: |
13 | 13 | # "Simpson Rule"
|
14 | 14 | # int(f) = delta_x/2 * (b-a)/3*(f1 + 4f2 + 2f_3 + ... + fn)
|
| 15 | + """ |
| 16 | + Calculate the definite integral of a function using Simpson's Rule. |
| 17 | + :param boundary: A list containing the lower and upper bounds of integration. |
| 18 | + :param steps: The number of steps or resolution for the integration. |
| 19 | + :return: The approximate integral value. |
| 20 | +
|
| 21 | + >>> round(method_2([0, 2, 4], 10), 10) |
| 22 | + 2.6666666667 |
| 23 | + >>> round(method_2([2, 0], 10), 10) |
| 24 | + -0.2666666667 |
| 25 | + >>> round(method_2([-2, -1], 10), 10) |
| 26 | + 2.172 |
| 27 | + >>> round(method_2([0, 1], 10), 10) |
| 28 | + 0.3333333333 |
| 29 | + >>> round(method_2([0, 2], 10), 10) |
| 30 | + 2.6666666667 |
| 31 | + >>> round(method_2([0, 2], 100), 10) |
| 32 | + 2.5621226667 |
| 33 | + >>> round(method_2([0, 1], 1000), 10) |
| 34 | + 0.3320026653 |
| 35 | + >>> round(method_2([0, 2], 0), 10) |
| 36 | + Traceback (most recent call last): |
| 37 | + ... |
| 38 | + ZeroDivisionError: Number of steps must be greater than zero |
| 39 | + >>> round(method_2([0, 2], -10), 10) |
| 40 | + Traceback (most recent call last): |
| 41 | + ... |
| 42 | + ZeroDivisionError: Number of steps must be greater than zero |
| 43 | + """ |
| 44 | + if steps <= 0: |
| 45 | + raise ZeroDivisionError("Number of steps must be greater than zero") |
| 46 | + |
15 | 47 | h = (boundary[1] - boundary[0]) / steps
|
16 | 48 | a = boundary[0]
|
17 | 49 | b = boundary[1]
|
@@ -41,11 +73,14 @@ def f(x): # enter your function here
|
41 | 73 | def main():
|
42 | 74 | a = 0.0 # Lower bound of integration
|
43 | 75 | b = 1.0 # Upper bound of integration
|
44 |
| - steps = 10.0 # define number of steps or resolution |
45 |
| - boundary = [a, b] # define boundary of integration |
| 76 | + steps = 10.0 # number of steps or resolution |
| 77 | + boundary = [a, b] # boundary of integration |
46 | 78 | y = method_2(boundary, steps)
|
47 | 79 | print(f"y = {y}")
|
48 | 80 |
|
49 | 81 |
|
50 | 82 | if __name__ == "__main__":
|
| 83 | + import doctest |
| 84 | + |
| 85 | + doctest.testmod() |
51 | 86 | main()
|
0 commit comments