Complex cube root

Introduction

File:0419polarDiagram.png

Complex number W = complex number w³.
Origin at point $(0,0)$.
$w_{real},W_{real}$ parallel to $X$ axis.
$w_{imag},W_{imag}$ parallel to $Y$ axis.
$w_{real}=1.2;w_{imag}=0.5;w_{mod}=\sqrt{1.2^2+0.5^2}=1.3.$
$W_{real}=0.828;;W_{imag}=2.035;$
$W_{mod}=\sqrt{0.828^2+2.035^2}=2.197.$
$W_{mod}=w_{mod}^3=1.3^3=2.197.$
$W_{\varphi }=3w_{\varphi }.$ By cosine triple angle formula:
$\cos W_{phi}=4\cdot (\frac{1.2}{1.3}{)}^3-3\cdot \frac{1.2}{1.3}=\frac{828}{2197}=\frac{W_{real}}{W_{mod}}.$
See "3 cube roots of W" in Gallery below.

Let complex numbers $W=$ a $+$ b$\cdot i$ and $w=$ p $+$ q$\cdot i.$ Let $W=w^3.$

When given $a,b,$ aim of this page is to calculate $p,q.$

In the diagram complex number $w=p+qi=w_{real}+i\cdot w_{imag}=w_{mod}(\cos w_{\varphi }+i\cdot \sin w_{\varphi }),$ where

• $w_{real},w_{imag}$ are the real and imaginary components of $w,$ the rectangular components.

• $w_{mod},w_{\varphi }$ are the modulus and phase of $w,$ the polar components.

Similarly, $W_{real},W_{imag},W_{mod},W_{\varphi }$ are the corresponding components of $W.$

$W=w^3=(w_{mod}(\cos w_{\varphi }+i\cdot \sin w_{\varphi }){)}^3$

$=w_{mod}^3(\cos (3w_{\varphi })+i\cdot \sin (3w_{\varphi }))$

$=W_{mod}(\cos (W_{\varphi })+i\cdot \sin (W_{\varphi }))$

$=W$

There are 2 significant calculations:

$w_{mod}=\sqrt[3]{W_{mod}}$ and

$\cos w_{\varphi }=\cos \frac{W_{\varphi }}{3}.$ Template:RoundBoxBottom

Implementation

Cos φ from cos 3φ

File:1008Cos3A.png

Graph of $\cos (3A)=4{\cos }^3(A)-3\cos (A)$
Formula and/or graph are used to calculate $\cos (A)$ if $\cos (3A)$ is given.

The cosine triple angle formula is: $\cos (3\theta )=4{\cos }^3\theta -3\cos \theta .$ This formula, of form $y=4x^3-3x,$ permits $\cos (3\theta )$ to be calculated if $\cos \theta$ is known.

If $\cos (3\theta )$ is known and the value of $\cos \theta$ is desired, this identity becomes: $4{\cos }^3\theta -3\cos \theta -\cos (3\theta )=0.$ $\cos \theta$ is the solution of this cubic equation.

In fact this equation has three solutions, the other two being $\cos (\theta \pm 120^{\circ }).$

$\cos (3(\theta \pm 120^{\circ }))=\cos (3\theta \pm 360^{\circ })=\cos (3\theta ).$

It is sufficient to calculate only $\cos \theta$ from $\cos 3\theta ,$ accomplished by the following code:

# python code

cosAfrom_cos3Adebug = 0

def cosAfrom_cos3A(cos3A) :
    cos3A = Decimal(str(cos3A))+0
    if 1 >= cos3A >= -1 : pass
    else :
        print ('cosAfrom_cos3A(cos3A) : cos3A not in valid range.')
        return None
    '''
    if cos3A == 0 :
        A = 90 and cos3A = cosA
    if cos3A == 1 :
        A = 0 and cos3A = cosA
    if cos3A == -1 :
        A = 180 and cos3A = cosA
'''
    if cos3A in (0,1,-1) : return cos3A

    # From the cosine triple angle formula:
    a,b,c,d = 4,0,-3,-cos3A

    # prepare for Newton's method.
    if d < 0 : x = Decimal(1)
    else : x = -Decimal(1)
    count = 31; L1 = [x]; almostZero = Decimal('1e-' + str(prec-2))

    # Newton's method:
    while 1 :
        count -= 1
        if count <= 0 :
            print ('cosAfrom_cos3A(cos3A): count expired.')
            break
        y = a*x*x*x + c*x + d
        if cosAfrom_cos3Adebug :
            print ('cosAfrom_cos3A(cos3A) : x,y =',x,y)
        if abs(y) < almostZero : break
        slope = 12*x*x + c
        delta_x = y/slope
        x -= delta_x
        if x in L1[-1:-5:-1] :
            # This value of x has been used previously.
            print ('cosAfrom_cos3A(cos3A): endless loop detected.')
            break
        L1 += [x]

    return x

