import numpy as np

def bisectionMethod(f, a, b, tol=10**-15):
    assert f(a) * f(b) < 0, "Can only find root if initial endpoints lead to different sign when evaluating f."

    m = (a + b) / 2
    if b - a < tol: return m

    if np.abs(f(m)) < tol:
        return m

    if f(a) * f(m) < 0:
        return bisectionMethod(f, a, m)
    return bisectionMethod(f, m, b)

def secantMethod(f, x0, x1, tol=10**-15, maxiter=100):
    n = 1
    while np.abs(x0 - x1) >= tol and n < maxiter:
        n = n + 1
        assert np.abs(f(x1) - f(x0)) > tol, "Derivative of function too flat for secant method."
        x1, x0 = x1 - f(x1) * (x1 - x0) / (f(x1) - f(x0)), x1

    return x1

def newton(f, Df, x0, tol=10**-15, maxiter=100):
    n = 1
    while np.linalg.norm(f(x0)) > tol and n < maxiter:
        n = n + 1
        dx = np.linalg.solve(Df(x0), -f(x0))
        x0 = x0 + dx

    return x0

if __name__ == "__main__":
    # Approximation des goldenen Schnitts: hier wissen wir ja schon, dass es
    # zwei Nullstellen gibt p und -1/p, also kann ein Punkt x0 = 0 gewählt
    # werden. Die Parabel ist nach oben geöffnet, also muss lediglich x1 >> 0.
    f = lambda x: x**2 - x - 1
    z = bisectionMethod(f, 0, 10)
    print("The golden ratio is the unique positive root of x^2 - x - 1 and is approximately {}.".format(z))

    # Die Quadratwurzel aus 2 ausrechnen:
    f = lambda x: x**2 - 2
    z = secantMethod(f, 2, 3/2)
    print("Approximating sqrt(2) without relying on derivates: {}.".format(z))

    # Beispiel aus der Vorlesung auf S. 139
    f = lambda z: np.array([z[0]**2 + z[1]**2 - 6, z[0]**3 - z[1]**2])
    Df = lambda z: np.array([[2*z[0], 2*z[1]], [3*z[0]**2, -2*z[1]]])

    z = newton(f, Df, np.array([1, 1]))
    print("Found root {} of nonlinear system using Newton's method. Error: {}.".format(z, np.linalg.norm(f(z))))