我想制作一个 Python 程序,该程序将运行二分法来确定以下内容的根:
f(x) = -26 + 85x - 91x2 +44x3 -8x4 + x5
二分法是一种用于估计多项式 f(x) 的根的数值方法。
是否有任何可用的伪代码、算法或库可以用来告诉我答案?
我想制作一个 Python 程序,该程序将运行二分法来确定以下内容的根:
f(x) = -26 + 85x - 91x2 +44x3 -8x4 + x5
二分法是一种用于估计多项式 f(x) 的根的数值方法。
是否有任何可用的伪代码、算法或库可以用来告诉我答案?
这是一些显示基本技术的代码:
>>> def samesign(a, b):
return a * b > 0
>>> def bisect(func, low, high):
'Find root of continuous function where f(low) and f(high) have opposite signs'
assert not samesign(func(low), func(high))
for i in range(54):
midpoint = (low + high) / 2.0
if samesign(func(low), func(midpoint)):
low = midpoint
else:
high = midpoint
return midpoint
>>> def f(x):
return -26 + 85*x - 91*x**2 +44*x**3 -8*x**4 + x**5
>>> x = bisect(f, 0, 1)
>>> print(x, f(x))
0.557025516287 3.74700270811e-16
要在达到给定容差时提前退出,请在循环末尾添加一个测试:
def bisect(func, low, high, tolerance=None):
assert not samesign(func(low), func(high))
for i in range(54):
midpoint = (low + high) / 2.0
if samesign(func(low), func(midpoint)):
low = midpoint
else:
high = midpoint
if tolerance is not None and abs(high - low) < tolerance:
break
return midpoint
您可以在此处使用scipy.optimize.bisect的早期 Stack Overflow 问题中看到解决方案。或者,如果您的目的是学习,则Wikipedia 条目中有关二分法的伪代码是在 Python 中进行自己的实现的一个很好的指南,正如前面问题的评论者所建议的那样。
我的实现比其他解决方案更通用但更简单:(和公共领域)
def solve(func, x = 0.0, step = 1e3, prec = 1e-10):
"""Find a root of func(x) using the bisection method.
The function may be rising or falling, or a boolean expression, as long as
the end points have differing signs or boolean values.
Examples:
solve(lambda x: x**3 > 1000) to calculate the cubic root of 1000.
solve(math.sin, x=6, step=1) to solve sin(x)=0 with x=[6,7).
"""
test = lambda x: func(x) > 0 # Convert into bool function
begin, end = test(x), test(x + step)
assert begin is not end # func(x) and func(x+step) must be on opposite sides
while abs(step) > prec:
step *= 0.5
if test(x + step) is not end: x += step
return x
# Defining Function
def f(x):
return x**3-5*x-9
# Implementing Bisection Method
def bisection(x0,x1,e):
step = 1
print('\n\n*** BISECTION METHOD IMPLEMENTATION ***')
condition = True
while condition:
x2 = (x0 + x1)/2
print('Iteration-%d, x2 = %0.6f and f(x2) = %0.6f' % (step, x2, f(x2)))
if f(x0) * f(x2) < 0:
x1 = x2
else:
x0 = x2
step = step + 1
condition = abs(f(x2)) > e
print('\nRequired Root is : %0.8f' % x2)
# Input Section
x0 = input('First Guess: ')
x1 = input('Second Guess: ')
e = input('Tolerable Error: ')
# Converting input to float
x0 = float(x0)
x1 = float(x1)
e = float(e)
#Note: You can combine above two section like this
# x0 = float(input('First Guess: '))
# x1 = float(input('Second Guess: '))
# e = float(input('Tolerable Error: '))
# Checking Correctness of initial guess values and bisecting
if f(x0) * f(x1) > 0.0:
print('Given guess values do not bracket the root.')
print('Try Again with different guess values.')
else:
bisection(x0,x1,e)
此外, codesansar.com/numerical-methods/拥有大量算法、伪代码和程序,它们使用不同的编程语言进行数值分析。
有公差:
# there is only one root
def fn(x):
return x**3 + 5*x - 9
# define bisection method
def bisection( eq, segment, app = 0.3 ):
a, b = segment['a'], segment['b']
Fa, Fb = eq(a), eq(b)
if Fa * Fb > 0:
raise Exception('No change of sign - bisection not possible')
while( b - a > app ):
x = ( a + b ) / 2.0
f = eq(x)
if f * Fa > 0: a = x
else: b = x
return x
#test it
print bisection(fn,{'a':0,'b':5}, 0.00003) # => 1.32974624634
我正在尝试进行一些修改:
这是我的代码:
import numpy as np
def bisection3(f,x0,x1,tol,max_iter):
c = (x0+x1)/2.0
x0 = c
Xs = 0.3604217029603
err_Abs = np.linalg.norm(x0-Xs)
itr = 0
f_x0 = f(x0)
f_x1 = f(x1)
xp = [] # sucesion que converge a la raiz
errores_Abs = []
errores_Aprox = []
fs = [] # esta sucecion debe converger a cero
while(((x1-x0)/2.0 > tol) and (itr< max_iter)):
if f(c) == 0:
return c
elif f(x0)*f(c) < 0:
x1 = c
else :
x0 = c
c = (x0+x1)/2.0
itr +=1
err_Abs = np.linalg.norm(x0-Xs)
err_Aprox = np.linalg.norm(x1 - x0)
fs.append(f(c))
xp.append(c)
errores_Abs.append(err_Abs)
errores_Aprox.append(err_Aprox)
return x0,errores_Abs, errores_Aprox,fs,xp
我有一个执行示例:
f = lambda x : 3*x + np.sin(x) - np.exp(x)
X0_r1 , err_Abs_r1,err_Aprox_r1, fs_r1 , xp_r1 = bisection3(f,0,0.5,1e-5,100)
可以修改上面的二分法,用一个公差作为止动器:
def samesign(a, b):
return a*b > 0
def bisect(func, low, high):
assert not samesign(func(low), func(high))
n = 20
e = 0.01 # the epsilon or error we justify for tolerance
for i in range(n):
if abs(func(low)-func(high)) > e:
midpoint = (low + high) / 2
print(i, midpoint, f(midpoint))
if samesign(func(low), func(midpoint)):
low = midpoint
else:
high = midpoint
else:
return round(midpoint, 2)
return round(midpoint, 2)
def f(x):
return -26 + 58 * x - 91 * x**2 + 44 * x**3 - 8 * x**4 + x**5
def relative_error(p, p_old):
return abs((p - p_old) / p)
def bisection(a, b, no_iterations, accuracy):
i = 1
f_a = f(a)
p_old = 0
while i <= no_iterations:
p = (a + b) / 2
f_p = f(p)
error = relative_error(p, p_old)
print('%2d %4f %4f %6f %6f %4.6f' % (i, a, b, p, f_p, error))
if error < accuracy or f_p == 0:
break
p_old = p
i += 1
if f_a * f_p > 0:
a = p
f_a = f_p
else:
b = p
bisection(1, 2, 20, 0.00001)