一元三次、四次方程求解

一、一元四次方程求解

一元四次方程求根公式，百科： https://baike.baidu.com/item/%E4%B8%80%E5%85%83%E5%9B%9B%E6%AC%A1%E6%96%B9%E7%A8%8B%E6%B1%82%E6%A0%B9%E5%85%AC%E5%BC%8F/10721996?fr=aladdin

1.网上搜到的求解代码

python求解代码，见：https://github.com/Larissa1990/Solve-cubic-and-quartic-equations-with-one-unknown/blob/master/Equations.py 经测试，有bug，求解不正确

说明，见：https://www.cnblogs.com/larissa-0464/p/11706131.html

2. 一元四次方程在线求解工具：https://www.osgeo.cn/app/s2085 可用来验证求解程序是否正确；

3. 一元四次方程，沈天珩简化求根公式, 具体公式详见百科链接；

import math
import cmath
import numpy as npdef cal_quartic_ik(args_list):a, b, c, d, e = args_listD = 3*pow(b,2) - 8*a*cE = -pow(b, 3) + 4*a*b*c - 8*pow(a, 2)*dF = 3*pow(b, 4) + 16*pow(a, 2)*pow(c, 2) - 16*a*pow(b, 2)*c + 16*pow(a, 2)*b*d - 64*pow(a, 3)*eA = D**2 - 3*FB = D*F - 9*pow(E, 2)C = F**2 - 3*D*pow(E, 2)delta = B**2 - 4*A*C  # 总判别式if (D == 0) & (E == 0) & (F == 0):""" 四重实根"""x = -b/(4*a)return 1, [x]if (A == 0) & (B == 0) & (C == 0) & (D*E*F != 0):""" 两个实根，其中一个三重实根"""x1 = (-b*D + 9*E)/(4*a*D)x234 = (-b * D - 3 * E) / (4 * a * D)return 2, [x1, x234]if (E == 0) & (F == 0) & (D != 0):""" 一对二重根"""if D>0:  # 根为实数x13 = (-b + math.sqrt(D))/(4*a)x24 = (-b - math.sqrt(D)) / (4 * a)return 2, [x13, x24]if D<0:  # 根为虚数# x13 = (-b + cmath.sqrt(D))/(4*a)# x24 = (-b - cmath.sqrt(D)) / (4 * a)return 0, 0if (A*B*C != 0) & (delta == 0):""" 一对二重实根 """x3 = (-b - np.sign(A*B*E)*math.sqrt( D - B/A))/(4*a)x4 = (-b - np.sign(A*B*E)*math.sqrt( D - B/A))/(4*a)if A*B>0 :  # 其余两根为不等实根x1 = (-b + np.sign(A*B*E)*math.sqrt( D - B/A) + math.sqrt( 2*B/A) )/(4*a)x2 = (-b + np.sign(A * B * E) * math.sqrt(D - B / A) - math.sqrt(2 * B / A)) / (4 * a)return 4, [x1, x2, x3, x4]if A*B < 0:  # 其余两根为共轭虚根# x1 = (-b + np.sign(A * B * E) * math.sqrt(D - B / A) + cmath.sqrt(2 * B / A)) / (4 * a)# x2 = (-b + np.sign(A * B * E) * math.sqrt(D - B / A) - cmath.sqrt(2 * B / A)) / (4 * a)return 2,  [x3, x4]if delta > 0:"""" 两个不等实根和一对共轭虚根"""z1 = A*D + 3*(( -B + math.sqrt(delta))/2.0)z2 = A * D + 3 * ((-B - math.sqrt(delta)) / 2.0)# print """ z1 =  """, z1# print """ z2 =  """, z2# print """ abs(z1) =  """, abs(z1)# print """ abs(z2) =  """, abs(z2)z = D**2 - D*(np.sign(z1)*pow(abs(z1), 1.0/3.0) + np.sign(z2)*pow(abs(z2), 1.0/3.0)) + \(np.sign(z1)*pow(abs(z1), 1.0/3.0) + np.sign(z2)*pow(abs(z2), 1.0/3.0))**2 - 3*Ax1 = (-b + np.sign(E)*math.sqrt((D + np.sign(z1)*pow(abs(z1), 1.0/3.0) + np.sign(z2)*pow(abs(z2), 1.0/3.0))/3.0)+ math.sqrt((2*D - np.sign(z1)*pow(abs(z1), 1.0/3.0) - np.sign(z2)*pow(abs(z2), 1.0/3.0)+ 2*math.sqrt(z))/3.0))/(4*a)x2 = (-b + np.sign(E)*math.sqrt((D + np.sign(z1)*pow(abs(z1), 1.0/3.0) + np.sign(z2)*pow(abs(z2), 1.0/3.0))/3.0)- math.sqrt((2*D - np.sign(z1)*pow(abs(z1), 1.0/3.0) - np.sign(z2)*pow(abs(z2), 1.0/3.0)+ 2*math.sqrt(z))/3.0))/(4*a)# 虚根忽略return 2, [x1, x2]if delta < 0:if E == 0:if (D>0) & (F>0) :""" 四个不等实根 """x1 = (-b + math.sqrt(D + 2*math.sqrt(F))) / (4 * a)x2 = (-b - math.sqrt(D + 2 * math.sqrt(F))) / (4 * a)x3 = (-b + math.sqrt(D - 2 * math.sqrt(F))) / (4 * a)x4 = (-b - math.sqrt(D - 2 * math.sqrt(F))) / (4 * a)return 4, [x1, x2, x3, x4]else:""" 两对不等共轭虚根 """# 虚根忽略print " 两对不等共轭虚根 "return 0, 0else:if (D > 0) & (F > 0):""" 四个不等实根 """theta = math.acos((3*B-2*A*D)/(2*A*math.sqrt(A)))y1 = (D - 2*math.sqrt(A)*math.cos(theta/3.0))/3.0y2 = (D + math.sqrt(A)*(math.cos(theta/3.0) + math.sqrt(3)*math.sin(theta/3.0)))/3.0y3 = (D + math.sqrt(A) * (math.cos(theta / 3.0) - math.sqrt(3) * math.sin(theta / 3.0))) / 3.0x1 = (-b + np.sign(E) * math.sqrt(y1) + ( math.sqrt(y2) + math.sqrt(y3)))/(4*a)x2 = (-b + np.sign(E) * math.sqrt(y1) - (math.sqrt(y2) + math.sqrt(y3))) / (4 * a)x3 = (-b - np.sign(E) * math.sqrt(y1) + (math.sqrt(y2) - math.sqrt(y3))) / (4 * a)x4 = (-b - np.sign(E) * math.sqrt(y1) - (math.sqrt(y2) - math.sqrt(y3))) / (4 * a)return 4, [x1, x2, x3, x4]else:""" 两对不等共轭虚根 """# 虚根忽略print " 两对不等共轭虚根 "return 0, 0

