STATUS: Draft

---

import numpy as np
import sympy as sp
import pickle
from IPython.display import HTML
import ipywidgets as widgets
import matplotlib as mpl
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt
mpl.rcParams['legend.fontsize'] = 10
import pandas as pd
import itertools

# function to print latex
def renderListToLatex(e):
    latex_rendering = []

    for i in range(len(e)):
        latex_rendering.append("$$" + sp.latex(e[i]) + "$$<br/>")
    
    return(HTML("".join(latex_rendering[0:])))

Solving Polynomial Equations (8)

Create needed variables

c_0, c_1, c_2, c_3, c_4, x, t, a_0, a_1, a_2, a_3, a_4, a_5, a_6, a_7, a_8, a_9, x_1, x_2, s_1, s_2 = sp.symbols('c_0, c_1, c_2, c_3, c_4, x, t, a_0, a_1, a_2, a_3, a_4, a_5, a_6, a_7, a_8, a_9, x_1, x_2, s_1, s_2')

Introduction

Recall from Solving Polynomials (7), the function $C$ that is equivalent using only $m_2$ and $m_3$.

$$C(m_2, m_3) \equiv(-1)^{m_3 + 1} \frac{(2 m_{2} + 3 m_{3})!}{(1 + m_{2} + 2 m_{3})!m_2!m_3!} \frac{c_0^{1 + m_{2} + 2 m_{3}} c_2^{m_2} c_3^{m_3} }{c_1^{2 m_{2} + 3 m_{3} + 1}}$$

Recall the formula to get values of this this function.

def C(m2, m3, returnCoefficientsOnly = False, returnCoefficientsOnlyWithoutSigns = False):
    c_0, c_1, c_2, c_3 = sp.symbols('c_0, c_1, c_2, c_3')
    s1 = (-1)**(m3 + 1)
    s2 = sp.factorial(2 * m2 + 3 * m3)
    s3 = sp.factorial(1 + m2 + 2 * m3) * sp.factorial(m2) * sp.factorial(m3)
    s4 = c_0**(1 + m2 + 2 * m3) * c_2**m2 *c_3**m3
    s5 = c_1**(2 * m2 + 3 * m3 + 1)
    
    if returnCoefficientsOnly:
        s6 = s1 * (s2 / s3)
    elif returnCoefficientsOnlyWithoutSigns:
        s6 = (s2 / s3)
    else:
        s6 = s1 * (s2 / s3) * (s4 / s5)
    return(s6)

Let g8 be a matrix of values.

g1 = np.arange(6)
g2 = [[C(j, i) for i in g1] for j in g1]
g3 = sp.Matrix(g2)
g3

$\displaystyle \left[\begin{matrix}- \frac{c_{0}}{c_{1}} & \frac{c_{0}^{3} c_{3}}{c_{1}^{4}} & - \frac{3 c_{0}^{5} c_{3}^{2}}{c_{1}^{7}} & \frac{12 c_{0}^{7} c_{3}^{3}}{c_{1}^{10}} & - \frac{55 c_{0}^{9} c_{3}^{4}}{c_{1}^{13}} & \frac{273 c_{0}^{11} c_{3}^{5}}{c_{1}^{16}}\\- \frac{c_{0}^{2} c_{2}}{c_{1}^{3}} & \frac{5 c_{0}^{4} c_{2} c_{3}}{c_{1}^{6}} & - \frac{28 c_{0}^{6} c_{2} c_{3}^{2}}{c_{1}^{9}} & \frac{165 c_{0}^{8} c_{2} c_{3}^{3}}{c_{1}^{12}} & - \frac{1001 c_{0}^{10} c_{2} c_{3}^{4}}{c_{1}^{15}} & \frac{6188 c_{0}^{12} c_{2} c_{3}^{5}}{c_{1}^{18}}\\- \frac{2 c_{0}^{3} c_{2}^{2}}{c_{1}^{5}} & \frac{21 c_{0}^{5} c_{2}^{2} c_{3}}{c_{1}^{8}} & - \frac{180 c_{0}^{7} c_{2}^{2} c_{3}^{2}}{c_{1}^{11}} & \frac{1430 c_{0}^{9} c_{2}^{2} c_{3}^{3}}{c_{1}^{14}} & - \frac{10920 c_{0}^{11} c_{2}^{2} c_{3}^{4}}{c_{1}^{17}} & \frac{81396 c_{0}^{13} c_{2}^{2} c_{3}^{5}}{c_{1}^{20}}\\- \frac{5 c_{0}^{4} c_{2}^{3}}{c_{1}^{7}} & \frac{84 c_{0}^{6} c_{2}^{3} c_{3}}{c_{1}^{10}} & - \frac{990 c_{0}^{8} c_{2}^{3} c_{3}^{2}}{c_{1}^{13}} & \frac{10010 c_{0}^{10} c_{2}^{3} c_{3}^{3}}{c_{1}^{16}} & - \frac{92820 c_{0}^{12} c_{2}^{3} c_{3}^{4}}{c_{1}^{19}} & \frac{813960 c_{0}^{14} c_{2}^{3} c_{3}^{5}}{c_{1}^{22}}\\- \frac{14 c_{0}^{5} c_{2}^{4}}{c_{1}^{9}} & \frac{330 c_{0}^{7} c_{2}^{4} c_{3}}{c_{1}^{12}} & - \frac{5005 c_{0}^{9} c_{2}^{4} c_{3}^{2}}{c_{1}^{15}} & \frac{61880 c_{0}^{11} c_{2}^{4} c_{3}^{3}}{c_{1}^{18}} & - \frac{678300 c_{0}^{13} c_{2}^{4} c_{3}^{4}}{c_{1}^{21}} & \frac{6864396 c_{0}^{15} c_{2}^{4} c_{3}^{5}}{c_{1}^{24}}\\- \frac{42 c_{0}^{6} c_{2}^{5}}{c_{1}^{11}} & \frac{1287 c_{0}^{8} c_{2}^{5} c_{3}}{c_{1}^{14}} & - \frac{24024 c_{0}^{10} c_{2}^{5} c_{3}^{2}}{c_{1}^{17}} & \frac{352716 c_{0}^{12} c_{2}^{5} c_{3}^{3}}{c_{1}^{20}} & - \frac{4476780 c_{0}^{14} c_{2}^{5} c_{3}^{4}}{c_{1}^{23}} & \frac{51482970 c_{0}^{16} c_{2}^{5} c_{3}^{5}}{c_{1}^{26}}\end{matrix}\right]$

