0

I am lookin for a solution to find cubic roots in Excel. I found the below code at this website.

http://www.mrexcel.com/forum/excel-questions/88804-solving-equations-excel.html

unfortunately, it doesn't work for me - I get #VALUE! when I run it and since I am only learning VBA, I have not had luck debugging it.

Sub QUBIC(P As Double, Q As Double, R As Double, ROOT() As Double)

' Q U B I C - Solves a cubic equation of the form:
' y^3 + Py^2 + Qy + R = 0 for real roots.
' Inputs:
' P,Q,R Coefficients of polynomial.

' Outputs:
' ROOT 3-vector containing only real roots.
' NROOTS The number of roots found. The real roots
' found will be in the first elements of ROOT.

' Method: Closed form employing trigonometric and Cardan
' methods as appropriate.

' Note: To translate and equation of the form:
' O'y^3 + P'y^2 + Q'y + R' = 0 into the form above,
' simply divide thru by O', i.e. P = P'/O', Q = Q'/O',
' etc.

Dim Z(3) As Double
Dim p2 As Double
Dim RMS As Double
Dim A As Double
Dim B As Double
Dim nRoots As Integer
Dim DISCR As Double
Dim t1 As Double
Dim t2 As Double
Dim RATIO As Double
Dim SUM As Double
Dim DIF As Double
Dim AD3 As Double
Dim E0 As Double
Dim CPhi As Double
Dim PhiD3 As Double
Dim PD3 As Double

Const DEG120 = 2.09439510239319
Const Tolerance = 0.00001
Const Tol2 = 1E-20

' ... Translate equation into the form Z^3 + aZ + b = 0

p2 = P ^ 2
A = Q - p2 / 3
B = P * (2 * p2 - 9 * Q) / 27 + R

RMS = Sqr(A ^ 2 + B ^ 2)
If RMS < Tol2 Then
' ... Three equal roots
nRoots = 3
ReDim ROOT(0 To nRoots)
For i = 1 To 3
ROOT(i) = -P / 3
Next i
Exit Sub
End If

DISCR = (A / 3) ^ 3 + (B / 2) ^ 2

If DISCR > 0 Then

t1 = -B / 2
t2 = Sqr(DISCR)
If t1 = 0 Then
RATIO = 1
Else
RATIO = t2 / t1
End If

If Abs(RATIO) < Tolerance Then
' ... Three real roots, two (2 and 3) equal.
nRoots = 3
Z(1) = 2 * QBRT(t1)
Z(2) = QBRT(-t1)
Z(3) = Z(2)
Else
' ... One real root, two complex. Solve using Cardan formula.
nRoots = 1
SUM = t1 + t2
DIF = t1 - t2
Z(1) = QBRT(SUM) + QBRT(DIF)
End If

Else