经测试，求解正确；

程序中，忽略了复数解，返回的解，均是实数解；可作为一基础求解工具；

一、一元三次方程求解

盛金公式： https://baike.baidu.com/item/%E7%9B%9B%E9%87%91%E5%85%AC%E5%BC%8F#9

def cal_cubic_ik(args_list):a, b, c, d = args_listA = b**2 - 3*a*cB = b*c - 9*a*dC = c**2 - 3*b*ddelta = B**2 - 4*A*C  # 总判别式if (A == 0) & (B == 0):""" 一个三重实根"""x = -b/(3*a)return 1, [x]if delta > 0:"""" 一个实根和一对共轭复根"""y1 = A*b + 3*a*(( -B + math.sqrt(delta))/2.0)y2 = A*b + 3*a*(( -B - math.sqrt(delta))/2.0)print """ y1 =  """, y1print """ y2 =  """, y2print """ abs(y2) =  """, abs(y2)# 虚根忽略x1 = (-b - (np.sign(y1)*pow(abs(y1), 1.0 / 3.0) + np.sign(y2)*pow(abs(y2), 1.0 / 3.0)) )/ (3 * a)  # 负数直接开立方根会报错return 1, [x1]if delta == 0:"""三个实根，其中有一个二重根 """K = B/float(A)x1 = -b/float(a) + Kx23 = -K/2.0return 2,  [x1, x23]if delta < 0:""" 三个不等实根 """T = (2*A*b - 3*a*B)/(2*math.sqrt(pow(A, 3)))theta = math.acos(T)x1 = (-b - 2 * math.sqrt(A) * math.cos(theta / 3.0)) / (3 * a)x2 = (-b + math.sqrt(A)*(math.cos(theta/3.0)) + math.sqrt(3)*math.sin(theta/3.0))/(3*a)x3 = (-b + math.sqrt(A) * (math.cos(theta / 3.0)) - math.sqrt(3) * math.sin(theta / 3.0)) / (3 * a)return 3, [x1, x2, x3]

程序中，忽略了复数解，返回的解，均是实数解；