# 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 )