mrootsort

Option Explicit

Public Function SortRootResult(ByVal RootA As Variant, _
                               ByVal Degree As Long, _
                               Optional ByVal Tol As Double = 0.0000000001) As Variant
    Dim i As Long, j As Long
    Dim tmpR As Variant, tmpI As Variant
    Dim nReal As Long
    Dim nRows As Long, nCols As Long

    On Error GoTo ReturnResult

    nRows = UBound(RootA, 1)
    nCols = UBound(RootA, 2)

    If Degree < 1 Then GoTo ReturnResult
    If nRows < Degree Then GoTo ReturnResult
    If nCols < 2 Then GoTo ReturnResult

    ' Treat very small imaginary parts as zero
    For i = 1 To Degree
        If Abs(DblOrZero(RootA(i, 2))) <= Tol Then
            RootA(i, 2) = 0#
        End If
    Next i

    ' Sort roots
    For i = 1 To Degree - 1
        For j = i + 1 To Degree
            If RootComesAfter(RootA(i, 1), RootA(i, 2), _
                              RootA(j, 1), RootA(j, 2), Tol) Then

                tmpR = RootA(i, 1)
                tmpI = RootA(i, 2)

                RootA(i, 1) = RootA(j, 1)
                RootA(i, 2) = RootA(j, 2)

                RootA(j, 1) = tmpR
                RootA(j, 2) = tmpI
            End If
        Next j
    Next i

    ' Recount real and complex roots
    nReal = 0
    For i = 1 To Degree
        If Abs(DblOrZero(RootA(i, 2))) <= Tol Then
            RootA(i, 2) = 0#
            nReal = nReal + 1
        End If
    Next i

    ' Count row is normally Degree + 1
    If nRows >= Degree + 1 Then
        RootA(Degree + 1, 1) = nReal
        RootA(Degree + 1, 2) = Degree - nReal
    End If

ReturnResult:
    SortRootResult = RootA
End Function

Private Function RootComesAfter(ByVal rA As Variant, ByVal iA As Variant, _
                                ByVal rB As Variant, ByVal iB As Variant, _
                                ByVal Tol As Double) As Boolean
    Dim ar As Double, ai As Double
    Dim br As Double, bi As Double
    Dim aIsReal As Boolean, bIsReal As Boolean

    ar = DblOrZero(rA)
    ai = DblOrZero(iA)
    br = DblOrZero(rB)
    bi = DblOrZero(iB)

    If Abs(ai) <= Tol Then ai = 0#
    If Abs(bi) <= Tol Then bi = 0#

    aIsReal = (ai = 0#)
    bIsReal = (bi = 0#)

    ' Real roots must come before complex roots
    If aIsReal <> bIsReal Then
        RootComesAfter = Not aIsReal
        Exit Function
    End If

    ' Sort by real part ascending
    If Abs(ar - br) > Tol Then
        RootComesAfter = (ar > br)
        Exit Function
    End If

    ' If complex roots have the same real part,
    ' positive imaginary part comes first
    If Not aIsReal Then
        If Abs(ai - bi) > Tol Then
            RootComesAfter = (ai < bi)
            Exit Function
        End If
    End If

    RootComesAfter = False
End Function

Private Function DblOrZero(ByVal v As Variant) As Double
    On Error GoTo NotNumeric

    If IsEmpty(v) Then GoTo NotNumeric
    If IsNull(v) Then GoTo NotNumeric

    If IsNumeric(v) Then
        DblOrZero = CDbl(v)
    Else
        DblOrZero = 0#
    End If

    Exit Function

NotNumeric:
    DblOrZero = 0#
End Function