64 lines
1.6 KiB
Python
64 lines
1.6 KiB
Python
#抛物线法解方程（主）
|
||
def MullerSolve(fx,x0,x1,x2,err1,err2,N):
|
||
count = 0
|
||
f0 = fx(x0)
|
||
f1 = fx(x1)
|
||
f2 = fx(x2)
|
||
q = (x2 - x1) / (x1 - x0)
|
||
p = 0
|
||
a = 0
|
||
b = 0
|
||
c = 0
|
||
while True:
|
||
p = (x2 - x0) / (x1 - x0)
|
||
|
||
a = q**2 * f0 - q*p*f1 + q*f2
|
||
b = q**2 *f0 - p**2 *f1 + (p + q)*f2
|
||
c = p*f2
|
||
h1 = 0
|
||
if b < 0:
|
||
h1 = -2 * c / (b - (b**2 - 4*a*c)**0.5)
|
||
else:
|
||
h1 = -2 * c / (b + (b**2 - 4*a*c)**0.5)
|
||
x3 = x2 + h1 * (x2 - x1)
|
||
f3 = fx(x3)
|
||
|
||
print(f"k={count}: x{count}={x0:.5}, f(x{count})={f0:.5}; x{count+1}={x1:.5}, f(x{count+1})={f1:.5}; x{count+2}={x2:.5}, f(x{count+2})={f2:.5}; x{count+3}={x3:.5}, f(x{count+3})={f3:.5}")
|
||
print(f"p={p:.5}, q={q:.5}, a={a:.5}, b={b:.5}, c={c:.5}, h={h1:.5}")
|
||
|
||
k = err1 + 1
|
||
if abs(f3) < 1:
|
||
k = abs(x3 - x2)
|
||
else:
|
||
k = abs(x3 - x2) / abs(f3)
|
||
|
||
if abs(f3) < err2 or k < err1:
|
||
return x3, 1
|
||
count += 1
|
||
if count > N:
|
||
return None, 0
|
||
x0 = x1
|
||
x1 = x2
|
||
x2 = x3
|
||
f0 = f1
|
||
f1 = f2
|
||
f2 = f3
|
||
q = h1
|
||
|
||
#精度要求#########################（152页）
|
||
if __name__ == "__main__":
|
||
err1 = 1e-5
|
||
err2 = 1e-5
|
||
N = 100
|
||
##把初始函数和初始值改成题干要求的##########################
|
||
x0 = 0.3
|
||
x1 = 0.5
|
||
x2 = 0.4
|
||
fx = lambda x: 8*x**4 - 8*x**2 + 1
|
||
|
||
result, status = MullerSolve(fx, x0, x1, x2, err1, err2, N)
|
||
if status == 1:
|
||
print(f"fx收敛 解为: {result}")
|
||
else:
|
||
print("fx不收敛")
|