Recall also that this value can be recovered as a sum and it is also one solution for a general cubic equation. Does it actually work in practice?

Let g4 be the equation where x is set to the sum of these values.

g4 = sp.Eq(x, np.sum(g3))
g4

$\displaystyle x = \frac{51482970 c_{0}^{16} c_{2}^{5} c_{3}^{5}}{c_{1}^{26}} + \frac{6864396 c_{0}^{15} c_{2}^{4} c_{3}^{5}}{c_{1}^{24}} + \frac{813960 c_{0}^{14} c_{2}^{3} c_{3}^{5}}{c_{1}^{22}} - \frac{4476780 c_{0}^{14} c_{2}^{5} c_{3}^{4}}{c_{1}^{23}} + \frac{81396 c_{0}^{13} c_{2}^{2} c_{3}^{5}}{c_{1}^{20}} - \frac{678300 c_{0}^{13} c_{2}^{4} c_{3}^{4}}{c_{1}^{21}} + \frac{6188 c_{0}^{12} c_{2} c_{3}^{5}}{c_{1}^{18}} - \frac{92820 c_{0}^{12} c_{2}^{3} c_{3}^{4}}{c_{1}^{19}} + \frac{352716 c_{0}^{12} c_{2}^{5} c_{3}^{3}}{c_{1}^{20}} + \frac{273 c_{0}^{11} c_{3}^{5}}{c_{1}^{16}} - \frac{10920 c_{0}^{11} c_{2}^{2} c_{3}^{4}}{c_{1}^{17}} + \frac{61880 c_{0}^{11} c_{2}^{4} c_{3}^{3}}{c_{1}^{18}} - \frac{1001 c_{0}^{10} c_{2} c_{3}^{4}}{c_{1}^{15}} + \frac{10010 c_{0}^{10} c_{2}^{3} c_{3}^{3}}{c_{1}^{16}} - \frac{24024 c_{0}^{10} c_{2}^{5} c_{3}^{2}}{c_{1}^{17}} - \frac{55 c_{0}^{9} c_{3}^{4}}{c_{1}^{13}} + \frac{1430 c_{0}^{9} c_{2}^{2} c_{3}^{3}}{c_{1}^{14}} - \frac{5005 c_{0}^{9} c_{2}^{4} c_{3}^{2}}{c_{1}^{15}} + \frac{165 c_{0}^{8} c_{2} c_{3}^{3}}{c_{1}^{12}} - \frac{990 c_{0}^{8} c_{2}^{3} c_{3}^{2}}{c_{1}^{13}} + \frac{1287 c_{0}^{8} c_{2}^{5} c_{3}}{c_{1}^{14}} + \frac{12 c_{0}^{7} c_{3}^{3}}{c_{1}^{10}} - \frac{180 c_{0}^{7} c_{2}^{2} c_{3}^{2}}{c_{1}^{11}} + \frac{330 c_{0}^{7} c_{2}^{4} c_{3}}{c_{1}^{12}} - \frac{28 c_{0}^{6} c_{2} c_{3}^{2}}{c_{1}^{9}} + \frac{84 c_{0}^{6} c_{2}^{3} c_{3}}{c_{1}^{10}} - \frac{42 c_{0}^{6} c_{2}^{5}}{c_{1}^{11}} - \frac{3 c_{0}^{5} c_{3}^{2}}{c_{1}^{7}} + \frac{21 c_{0}^{5} c_{2}^{2} c_{3}}{c_{1}^{8}} - \frac{14 c_{0}^{5} c_{2}^{4}}{c_{1}^{9}} + \frac{5 c_{0}^{4} c_{2} c_{3}}{c_{1}^{6}} - \frac{5 c_{0}^{4} c_{2}^{3}}{c_{1}^{7}} + \frac{c_{0}^{3} c_{3}}{c_{1}^{4}} - \frac{2 c_{0}^{3} c_{2}^{2}}{c_{1}^{5}} - \frac{c_{0}^{2} c_{2}}{c_{1}^{3}} - \frac{c_{0}}{c_{1}}$