Template:RoundBoxTop When $x==1,$ slope $=12x^2-3=9.$ Within area of interest, maximum absolute value of slope $=9,$ a rather small value for slope.

Consequently, with only 9 passes through loop, Newton's method produces a result accurate to 200 places of decimals . Template:RoundBoxBottom Template:RoundBoxTop There are 3 conditions, any 1 of which terminates the loop:

• abs(y) very close to 0 (normal termination).

• count expired.

• endless loop detected with the same value of x repeated.

Calculation of cube roots of W

Calculation of 1 root

# python code.

def complexCubeRoot (a,b) :
    '''
p,q = complexCubeRoot (a,b) where:
* a,b are rectangular coordinates of complex number a + bi.
* (p + qi)^3 = a + bi.
'''

    a = Decimal(str(a))
    b = Decimal(str(b))

    if a == b == 0 : return (Decimal(0), Decimal(0))
    if a == 0 : return (Decimal(0), -simpleCubeRoot(b))
    if b == 0 : return (simpleCubeRoot(a), Decimal(0))

    Wmod = (a*a + b*b).sqrt()
    wmod = simpleCubeRoot (Wmod)
    cosWφ = a/Wmod
    coswφ = cosAfrom_cos3A(cosWφ)
    p = coswφ * wmod # wreal
    q = (wmod*wmod - p*p).sqrt() # wimag
    # Resolve ambiguity of q, + or -.
    v1 = - 3*p*p*q + q*q*q
    qp = abs(b + v1 )
    qn = abs(b - v1)
    if qp > qn : q *= -1

    return p,q

For function simpleCubeRoot() see Cube_root#Implementation Template:RoundBoxBottom

Calculation of 3 roots

Template:RoundBoxTop See Cube roots of unity.

The cube roots of unity are : $1,\frac{-1\pm i\sqrt{3}}{2}.$

When $r_0=\sqrt[3]{W}$ is known, the other 2 cube roots are:

• $r_1=r_0\cdot \frac{-1+i\sqrt{3}}{2}$

• $r_2=r_0\cdot \frac{-1-i\sqrt{3}}{2}$

# python code
def complexCubeRoots (a,b) :
    '''
r0,r1,r2 = complexCubeRoots (a,b) where:
* a,b are rectangular coordnates of complex number a + bi.
* (p + qi)^3 = a + bi.
* r0 = (p0,q0)
* r1 = (p1,q1)
* r2 = (p2,q2)
'''
    p,q = complexCubeRoot (a,b)
    r3 = Decimal(3).sqrt() # Square root of 3.
    pr3,qr3 = p*r3, q*r3
#  r1 = ((-p-q*r)/2, (p*r - q)/2)
#  r2 = ((-p+q*r)/2, (-p*r - q)/2)
    r0 = (p,q)
    r1 = ((-p-qr3)/2, (pr3 - q)/2)
    r2 = ((-p+qr3)/2, (-pr3 - q)/2)
    return r0,r1,r2

Testing results

Template:RoundBoxTop $(p+q\cdot i{)}^3=p^3-3p{q}^2+(3p^{2}q-q^3)\cdot i=a+b\cdot i.$

$a=p^3-3p{q}^2;b=3p^{2}q-q^3.$

# Python code used to test results of code above,

