cs107-lecture-examples

Example codes used during Harvard CS107 lectures
git clone https://git.0xfab.ch/cs107-lecture-examples.git
Log | Files | Refs | README | LICENSE

newton.py (1467B)


      1 #!/usr/bin/env python3
      2 # File       : newton.py
      3 # Description: Newton's method with exact and FD approximated Jacobian
      4 # Copyright 2022 Harvard University. All Rights Reserved.
      5 import cProfile
      6 import numpy as np
      7 
      8 f = lambda x: x - np.exp(-2.0 * np.sin(4.0 * x) * np.sin(4.0 * x))
      9 J = lambda x: 1.0 + 16.0 * np.exp(-2.0 * np.sin(4.0 * x)**2) * np.sin(4.0 * x) * np.cos(4.0 * x)
     10 
     11 def J_fd(eps):
     12     """Returns finite-difference function object for given discretization `eps`."""
     13     return lambda x: (f(x + eps) - f(x)) / eps
     14 
     15 
     16 def newton(f, J, x_k, tol=1.0e-8, max_it=100):
     17     """Newton-Raphson method."""
     18     root = None
     19     for k in range(max_it):
     20         dx_k = -f(x_k) / J(x_k)
     21         if abs(dx_k) < tol:
     22             root = x_k + dx_k
     23             print(f"Found root {root:e} at iteration {k+1}")
     24             break
     25         print(f"Iteration {k+1}: Delta x = {dx_k:e}")
     26         x_k += dx_k
     27     return root
     28 
     29 
     30 def main(J):
     31     """
     32     Main function with specific Jacobian `J`.
     33 
     34     Parameters
     35     ----------
     36     J : function
     37         Jacobian function object to be used in the Newton-Raphson method.
     38 
     39     """
     40     newton(f, J, 0.1, 1.0e-8, 100)
     41 
     42 
     43 if __name__ == "__main__":
     44     # TODO: implement the three cProfile calls here:
     45     # 1. Profile the exact Jacobian and save it in file `exact`
     46     # 2. Profile the FD Jacobian with eps=1.0e-1 and save it in file `approx_coarse`
     47     # 3. Profile the FD Jacobian with eps=1.0e-8 and save it in file `approx_fine`
     48     pass