# Some Python code for use at https://sagecell.sagemath.org/
**********
# Numerical Limits
def f(x):
return (x**2-1)/(x-1)
c = 1.0
for i in range(5):
e = (10.0)**(-i-1)
x1 = c + e
print( "x = ", x1, "f(x) = ", f(x1) )
x2 = c - e
print( "x = ", x2, "f(x) = ", f(x2) )
print()
**********
# Newton's Method
def f(x):
return x - cos(x)
def fp(x):
return 1.0 + sin(x)
x = 1.0
print(x)
for i in range(10):
x = x - f(x)/fp(x)
print(x)
**********
# Simple bisection
def f(x):
return x - cos(x)
a = 0.0 # f(a) and f(b) must be nonzero and have different signs
b = 1.0
for i in range(20):
c = ( a + b ) / 2
print(c)
x = f(c) * f(a)
if x == 0:
print( "Found it." )
break
elif x < 0:
b = c
else:
a = c
**********
# Trapzoid Rule
def f(x):
return 1.0/x
a = 1.0
b = 2.0
N = 8
h = ( b - a ) / (1.0*N)
x = float(a)
s = ( f(a) + f(b) ) / 2.0
for i in range(N-1):
x = x + h
s = s + f(x)
print( "Trapezoid rule approximation:", h * s)
**********
# Simpson's rule
def f(x):
return x**2 + x
a = 0.0
b = 3.0
N = 30 # Must be even
h = ( b - a ) / (1.0*N)
h2 = 2 * h
s0 = f(a) + f(b)
x = a - h
s1 = 0.0
for i in range(0,N,2):
x = x + h2
s1 = s1 + f(x)
s1 = 4 * s1
x = float(a)
s2 = 0.0
for i in range(1,N-1,2):
x = x + h2
s2 = s2 + f(x)
s2 = 2 * s2
s = h * ( s0 + s1 + s2 ) / 3.0
print( "Simpson's rule approximation:", s )