def complexCubeRootsTest (a,b) :
    print ('\n++++++++++++++++++++')
    print ('a,b =', a,b)
    almostZero = Decimal('1e-' + str(prec-5))
    r0,r1,r2 = complexCubeRoots (a,b)
    for root in (r0,r1,r2) :
        p,q = root
        print ('  pq =',(p), (q))
        a_,b_ = (p*p*p - 3*p*q*q), (3*p*p*q - q*q*q)
        if a :
            v1 = abs((a_-a)/a)
            if v1 > almostZero : print ('error *',a_,a,v1)
        else :
            v1 = abs(a_)
            if v1 > almostZero : print ('error !',a_,a,v1)
        if b :
            v1 = abs((b_-b)/b)
            if v1 > almostZero : print ('error &',b_,b,v1)
        else :
            v1 = abs(b_)
            if v1 > almostZero : print ('error %',b_,b,v1)
    return r0,r1,r2
    
import decimal
Decimal = D = decimal.Decimal
prec = decimal.getcontext().prec # precision

cosAfrom_cos3Adebug = 1

for p in range (-10,11,1) :
    for q in range (-10,11,1) :
        a = p*p*p - 3*p*q*q
        b = 3*p*p*q - q*q*q
        complexCubeRootsTest(a,b)

Gallery

• Cube roots of 8i.
  Cube roots of 8i.
• Cube roots of -8.
  Cube roots of -8.
• 3 cube roots of W.
  3 cube roots of W.

Method #2 (Graphical)

Introduction

File:0406 2curves01a.png

Graphs of 2 curves showing complex cube roots as points of intersection of the 2 curves.

Let complex number $w=p+qi.$

Then $W=w^3=(p+qi{)}^3=p^3-3p{q}^2+(3p^{2}q-q^3)i.$

Let $W=a+bi.$

Then:

$a=p^3-3p{q}^2$ and

$b=3p^{2}q-q^3.$

When $W$ is given and $w$ is desired, $w$ may be calculated from the solutions of 2 simultaneous equations:

$p^3-3p{q}^2-a=0\dots (1)$ and

$3p^{2}q-q^3-b=0\dots (2).$

For example, let $W=(39582+3799i).$

Then equations $(1)$ and $(2)$ become (for graphical purposes):

$x^3-3x{y}^2-39582=0\dots (3),$ black curve in diagram, and

$3x^{2}y-y^3-3799=0\dots (4),$ red curve in diagram.

Three points of intersection of red and black curves are:

$(-18,29),$

$(34.11473670974872,1.0884572681198943),$ and

$(-16.11473670974872,-30.088457268119896),$

interpreted as the three complex cube roots of $W,$ namely:

$w_0=(-18+29i),$

$w_1=(34.11473670974872+1.0884572681198943i)$ and

$w_2=(-16.11473670974872-30.088457268119896i).$ Template:RoundBoxTop Proof:

# python code
# Each cube root cubed.
>>> w0 = (-18 + 29j)
>>> w1 = (34.11473670974872 + 1.0884572681198943j)
>>> w2 = (-16.11473670974872 - 30.088457268119896j)
>>> for w in (w0,w1,w2) : w**3
... 
(39582+3799j)
(39582+3799j)
(39582+3799j)
# The moduli of all 3 cube roots.
>>> for w in (w0,w1,w2) : (w.real**2 + w.imag**2) ** 0.5
... 
34.132096331752024
34.132096331752024
34.132096331752024

Preparation

This method depends upon selection of most appropriate quadrant.

In the example above, quadrant $2$ is chosen because any non-zero positive value of $y$ intersects red curve and any non-zero negative value of $x$ intersects black curve.

Figures 1-4 below show all possibilities of $\pm a$ and $\pm b.$ Template:RoundBoxTop

• Figure 1. When both a,b are positive, quadrant 2 is chosen. (−18+29i)3=(39582+3799i)
  Figure 1. When both $a,b$ are positive, quadrant 2 is chosen.
  $(-18+29i{)}^3$$=(39582+3799i)$
• Figure 2. When a is positive and b is negative, quadrant 3 is chosen. (−18−29i)3=(39582−3799i)
  Figure 2. When $a$ is positive and $b$ is negative, quadrant 3 is chosen.
  $(-18-29i{)}^3$$=(39582-3799i)$
• Figure 3. When a is negative and b is positive, quadrant 1 is chosen. (18+29i)3=(−39582+3799i)
  Figure 3. When $a$ is negative and $b$ is positive, quadrant 1 is chosen.
  $(18+29i{)}^3$$=(-39582+3799i)$
