source: anuga_work/development/sudi/sw_1d/periodic_waves/cg/bisect.py @ 7837

Last change on this file since 7837 was 7837, checked in by mungkasi, 12 years ago

Again, adding some codes for 1d problems on debris avalanche and periodic waves.

File size: 902 bytes
Line 
1## module bisect
2''' root = bisect(f,x1,x2,switch=0,tol=1.0e-9).
3    Finds a root of f(x) = 0 by bisection.
4    The root must be bracketed in (x1,x2).
5    Setting switch = 1 returns root = None if
6    f(x) increases as a result of a bisection.
7'''   
8from math import log,ceil
9import error
10
11def bisect(f,x1,x2,switch=0,epsilon=1.0e-15):#9):
12    f1 = f(x1)
13    if f1 == 0.0: return x1
14    f2 = f(x2)
15    if f2 == 0.0: return x2
16    if f1*f2 > 0.0: error.err('Root is not bracketed')
17    n = ceil(log(abs(x2 - x1)/epsilon)/log(2.0))
18    for i in range(n):
19        x3 = 0.5*(x1 + x2); f3 = f(x3)
20        if (switch == 1) and (abs(f3) >abs(f1)) \
21                         and (abs(f3) > abs(f2)):
22            return None   
23        if f3 == 0.0: return x3
24        if f2*f3 < 0.0:
25            x1 = x3; f1 = f3
26        else:
27            x2 = x3; f2 = f3
28    return (x1 + x2)/2.0
Note: See TracBrowser for help on using the repository browser.