Consider:
Find root of $f(x)=(x-2)^2 - 1$.
import numpy as np
def bisection(f, a, b, tol=1E-5, maxit=10000):
fa = f(a)
fb = f(b)
if fa==0: return a, 0
if fb==0: return b, 0
if (fa>=0 and fb>=0) or (fa<=0 and fb<=0) :
print(f"Error, a and b don't bracket the root")
return np.nan, -1
for nit in range(1,maxit+1):
c = 0.5*(a+b)
if np.abs(a-b) <= tol :
return c, nit
fc = f(c)
if fc == 0: return c, nit
if (fc > 0 and fa < 0) or (fc < 0 and fa > 0):
b = c
fb = fc
else :
a = c
fa = fc
print(f"Warning, no convergence in {maxit} iterations")
return c, maxit
#--------------------
def func1(x):
return (x-2)**2 - 1
#--------------------
x, nit = bisection(func1, 2., 3.5, tol=1E-8, maxit=100)
print(f"Solution x = {x} found in {nit} iterations")
print(f"f(x) = {func1(x)}")
Solution x = 2.9999999990686774 found in 29 iterations f(x) = -1.862645149230957e-09
1/2**27 / 2, (3-2.9999999962747097)
fsolve
¶from scipy.optimize import fsolve
xguess = 2
x = fsolve(func1, xguess)[0]
print(f"Solution x = {x}")
Solution x = 3.0