• Figure 4. When both a,b are negative, quadrant 4 is chosen. (18−29i)3=(−39582−3799i)
  Figure 4. When both $a,b$ are negative, quadrant 4 is chosen.
  $(18-29i{)}^3$$=(-39582-3799i)$

Description of method

Four points

File:0406 4points.png

Graphs of 2 curves showing 4 points that enclose one of the complex cube roots.
Area enclosed by the 4 points becomes progressively smaller and smaller until the point of intersection is identified.

Assume that $W=-39582-3799i,$ in which case both $a,b$ are negative and quadrant $4$ is chosen.

In quadrant $4$ any non-zero positive value of x intersects black curve and any non-zero negative value of y intersects red curve.

Choose any convenient negative, non-zero value of $y.$

Let $y=-18.$

Using this value of $y,$ calculate coordinates of point $A$ on red curve.

Using $x$ coordinate of point $A,$ calculate coordinates of point $B$ on black curve.

Using $y$ coordinate of point $B,$ calculate coordinates of point $C$ on red curve.

Using $x$ coordinate of point $C,$ calculate coordinates of point $D$ on black curve.

Points $A,B,C,D$ enclose the point of intersection of the 2 curves. Template:RoundBoxBottom

Point of intersection

File:0406intersection.png

Graphs of 2 curves showing complex cube root enclosed by four points $A,B,C,D$.
Point $E,$ intersection of lines $AC,BD$ is close to complex cube root, and is starting point for next iteration.

Calculate equations of lines $AC,BD.$

Calculate coordinates of point $E,$ intersection of lines $AC,BD.$

Point $E$ is used as starting point for next iteration.

Area of quadrilateral $ABCD$ becomes smaller and smaller until complex cube root, intersection of red and black curves, is identified. Template:RoundBoxBottom

Implementation

Initialization

# python code

import decimal
D = decimal.Decimal

precision = decimal.getcontext().prec = 100

useDecimal = 1

Tolerance = 1e-14
if useDecimal :
    Tolerance = D('1e-' + str(decimal.getcontext().prec - 2))


def line_mc (point1,point2) :
    '''
m,c = line_mc (point1,point2)
where y = mx + c.
'''
    x1,y1 = point1
    x2,y2 = point2
    m = (y2-y1) / (x2-x1)
    # y = mx + c
    c = y1 - m*x1
    return m,c


def intersectionOf2Lines (line1, line2) :
    m1,c1 = line1
    m2,c2 = line2
    # y = m1x + c1
    # y = m2x + c2
    # m1x + c1 = m2x + c2
    # m1x - m2x = c2 - c1
    x = (c2-c1)/(m1-m2)
    y = m1*x + c1
    return x,y

def almostEqual (v1,v2,tolerance = Tolerance) :
    '''
status = almostEqual (v1,v2)

For floats, tolerance is 1e-14.

1234567.8901234567 and
1234567.8901234893
are not almostEqual.

1234567.8901234567 and
1234567.8901234593
are almostEqual.
'''
    return abs(v1-v2) < tolerance*abs((v1+v2)/2)

Function two_points()

File:0406 4points.png

Graphs of 2 curves showing 4 points that enclose one of the complex cube roots.
On first invocation function two_points() returns points $A,B.$
On second invocation function two_points() returns points $C,D.$
If distance $AB$ or distance $CD$ is very small, function two_points() returns status True.

def two_points (y,a,b,quadrant) :
    '''
[point1,point2],status = two_points (y,a,b)
'''

    L1 = [] ; yInput = y
    if 0 :
        print ('two_points():')
        s1 = '  a,b,y' ; print (s1, eval(s1))

# a = ppp - pqq - 2pqq = ppp - 3pqq (1)
# b = 2ppq + ppq - qqq = 3ppq - qqq (2)
#
# Using (2)
# 3xxy - yyy = b
#
#     yyy + b
# X = ----------
#        3y
    X = (y**3 + b) / (3*y)

    if isinstance(X,D) : x = X.sqrt()
    else : x = X ** 0.5
    if quadrant in (2,3) : x *= -1
    L1 += [(x,y)]

# Using (1)
# xxx - 3xyy = a
#
#     xxx - a
# Y = -----------
#        3x
#
    Y = (x**3 - a) / (3*x)

    if isinstance(Y,D) : y = Y.sqrt()
    else : y = Y ** 0.5
    if quadrant in (3,4) : y *= -1

    L1 += [(x,y)]

    return L1, almostEqual(y, yInput)