' ... Three real unequal roots. Solve using trigonometric method.
nRoots = 3
AD3 = A / 3#
E0 = 2# * Sqr(-AD3)
CPhi = -B / (2# * Sqr(-AD3 ^ 3))
PhiD3 = Acos(CPhi) / 3#
Z(1) = E0 * Cos(PhiD3)
Z(2) = E0 * Cos(PhiD3 + DEG120)
Z(3) = E0 * Cos(PhiD3 - DEG120)

End If

' ... Now translate back to roots of original equation
PD3 = P / 3

ReDim ROOT(0 To nRoots)

For i = 1 To nRoots
ROOT(i) = Z(i) - PD3
Next i

End Sub

Function QBRT(X As Double) As Double

' Signed cube root function. Used by Qubic procedure.

QBRT = Abs(X) ^ (1 / 3) * Sgn(X)

End Function

Can anyone please guide me on how to fix it, so I can run it. Thanks.

EDIT: This is how I am running it in Excel (I changed Qubic to be a function instead of sub) cells A1:A3 contain p,q, r respectively cells B1:B3 contain Roots() cells C1:C3 contain array for the output of Qubic

A1:1 A2:1 A3:1

B1:0.1 B2:0.1 B3:0.1

C1: C2: C3: {=QUBIC(A1,A2,A3,B1:B3)}

ADD: now that it works with the fix from @assylias, I am trying the following from another sheet:

Function ParamAlpha(p,q,r) as Double
Dim p as Double
Dim q as Double 
Dim r as Double
p=-5
q=-2
r=24
    Dim Alpha as Double
    Dim AlphaVector() as Double
    AlphaVector=QubicFunction(p,q,r)
    Alpha=FindMinPositiveValue(AlphaVector)
End Function

Function FindMinPositiveValue(AlphaVector) As Double
Dim N As Integer, i As Integer
N = AlphaVector.Cells.Count
Dim Alpha() As Double
ReDim Alpha(N) As Double
For i = 1 To N
    If AlphaVector(i) > 0 Then
        Alpha(i) = AlphaVector(i)
    Else
        Alpha(i) = 100000000000#
    End If
Next i
FindMinPositiveValue = Application.Min(Alpha)
End Function

In Excel, I call =ParamAlpha(-5,-2,24) and it returns #VALUE!

Fionnuala
  • 90,370
  • 7
  • 114
  • 152
user1155299
  • 877
  • 5
  • 20
  • 29

1 Answers1

2

If you add the following procedure, it will show the results in a message box. You can then modify it to do something else as you require:

Public Sub test()

  Dim p As Double
  Dim q As Double
  Dim r As Double
  Dim roots() As Double

  p = 1
  q = 1
  r = 1

  QUBIC p, q, r, roots

  Dim i As Long
  Dim result As String

  result = "("
  For i = LBound(roots, 1) To UBound(roots, 1)
    result = result & roots(i) & ","
  Next i

  result = Left(result, Len(result) - 1) & ")"

  MsgBox "Roots of y^3 + " & p & ".y^2 + " & r & ".y + " & r & " = 0 has the following roots: " & result

End Sub

Alternatively, if you want the result in the form of a fomula array directly in a spreadsheet, you can add the following function in the same module:

Public Function QubicFunction(p As Double, q As Double, r As Double) As Double()

  Dim roots() As Double
  QUBIC p, q, r, roots
  QubicFunction = roots

End Function

You then call it from Excel by selecting a few cells (horizontally, for example A1:B1) and press CTRL+SHIFT+ENTER:

=QubicFunction(1, 1, 1)
assylias
  • 321,522
  • 82
  • 660
  • 783
  • thanks for the code. I added it but it doesn't show any messages. I also changed QUBIC to function from Sub (should that make a difference?) and it didn't work. I am updating my post with how I am running it in Excel. Please let me know if I am running it incorrectly. thanks! – user1155299 Nov 19 '12 at 16:43
  • @user1155299 Open a new workbook, open the VBA editor (ALT+F11), right click on the new workbook's project > Insert > New Module and paste your code (and mine) in that new module. – assylias Nov 19 '12 at 16:46
  • did that, and I still get #VALUE!. I updated my post with how I am running it in Excel. Is how I am doing it correct? – user1155299 Nov 19 '12 at 16:51
  • 1
    @user1155299 What do you mean I get `#VALUE!`. How do you run the macro? From the excel sheet, press ALT+F8, select the `Test` macro and run it. – assylias Nov 19 '12 at 16:53
  • @user1155299 I see that you call it as a function from Excel. I have added something at the bottom of my answer to do that. – assylias Nov 19 '12 at 16:59
  • oh ok. it works, thanks. the first window show P=1, and then the last message box `Roots of....following roots:(0,-1)`. So, the question is: why do I get a #VALUE! when I enter the parameters into QUBIC in Excel and run it as {=QUBIC(p,q,r,roots)}. where am I making a mistake in that? and how do I correct it so I can access it as a function in any spreadsheet? oops, I just saw your post, let me run that – user1155299 Nov 19 '12 at 17:05
  • I just updated my post to get the minimum positive value of the cubic root. The function `FindMinPositiveValue` works fine stand-alone, but from Excel when I call `=ParamAlpha(p,q,r)`, then it returns `#VALUE!`. Why so? The individual functions work fine, but I think I am making a mistake in calling it. Can you please tell me how to fix it. Thanks. – user1155299 Nov 19 '12 at 20:02
  • 1
    With `N = UBound(AlphaVector, 1)` it should work better (once you remove the p, q and r declaration and assignment in ParamAlpha) – assylias Nov 19 '12 at 20:23
  • I made the change, and it works. Thanks again. I should read up on UBound to understand what happened there ;-) – user1155299 Nov 19 '12 at 20:48