50 lines
1.7 KiB
Python
50 lines
1.7 KiB
Python
import math
|
||
|
||
def NewtonDownHillSolve(fx, dfx, x0, err1,err2, N0,min_t):
|
||
count = 0
|
||
print(f"k={count}, x0={x0}\n")
|
||
x1 = x0 + 1 + err1
|
||
while abs(x1 - x0) > err1 or abs(fx(x1)) > err2:
|
||
t = 1
|
||
if abs(dfx(x1)) < 1e-10:
|
||
print("导数为0，无法下山")
|
||
return None, 0
|
||
print(f"当前点: x0={x0}")
|
||
while t >= min_t:
|
||
x1 = x0 - t * fx(x0) / dfx(x0)
|
||
print(f"下山: t={t}, x1={x1}, fx(x{count+1})={fx(x1)}, fx(x{count})={fx(x0)}")
|
||
if abs(fx(x1)) < abs(fx(x0)):
|
||
break
|
||
t *= 0.5
|
||
if t < min_t:
|
||
print("达到最小t，下山失败")
|
||
return None, -2
|
||
# x1 = x0 - fx(x0) / dfx(x0)
|
||
count += 1
|
||
print(f"k={count}, x{count}={x1},x1-x0={abs(x1-x0)}\n")
|
||
if count > N0:
|
||
return None, -1
|
||
x0 = x1
|
||
print(f"收敛: x1={x1}, fx(x1)={fx(x1)}")
|
||
return x1, 1
|
||
|
||
|
||
if __name__ == "__main__":
|
||
##############################################################################################################
|
||
err1 = 1e-5 # 根的误差限 见P147
|
||
err2 = 1e-5 # 残量精度 见P147
|
||
N0 = 100 # 最大迭代次数
|
||
min_t = 1e-10 # 最小t值
|
||
x0 = 0.6 # 初始值
|
||
fx = lambda x: x**3 - x - 1 # 原函数
|
||
dfx = lambda x: 3*x**2 - 1 # 导函数
|
||
|
||
result, status = NewtonDownHillSolve(fx, dfx, x0, err1, err2, N0, min_t)
|
||
if status == 1:
|
||
print(f"收敛 解为: {result}")
|
||
elif status == -1:
|
||
print("不收敛")
|
||
elif status == -2:
|
||
print("下山失败")
|
||
else:
|
||
print("导数为0，无法收敛") |