Function pointOfIntersection()

File:0406intersection.png

Graphs of 2 curves showing complex cube root enclosed by four points $A,B,C,D$.
From points $A,B,C,D$ function pointOfIntersection() calculates coordinates of point $E.$
If distance $AB$ is very small, point $A$ is returned as equivalent to intersection of red and black curves.
If distance $CD$ is very small, point $C$ is returned as equivalent to intersection of red and black curves.

def pointOfIntersection(y,a,b,quadrant) :
    '''
pointE, status = pointOfIntersection(y,a,b)
y is Y coordinate of point A.
'''
#    print('\npointOfIntersection()')
    t1,status = two_points (y,a,b,quadrant)
    ptA,ptB = t1
    if status :
        # Distance between ptA and ptB is very small.
        # ptA is considered equivalent to the complex cube root.
        return ptA,status

    t2,status = two_points (ptB[1],a,b,quadrant)
    ptC,ptD = t2
    if status :
        # Distance between ptC and ptD is very small.
        # ptC is considered equivalent to the complex cube root.
        return ptC,status

    lineAC = line_mc (ptA,ptC)
    lineBD = line_mc (ptB,ptD)
    pointE = intersectionOf2Lines (lineAC, lineBD)
    return pointE,False

Execution

def CheckMake_pq(a,b,x,y) :
    '''
p,q = CheckMake_pq(a,b,x,y)
p = x and q = y.
p,q are checked as valid within tolerance, and are reformatted slightly to improve appearance.
'''
#    print ('\nCheckMake_pq(a,b,x,y)')
#    s1 = '(a,b)' ; print (s1,eval(s1))
#    s1 = '   x'  ; print (s1,eval(s1))
#    s1 = '   y'  ; print (s1,eval(s1))
    P = x**2 ; Q = y**2 ; p=q=-1
    for p in (x,) :
        # a = ppp - 3pqq; a + 3pqq should equal ppp.
        if not almostEqual (P*p, a + 3*p*Q, Tolerance*100) : continue
        for q in  (y,) :
            # b = 3ppq - qqq; b + qqq should equal 3ppq
            if not almostEqual (3*P*q, b + q*Q) : continue
            # Following 2 lines improve appearance of p,q
            if useDecimal :
                # 293.00000000000000000000000000000000000000034 becomes 293
                p,q = [ decimal.Context(prec=precision-3).create_decimal(s).normalize()
                        for s in (p,q) ]
            else :
                # 123.99999999999923 becomes 124.0
                p,q = [ float(decimal.Context(prec=14).create_decimal(s)) for s in (p,q) ]
            return p,q
    # If code gets to here there is internal error.
    s1 = '   p'  ; print (s1,eval(s1))
    s1 = '   q'  ; print (s1,eval(s1))
    s1 = 'P*p, a + 3*p*Q' ; print (s1,eval(s1))
    s1 = '3*P*q, b + q*Q' ; print (s1,eval(s1))
    1/0


def ComplexCubeRoot (a,b, y = 100, count_ = 20) :
    '''
p,q = ComplexCubeRoot (a,b, y, count_)
(p+qi)**3 = (a+bi)
'''
    print ('\nComplexCubeRoot(): a,b,y,count_ =',a,b,y,count_)

    if useDecimal : y,a,b = [ D(str(v)) for v in (y,a,b) ]

    if a == 0 :
        if b == 0 : return 0,0
        if useDecimal : cubeRoot = abs(b) ** (D(1)/3)
        else : cubeRoot = abs(b) ** (1/3)
        if b > 0 : return 0,-cubeRoot
        return 0,cubeRoot
    if b == 0 :
        if useDecimal : cubeRoot = abs(a) ** (D(1)/3)
        else : cubeRoot = abs(a) ** (1/3)
        if a > 0 : return cubeRoot,0
        return -cubeRoot,0

    # Select most appropriate quadrant.
    if a > 0: setx = {2,3}
    else: setx = {1,4}

    if b > 0: sety = {1,2}
    else: sety = {3,4}

    quadrant, = setx & sety
    # Make sign of y appropriate for this quadrant.
    if quadrant in (1,2) : y = abs(y)
    else : y = -abs(y)

    s1 = '    quadrant' ;    print (s1, eval(s1))

    for count in range (0,count_) :
        pointE,status = pointOfIntersection(y,a,b,quadrant)
        s1 = 'count,status' ; print (s1, eval(s1))
        s1 = '    pointE[0]' ;    print (s1, eval(s1))
        s1 = '    pointE[1]' ;    print (s1, eval(s1))
        x,y = pointE
        if status : break

    p,q = CheckMake_pq(a,b,x,y)

    return p,q