Evaluate a cubic equation using a general

z = sp.symbols('z')
g5 = sp.Eq(1 + 5. * z - z**2 + 3 * z**3, 0)
g5

$\displaystyle 3 z^{3} - z^{2} + 5.0 z + 1 = 0$

Let g6 be the solutions for z in g5.

g6 = sp.solve(g5, z)
g6

[-0.188828952032499,
 0.261081142682916 - 1.3027289283418*I,
 0.261081142682916 + 1.3027289283418*I]

Let g7 be the solutions using the derived formula, $C$.

g7 = {c_0:1, c_1:5, c_2: -1, c_3:3}
g8 = sp.N(g4.rhs.subs(g7))
g8

$\displaystyle -0.188828892735343$

Note the answers are similiar.

Check another example

g9 = sp.Eq(1 + 15. * z - z**2 + 3 * z**3, 0)
g9

$\displaystyle 3 z^{3} - z^{2} + 15.0 z + 1 = 0$

g10 = sp.solve(g9, z)
g10

[-0.0663151597975300,
 0.199824246565432 - 2.23306359608708*I,
 0.199824246565432 + 2.23306359608708*I]

Let g12 be an answer obtained from the formula $C$

g11 = {c_0:1, c_1:15, c_2: -1, c_3:3}
g12 = sp.N(g4.rhs.subs(g11))
g12

$\displaystyle -0.0663151597974653$

Note they are again similiar

# demonstrating a viable technology in some cases - we can solve a cubic 
# equation with high school algebra - you don't need cubed roots / square roots 
# note we don't get final solution

# variant

# try another examle
# note we are looking at non complex solutions

g11 = sp.Eq(1 + 3. * z - 4 * z**2 + 5 * z**3, 0)
g11
g12 = sp.solve(g11, z)
g12

[-0.236610107250854,
 0.518305053625427 - 0.75936308841047*I,
 0.518305053625427 + 0.75936308841047*I]

g13 = {c_0:1, c_1:3, c_2: -4, c_3:5}
g14 = sp.N(g4.rhs.subs(g13))
g14

$\displaystyle -33.9803439976563$

# this is clearly incorrect 

# why is this the case - consider the terms - denomitor has big power
# but it is a smaller number to large power makes it difficult - nature 
# of certain choice - c_1 needs to be bigger than other 1

# try another varint

g15 = sp.Eq(2 + 3. * z - 7 * z**2 + z**3, 0)
g15
g16 = sp.solve(g15, z)
g16

[-0.355967716679387 + 0.e-20*I,
 0.865675512845961 - 0.e-22*I,
 6.49029220383343 - 0.e-21*I]

g17 = {c_0:2, c_1:3, c_2: -7, c_3:1}
g18 = sp.N(g4.rhs.subs(g17))
g18

$\displaystyle -12825.756404674$

# so very incorrect but....reconsider our method words well when coefficnt 
# of z is bigger than others
# so perhaps we can transform uor equation to make this the case 
# one more likely to produce solution

# let z = 1 / u
g15

$\displaystyle z^{3} - 7 z^{2} + 3.0 z + 2 = 0$

# introduce new version  whre z = 1 / u. Then solve for u
u = sp.symbols('u')
g19 = g15.subs(z, 1/u)
g20 = sp.Eq(sp.simplify(g19.lhs * u**3), 0)
g20

$\displaystyle 2 u^{3} + 3.0 u^{2} - 7 u + 1 = 0$

g21 = sp.solve(g20, u)
g21

[-2.80924351603682 + 0.e-23*I,
 0.154076267846517 + 0.e-22*I,
 1.15516724819031 - 0.e-23*I]

# and now wolve in our way

g22 = {c_0:1, c_1:-7, c_2: 3, c_3:2}
g23 = sp.N(g4.rhs.subs(g22))
g23

$\displaystyle 0.154074694041552$

# so now similair. so 

# now recover z

g24 = 1 / g23
g24

$\displaystyle 6.49035849930239$

# this is a solution to original one. So yes our current solution doens't cover all 
# solutions, but if we are flexible, we can move into position where formula
# might be ale to use this 

# we don't yet know conditions as yet, but this technology does 
# solve cubic equations 

# need to turn back to coefficients - we know it work, time to move on