Решение кубического уравнения в python без numpy и solve

Кто-нибудь знает как сделать решение кубического уравнения в python без numpy и solve ?


Ответы (2 шт):

Автор решения: Сергей

Вот. Честно найдено на: https://www.cyberforum.ru/python-tasks/thread2874633.html за 30 секунд через поисковый запрос в Google "решение кубического уравнения python" (ответ был первым в выдаче). Теперь у нас в базе знаний есть решение кубического уравнения без numpy и solve...

def sign(x):  # вроде бы есть в numpy, но и написать несложно...
    if x>0:
        return 1
    elif x<0:
        return -1
    else:
        return 0
        
def dicho(f,a,b,eps=1.0e-14):  # Дихотомия. Один вещественный корень у  
    fa=f(a)                    # кубического полиномаесть всегда. Его и найдем...
    fb=f(b)
    while True:
        c=0.5*(a+b)
        if abs(b-a)<eps:
            return c
        fc=f(c)
        if abs(fc)<=eps:
            return c
        if sign(fa)*sign(fc)<0:
            b=c
            fb=fc
        else:
            a=c
            fa=fc
            
def div_poly(p,a):  # делим исходный полином на (x-a) и получаем квадратный
    r=[0,0,0]          # трехчлен, поставляющий два недостающих корня
    r[2]=p[3]
    r[1]=p[2]+a*p[3]
    r[0]=(p[1]+a*(p[2]+a*p[3]))
    return tuple(r)
    
def solve_qube(p):  # парадная функция
    q=max(p)
    left=-abs(q)/abs(p[3])  # границы вещественного корня
    right=-left
    x1=dicho(lambda x: p[3]*x**3+p[2]*x**2+p[1]*x+p[0],left,right)
    (c,b,a)=div_poly(p,x1)
    d=b**2-4*a*c
    x2=(-b+d**0.5)/(2*a)
    x3=(-b-d**0.5)/(2*a)
    return (x1,x2,x3)
    
print(solve_qube([12,4,8,5]))

Вывод:

(-1.8623879412594293, (0.13119397062971466+1.1275887015236796j), (0.13119397062971452-1.1275887015236796j))

→ Ссылка
Автор решения: MIkhail

вот решение кубического уравнения в общем виде

from cmath import *

solutions = set()


def cbrt(polynomial):
    solution = set()
    root1 = polynomial ** (1 / 3)
    root2 = (polynomial ** (1 / 3)) * (-1 / 2 + (sqrt(3) * 1j) / 2)
    root3 = (polynomial ** (1 / 3)) * (-1 / 2 - (sqrt(3) * 1j) / 2)
    solution.update({root1, root2, root3})
    return solution


def cardano(a, b, c, d):
    p = (3 * a * c - b ** 2) / (3 * a ** 2)
    q = (2 * b ** 3 - 9 * a * b * c + 27 * a ** 2 * d) / (27 * a ** 3)
    alpha = cbrt(-q / 2 + sqrt((q / 2) ** 2 + (p / 3) ** 3))
    beta = cbrt(-q / 2 - sqrt((q / 2) ** 2 + (p / 3) ** 3))
    for i in alpha:
        for j in beta:
            if abs((i * j) + p / 3) <= 0.00001:
                x = i + j - b / (3 * a)
                solutions.add(x)


def quadratic(a, b, c):
    D = b ** 2 - 4 * a * c
    x1 = (-b + sqrt(D)) / 2 * a
    x2 = (-b - sqrt(D)) / 2 * a
    solutions.update({x1, x2})


def linear(a, b):
    if a == 0 and b == 0:
        solutions.add("True")

    if a == 0 and b != 0:
        solutions.add("False")

    if a != 0:
        solutions.add(-b / a)


print('ax^3+bx^2+cx+d=0')
a, b, c, d = map(complex, input().split())

if a != 0:
    cardano(a, b, c, d)

elif b != 0:
    quadratic(b, c, d)

else:
    linear(c, d)

print(solutions)
→ Ссылка