An Example

p,q = 18,-29
w0 = p + q*1j
W = w0**3
a,b = W.real, W.imag
s1 = '\na,b' ; print (s1, eval(s1))
print ('Calculate one cube root of W =', W)

p,q = ComplexCubeRoot (a, b, -18)
s1 = '\np,q' ; print (s1, eval(s1))

sign = ' + '
if q < 0 : sign = ' - '
print ('w0 = ', str(p) ,sign, str(abs(q)),'i',sep='')

a,b (-39582.0, -3799.0)
Calculate one cube root of W = (-39582-3799j)

ComplexCubeRoot(): a,b,y,count_ = -39582.0 -3799.0 -18 20
    quadrant 4
count,status (0, False)
    pointE[0] 18.29530595866769796981147594954794157453427770441979517949705190002312860122683512802090262517713985
    pointE[1] -29.17605851188829785804056660826030733025475591125914094664311767722817017040959484906193571713185723
count,status (1, False)
    pointE[0] 18.00005338608833140244623119091867294731079673210031643698475740639013316464725974710029414039192178
    pointE[1] -29.00004871113589733281025047965490410760310240487028781197782902118493613318468122514940700337153822
count,status (2, False)
    pointE[0] 18.00000000000411724901243622639913339568271402799883998577879934397861645683113192688877413853607310
    pointE[1] -29.00000000000374389488957142528693701977412078931714190614174861872117809632376941851192785888688849
count,status (3, False)
    pointE[0] 18.00000000000000000000000002432197332441306371168312765664226520285586754485200564671348858119116700
    pointE[1] -29.00000000000000000000000002211642401405406173668734741733917293863102998696594080783163691859719835
count,status (4, False)
    pointE[0] 18.00000000000000000000000000000000000000000000000000000084875301166228708732099292145747224660296019
    pointE[1] -29.00000000000000000000000000000000000000000000000000000077178694502906372553368739653233322951192943
count,status (5, False)
    pointE[0] 18.00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
    pointE[1] -29.00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
count,status (6, True)
    pointE[0] 18
    pointE[1] -29.00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

p,q (Decimal('18'), Decimal('-29'))
w0 = 18 - 29i

Method #3 (Algebraic)

Introduction

Let complex number $w=p+qi.$

Then $W=w^3=(p+qi{)}^3=p^3-3p{q}^2+(3p^{2}q-q^3)i.$

Let $W=a+bi.$

Then:

$a=p^3-3p{q}^2\dots (1)$ and

$b=3p^{2}q-q^3\dots (2).$

When $W$ is given and $w$ is desired, $w$ may be calculated from the solutions of 2 simultaneous equations $(1)$ and $(2),$ where $a,b$ are known values, and $p,q$ are desired.

Implementation

$(1)$ squared: $a^2=p^6-6p^4{q}^2+9p^2{q}^4\dots (1a)$

From $(2):3p^{2}q=b+q^3\dots (2a)$

$(1a)*27q^3:$

$27q^3{a}^2=27q^3{p}^6-27q^3(6)p^4{q}^2+27q^3(9)p^2{q}^4$

$27q^3{a}^2=27(p^{2}q{)}^3-27(6)p^4{q}^5+27(9)p^2{q}^7$

$27q^3{a}^2=(3p^{2}q{)}^3-27(6)p^4{q}^2{q}^3+27(9)p^{2}q{q}^6$

$27q^3{a}^2=(3p^{2}q{)}^3-3(6)(9p^4{q}^2)q^3+27(3)(3p^{2}q)q^6$

$27q^3{a}^2=(3p^{2}q{)}^3-3(6)((3p^{2}q{)}^2)q^3+27(3)(3p^{2}q)q^6$

Let $Q=q^3:$

$27Q{a}^2=(3p^{2}q{)}^3-3(6)((3p^{2}q{)}^2)Q+27(3)(3p^{2}q)Q^2\dots (1b)$

