scientific-programming-exer.../ex_46.py

63 lines
1.0 KiB
Python
Raw Normal View History

2019-02-25 12:07:52 +00:00
import numpy as np
2019-01-23 14:27:28 +00:00
def bisec(f, a, b, eps, nmax):
"""
Takes f:[a,b] -> |R and tries to compute the null point
using bisection.
"""
if(f(a)*f(b) >= 0):
raise ValueError("f(a)*f(b) >= 0")
x_minus = a
x_plus = b
if(a > b):
x_minus = b
x_plus = a
for i in range(nmax):
x_minus, x_plus, err = bisec_one(x_minus, x_plus, f, eps)
if(err < eps):
break
if(err < eps):
return x_minus, err
raise ValueError("bisection hit nmax")
def bisec_one(x_minus, x_plus, f, eps):
a = x_minus
b = x_plus
x_m = (b + a) / 2
y_m = f(x_m)
if(y_m < 0):
if(y_m > -1*eps):
return x_m, x_m, -y_m
if(f(x_m)*f(b) >= 0):
return a, x_m, -y_m
return x_m, b, -y_m
if(y_m < eps):
return x_m, x_m, y_m
if(f(a)*f(x_m) >= 0):
return x_m, b, y_m
return a, x_m, y_m
if( __name__ == "__main__"):
f1 = lambda x: x
f2 = lambda x: x**3
f3 = lambda x: -x + 1
f4 = lambda x: -x**3 + 1
2019-02-25 12:07:52 +00:00
f4 = lambda x: (x - 2)*np.exp(-x**2)
2019-01-23 14:27:28 +00:00
fs = [f1, f2, f3, f4]
for f in fs:
2019-02-25 12:07:52 +00:00
print(bisec(f, -12, 10, 0.0000001, 100))
print(bisec(f4, 1.2, 2.4, 0.001, 100))