meyers hw 15.ipynb

Kernel: SageMath 8.2

#Chapter 5.2 #1

import numpy as np

def f1(x):
    return (x/np.sqrt(x**2+9))

x = np.array([2])
print('1 panel approximate answer is: ')
print((f1(4)-f1(0)/2*4))
print('2 panel approximate answer is: ')
x = np.linspace(0,4,2)
print(np.trapz(f1(x), x))
print('4 panel approximate answer is: ')
x = np.linspace(0,4,4)
print(np.trapz(f1(x), x))

1 panel approximate answer is: 
0.8
2 panel approximate answer is: 
1.6
4 panel approximate answer is: 
1.96066973984

## Chapter 5.2 #2 Composite Midpoint

import numpy as np

f1 = lambda x: (x/np.sqrt(x**2+9))

def midpoint(f, a, b, n):
    h = float(b-a)/n
    result = 0
    for i in range(n):
        result += f((a + h/2.0) + i*h)
    result *= h
    print(result)
midpoint(f1,0,4,1)
midpoint(f1,0,4,2)
midpoint(f1,0,4,4)

2.2188007849
2.04666909441
2.01105358484

from math import *
def f1(x):
    return (x/np.sqrt(x**2+9))

def Simpson(f, a, b, n):
    if a > b:
        print 'Incorrect bounds'
        return None
    if n%2 != 0: # also an 'if' because both tests are NOT
        # mutually exclusive
        print 'Invalid choice of n'
        return None
    else:
    	h = (b - a)/float(n) # need to cast 'n' as float in order to avoid
        # integer division
    	sum1 = 0
    	for i in range(1, n/2 + 1):
        	sum1 += f(a + (2*i - 1)*h)
    	sum1 *= 4
    	sum2 = 0
    	for i in range(1, n/2): # range must be ints: range() integer 
            #end argument expected, got float.
        	sum2 += f(a + 2*i*h)
    	sum2 *= 2
    	approx = (b - a)/(3.0*n)*(f(a) + f(b) + sum1 + sum2)
    	return approx
    
a=Simpson(f1,0,4,2)
b=Simpson(f1,0,4,4)

print(a)
print(b)

2.0125338566
2.00091286042

a-2 #error from 2 panel

0.012533856600610793

b-2 #error from 4 panel

0.00091286042133287282