[NumCS] Add code (Zeros)

semester3/numcs/parts/03_zeros/03_bisection-method.tex

@@ -4,10 +4,32 @@ Die Idee hier ist, das Intervall immer weiter zu halbieren und ein bekannterer N
 
 \innumpy haben wir \texttt{scipy.optimize.bisect} und \texttt{scipy.optimize.fsolve}, wobei \texttt{fsolve} ein alter Algorithmus ist.
 
-Im Skript auf Seiten 206 - 207 findet sich eine manuelle implementation des Bisektionsverfahren.
-Der Code ist jedoch (at the time of writing) nicht ausführbar aufgrund von \texttt{IndentationErrors}
-
 Das Bisektionsverfahren konvergiert linear und kann nur für Funktionen verwenden, bei welchen die Nullstellen auf beiden Seiten jeweils ungleiche Vorzeichen haben.
 
 Für jeden Iterationsschritt ermitteln wir die Mitte des Intervalls und berechnen die Funktionswerte an den Rändern, wie auch dem Mittelpunkt.
 Dann ersetzen wir den Rand des Intervalls, dessen Funktionswert dasselbe Vorzeichen hat, wie der Funktionswert des Mittelpunkts.
+
+\innumpy Lässt sich das Intervallhalbierungsverfahren (Bisektion) leicht direkt implementieren:
+
+\begin{code}{python}
+def bisection_method(f, a: float, b: float, tol=1e-12, maxIter=100):
+    """ Bisection method on f using initial bounds a,b """
+    if a > b: a, b = b, a
+    fa = f(a); fb = f(b)
+    if fa*fb > 0: raise ValueError("f(a) & f(b) must have different signs for bisection")
+    
+    sgn_fb = 1
+    if fa > 0: sgn_fb = -1
+    x = 0.5 * (a + b)
+    
+    iter = 1
+    while (b-a > tol and a<x and x<b and iter<maxIter):
+        # Check on which side f(x) is, update bounds
+        if sgn_fb*f(x) > 0: b = x
+        else:               a = x
+
+        x = 0.5 * (a + b)
+        iter += 1
+        
+    return x, iter
+\end{code}

semester3/numcs/parts/03_zeros/04_newton-one-d.tex

@@ -15,3 +15,17 @@ Beim Newtonverfahren verwendet man für jeden Iterationsschritt die lineare Funk
 
 \drmvspace
 falls $F(x^*) = 0$ und $F^(x^*) \neq 0$
+
+\newpage
+
+\innumpy Ist das Newton-Verfahren mit sehr wenig code implementierbar:
+
+\begin{code}{python}
+def newton_method(f, df, x: float, tol=1e-12, maxIter=100):
+    """ Use Newton's method to find zeros, requires derivative of f """
+    iter = 0
+    while (np.abs(df(x)) > tol and iter < maxIter):
+        x -= f(x)/df(x); iter += 1
+    
+    return x, iter
+\end{code}

semester3/numcs/parts/03_zeros/05_sectant-method.tex

@@ -6,3 +6,20 @@ Dann ist ein Schritt:
 \begin{align*}
     \tilde{F}(x) = F(x^{(k)}) + q^{(k)} (x - x^{(k)}) \Longrightarrow x^{(k + 1)} := x^{(k)} - \frac{F(x^{(k)})}{q^{(k)}}, \smallhspace \text{ falls } q^{(k)} \neq 0
 \end{align*}
+
+\innumpy Das Sekanten-Verfahren lässt sich im Prinzip ähnlich implementieren wie Newton:
+
+\begin{code}{python}
+def secant_method(f, x0: float, x1: float, tol=1e-12, maxIter=100):
+    """ Use secant method, which approximates the derivative """
+    f0 = f(x0)
+    for i in range(maxIter):
+        fn = f(x1)
+        secant = fn * (x1-x0) / (fn - f0)   # Approximate derivative via secant
+        x0 = x1; x1 -= secant
+        
+        if np.abs(secant) < tol: return x1, i
+        else: f0 = fn
+        
+    return None, maxIter
+\end{code}