Seite - 144 - in Programming for Computations – Python - A Gentle Introduction to Numerical Simulations with Python 3.6, Band Second Edition

144 6 ComputingIntegralsandTestingCode 6.3.2 AGeneral Implementation Wefollowtheadviceandlessonslearnedfromtheimplementationofthetrapezoidal method. Thus, we make a module midpoint.pywith a general implementation of (6.20) and a function application, just like we did with the trapezoidal function: def midpoint(f, a, b, n): h = (b-a)/n f_sum = 0 for i in range(0, n, 1): x = (a + h/2.0) + i*h f_sum = f_sum + f(x) return h*f_sum def application(): from math import exp v = lambda t: 3*(t**2)*exp(t**3) n = int(input(’n: ’)) numerical = midpoint(v, 0, 1, n) # Compare with exact result V = lambda t: exp(t**3) exact = V(1) - V(0) error = abs(exact - numerical) print(’n={:d}: {:.16f}, error: {:g}’.format(n, numerical, error)) if __name__ == ’__main__’: application() Inmidpoint, observe how thex values in the loop start out atx = a+ h2 (when i is 0), which is in the middle of the first rectangle. Thex values then increase by h foreach iteration,meaning thatwe repeatedly“jump” to the midpointof the next rectangleas i increases. This is consistent with the formula in (6.20), as is the final x value ofx= a+ h2 + (n−1)h. To convince yourself that the first, intermediate andfinalxvaluesarecorrect, lookatacasewithonly threerectangles, forexample. Whenapplication is called, the particular problem ∫1 0 3t 2et 3 dt is computed, i.e., the same integral that we handled with the trapezoidal method. Running the programwithn = 4gives the output n=4: 1.6189751378083810, error: 0.0993067 The magnitude of this error is now about 0.1 in contrast to 0.2, which we got with the trapezoidalrule.This is infactnotaccidental:onecanshowmathematically that the error of the midpoint method is a bit smaller than for the trapezoidal method. The differences are seldom of any practical importance, and on a laptop we can easily usen= 106 and get the answer with an error of about 10−12 in a couple of seconds.