For $(3p^{2}q)$ in $(1b)$ substitute $(b+Q):$

$27Q{a}^2=(b+Q{)}^3-3(6)((b+Q{)}^2)Q+27(3)(b+Q)Q^2\dots (1c)$

Expand $(1c),$ simplify, gather like terms and result is:

$f(Q)=s{Q}^3+t{Q}^2+uQ+v\dots (3)$ where:

$Q=q^3$

$s=64$

$t=48b$

$u=-(15b^2+27a^2)$

$v=b^3$

Calculate one real root of $(3):Q_1$

$q_1=\sqrt[3]{Q_1}$

From $(2a):p_1=\sqrt{\frac{b+Q_1}{3q_1}}$

Check $p_1$ against $(1)$ to resolve ambiguity of sign of $p_1.$

$p_1+q_{1}i$ is one cube root of $W.$

Template:RoundBoxTop $p$ may be calculated without ambiguity as follows:

$a=p^3-3p{q}^2\dots (1)$ and

$b=3p^{2}q-q^3\dots (2).$

From $(1):p^3-3q^{2}p-a=0\dots (3)$

From $(2):3q{p}^2-(Q+b)=0\dots (4)$

Let:

$A=-3q^2$

$B=-a$

$C=3q$

$D=-(Q+b)$

Then $(3),(4)$ become:

$p^3+Ap+B=0\dots (5)$

$C{p}^2+D=0\dots (6)$

$(5)*D:D{p}^3+DAp+DB=0\dots (7)$

$(6)*B:BC{p}^2+BD=0\dots (8)$

$(7)-(8):D{p}^3-BC{p}^2+DAp=0\dots (9)$

Simplify $(9):D{p}^2-BCp+DA=0\dots (10)$

Repeat $(6):C{p}^2+D=0\dots (6)$

$(10)*C:CD{p}^2-BCCp+DAC=0\dots (11)$

$(6)*D:CD{p}^2+DD=0\dots (12)$

$(11)-(12):-BCCp+DAC-DD=0\dots (13)$

From $(13):p=\frac{DAC-D^2}{B{C}^2}=\frac{D(AC-D)}{B{C}^2}$ Template:RoundBoxBottom

An Example

Calculate cube roots of complex number $W=39582+3799i.$ Template:RoundBoxTop

File:1219cubic01.png

Graph of $f(Q)$ shown as graph of $f(x)$ and showing three values of $Q:Q_1,Q_2,Q_3$.
$Y$ axis compressed for clarity.

# python code:

a,b = 39582,3799

s = 64
t = 48*b
u = -(15*b**2 + 27*a**2)
v = b**3

s,t,u,v

(64, 182352, -42518323563, 54828691399)

Calculate roots of cubic function: $y=f(x)$$=64x^3$$+182352x^2$$-42518323563x$$+54828691399.$

Three roots are: $Q_1,Q_2,Q_3=-27239.53953801976,1.2895380197588122,24389$ Template:RoundBoxBottom

# python code:

values_of_Q = Q1,Q2,Q3 = -27239.53953801976, 1.2895380197588122, 24389
q1,q2,q3 = [ abs(Q)**(1/3) for Q in values_of_Q ]
q1 *= -1
values_of_q = q1,q2,q3
s1 = 'values_of_q' ; print(s1, eval(s1))

for m in zip(values_of_q, values_of_Q) :
    q,Q = m
    p = ((b + Q)/(3*q)) ** .5
    sum = p**3 - 3*p*q**2 - a
    if abs(sum) > 1e-10 : p *= -1
    w = p + q*1j
    s1 = 'w, w**3' ; print (s1, eval(s1))

values_of_q (-30.08845726811989, 1.0884572681198943, 29)
w, w**3 ((-16.114736709748723 - 30.08845726811989j  ), (39582+3799j))
w, w**3 (( 34.11473670974874  +  1.0884572681198943j), (39582+3799j))
w, w**3 ((-18.0               + 29.0j               ), (39582+3799j))

Three cube roots of $W=39582+3799i$ are:

$w_1=(-16.114736709748723-30.08845726811989i)$

$w_2=(34.11473670974874+1.0884572681198943i)$

$w_3=(-18.0+29.0i)$