The present code snippet just does an especific task. This is, given an ordered sequence of natural numbers, starting at zero, generate all the different permutations there can be. When we talk about permutations, the order of the sequence does matter. So, the sequence 0, 1, 2 differs from 2, 0, 1.
For example, think about a 4-7-2 lock combination sequence. Changing the sequence to 2-7-4 would not open the lock, of course. Properly speaking in mathematics, we are talking about permutations, because the order does matter. On the other hand, if the order did not matter, we would be talking about combinations.
So, what does present Console application do? First, asks for how large the sequence must be. Possible values of N must be greater or equal to 3 and less or equal to 10 (3<=N<=10). Then asks if the sequence generation must be saved into a text file and/or displayed. After the sequence is generated, requests if the generated N! permutations should be verified. All of them should be different, this is, there must be no 2 or more equal sequences. If you play around with the application, changing just one of the numbers of the variable vTable()(), for instance, changing New Int32() {1, 2, 3} into New Int32() {1, 1, 3} implies failing any generation N>=4. Any change to the next array New Int32() {3, 1, 3, 1} would imply failing for N>=5.
The table is taken from http://comjnl.oxfordjournals.org/content/6/3/293.full.pdf and it is linked at Wikipedia http://en.wikipedia.org/wiki/Heap's_algorithm
I have to be honest, because I did not follow too much the explanations and flow diagrams. So I had to think by myself for a way to make the algorithm happen and work. If you find any failure, please let me know.
As you may see, the algorithm is not recursive. There is an iteration process and exits when i>= N, i.e. when the right-most number is going to be interchanged for the N-th time. If you like to dig into the code, first take a look to the above links.
Just let me comment a bit the array vIndex(). When a new sequence is generated vIndex(1) is incremented from 0 to 1. The next time from 1 to 2 but because of the tailing ...mod (i+1), here mod(2), vIndex(1) will be 0. So i is incremented to 2 and so vIndex(2)=1. The process continues and when vIndex(2) equals 3, a zero is assigned to vIndex(2) and vIndex(3) is incremented by one (i=3). Then the "if i>2 then" instruction branches so that the positions interchanged follow the rules given by the table inside vTable()() (see Table 1 at Heap's .pdf document).

9a92ed49a9e81b164b3edf1e804fb9a6

Imports System.IO

Module Module1

    Sub Main()
        Try
            Console.WriteLine("Numbers to permutate (3<= N <= 10), enter N= ? ")
            Dim e1 As String = Console.ReadLine
            Dim N As Int32
            If Not Int32.TryParse(e1, N) Then
                Console.WriteLine("Could not parse N.")
                Exit Try
            End If
            If N < 3 OrElse N > 10 Then
                Console.WriteLine("N is out of bounds.")
                Exit Try
            End If
            Console.WriteLine("Save to file (Y/N)?: ")
            e1 = Console.ReadLine
            Dim bSave As Boolean = False
            If Len(e1) AndAlso LCase(e1.Chars(0) = "y") Then
                bSave = True
            End If
            Console.WriteLine("Should permutations be displayed (Y/N)?: ")
            e1 = Console.ReadLine
            Dim bDsp As Boolean = False
            If Len(e1) AndAlso LCase(e1.Chars(0)) = "y" Then
                bDsp = True
            End If
            Console.WriteLine("")

            Dim sFile As String = String.Empty
            Dim perm()() As Int32 = _
                HeapsPermByInterchange(N, bDsp, bSave, sFile)
            Console.WriteLine(vbCrLf + "Done!" + vbCrLf)
            If bSave Then
                Console.WriteLine("Saved to file: " + _
                        vbCrLf + sFile + vbCrLf)
            End If


            Console.WriteLine("Verify (Y/N)?: ")
            e1 = Console.ReadLine
            If Len(e1) AndAlso LCase(e1.Chars(0)) <> "y" Then
                Exit Try ' no verification
            End If
            Dim msg As String = String.Empty
            If Not verifyPermutations(perm, N, msg) Then
                Console.WriteLine(msg)
            Else
                Console.WriteLine( _
                vbCrLf + "Permutations have been successfully verified!" _
                + vbCrLf)
            End If
        Catch ex As Exception
            Console.WriteLine(ex.ToString)
        End Try
        Console.WriteLine("Press Enter to exit.")
        Console.ReadLine()
    End Sub
    Function HeapsPermByInterchange(n As Int32, _
                    Optional bDsp As Boolean = False, _
                    Optional bSave As Boolean = False, _
                    Optional ByRef filePath As String = "") As Int32()()

        ' See http://en.wikipedia.org/wiki/Heap's_algorithm
        ' and also http://comjnl.oxfordjournals.org/content/6/3/293.full.pdf

        Static vTable()() As Int32 = { _
            New Int32() {1, 2, 3}, _
            New Int32() {3, 1, 3, 1}, _
            New Int32() {3, 4, 3, 2, 3}, _
            New Int32() {5, 3, 1, 5, 3, 1}, _
            New Int32() {5, 2, 7, 2, 1, 2, 3}, _
            New Int32() {7, 1, 5, 5, 3, 3, 7, 1}, _
            New Int32() {7, 8, 1, 6, 5, 4, 9, 2, 3}, _
            New Int32() {9, 7, 5, 3, 1, 9, 7, 5, 3, 1}, _
            New Int32() {9, 6, 3, 10, 9, 4, 3, 8, 9, 2, 3}}
        Dim ret()() As Int32 = Nothing
        Dim fs As FileStream = Nothing
        Dim sw As StreamWriter = Nothing
        Try
            Dim vIndex(n) As Int32
            Dim c As Int32 = 0
            Dim aux, iAux, viAux(n) As Int32
            Dim i As Int32
            Dim fact As Int32 = n
            For i = n - 1 To 2 Step -1
                fact *= i
            Next
            Dim nDigits As Int32 = Math.Floor(Math.Log10(fact)) + 1
            Dim sFormat As String = "{0:" _
                + Microsoft.VisualBasic.StrDup(nDigits, "0") + "}:  "
            ReDim ret(fact - 1)
            Dim cI(n - 1) As Int32
            For i = 0 To fact - 1
                If i < n Then cI(i) = i
                ReDim ret(i)(n - 1)
            Next
            If bSave Then
                filePath = Microsoft.VisualBasic.CurDir _
                    + "\perm" + n.ToString + ".txt"
                fs = New FileStream(filePath, FileMode.Create)
                sw = New StreamWriter(fs)
            End If
            Dim iOld As Int32
            Do
                If bDsp Then
                    Console.Write(String.Format(sFormat, c + 1))
                End If
                For i = 0 To cI.Length - 1
                    If bDsp Then
                        Console.Write(cI(i).ToString + " ")
                    End If
                    If bSave Then
                        sw.Write(cI(i).ToString + " ")
                    End If
                    ret(c)(i) = cI(i)
                Next
                If bDsp Then
                    If c Then
                        Console.Write(" --- ")
                        If iOld > 2 Then
                            Console.WriteLine((iOld + 1).ToString + _
                               "," + (iAux + 1).ToString)
                        Else
                            Console.WriteLine("1," _
                                    + (iOld + 1).ToString)
                        End If
                    Else
                        Console.WriteLine("")
                    End If
                End If
                If bSave Then
                    sw.WriteLine("")
                End If
                c += 1
                i = 1
                Do
                    vIndex(i) = (vIndex(i) + 1) Mod (i + 1)
                    If vIndex(i) Then Exit Do
                    i += 1
                    If i >= n Then Exit Try
                Loop
                iOld = i
                If i > 2 Then
                    iAux = viAux(i)
                    iAux = vTable(i - 3)(iAux) - 1
                    aux = cI(iAux)
                    cI(iAux) = cI(i) : cI(i) = aux
                    viAux(i) = (viAux(i) + 1) Mod i
                Else
                    aux = cI(0)
                    cI(0) = cI(i) : cI(i) = aux
                End If
            Loop
        Catch ex As Exception
            Throw ex
        Finally
            If sw IsNot Nothing Then
                sw.Close()
            End If
            If fs IsNot Nothing Then
                fs.Close()
            End If
        End Try
        Return ret
    End Function
    Public Function verifyPermutations( _
        perm()() As Int32, N As Int32, _
        Optional ByRef msg As String = "") As Boolean
        Try
            Dim i, j As Int32
            Dim fact As Int32 = N
            For i = N - 1 To 2 Step -1
                fact *= i
            Next
            Dim vStr(fact - 1) As String
            For i = 0 To fact - 1
                For j = 0 To N - 1
                    vStr(i) += Chr(&H30 + perm(i)(j))
                Next
            Next
            Array.Sort(vStr)
            For i = fact - 1 To 2 Step -1
                Dim pos As Int32 = Array.BinarySearch(vStr, vStr(i))
                If pos >= 0 AndAlso pos < i Then
                    Exit For
                End If
            Next
            If i > 1 Then
                j = Array.IndexOf(vStr, vStr(i))
                msg = String.Format( _
                    "rows {0} and {1} are repeated: ", _
                    i + 1, j + 1)
                msg += " ("
                For j = 0 To N - 1
                    msg += perm(i)(j).ToString
                    If j < N - 1 Then
                        msg += ", "
                    End If
                Next
                msg += ")"
                Return False
            End If
        Catch ex As Exception
            Throw ex
        End Try
        Return True
    End Function
End Module

I think this is a case where you'll find the recursive approach much more clear and concise.

Module Module1

    Sub Main()

        Permute({1, 2, 3, 4})
        Console.ReadLine()

    End Sub

    Sub Permute(nums() As Integer, Optional start As Integer = 0)

        If start = nums.Length Then
            OutputArray(nums)
        Else
            For i As Integer = start To nums.Length - 1
                Swap(nums, start, i)
                Permute(nums, start + 1)
                Swap(nums, start, i)
            Next
        End If

    End Sub

    Sub Swap(ByRef nums() As Integer, i As Integer, j As Integer)

        Dim temp As Integer
        temp = nums(i)
        nums(i) = nums(j)
        nums(j) = temp

    End Sub

    Sub OutputArray(nums() As Integer)

        For Each num As Integer In nums
            Console.Write(num.ToString & " ")
        Next
        Console.WriteLine()

    End Sub

End Module

Edited 2 Years Ago by Reverend Jim

Here is the reason I was working with permutations: to obtain the determinant and the inverse of a matrix. The calculator already had a method to retrieve the inverse (Gauss-Jordan through absolute pivot), but it was extremely slow. For a 5x5 matrix it could take about one minute to resolve. Now, for the input:
roots(det(A))@A=(3,a,a,a,a|a,3,a,a,a|a,a,3,a,a|a,a,a,3,a|a,a,a,a,3)
it takes 203ms, more than 300 times faster! and, now a days, one minute is much more.
I'll explain the input:
1) Obtain the determinant of matrix A, i.e. a polynomial P(a) =(4a+3)(a-3)^4.
(corresponds to function det(A) in the input).
2) Obtain the roots of P(a) (specified in function roots(...))
3) Matrix A is a 5x5 matrix (corresponds to @A=.... )
3,a,a,a,a
a,3,a,a,a
a,a,3,a,a
a,a,a,3,a
a,a,a,a,3
The user's interface looks like this:

5539f93ebe10d791883b03ae26530d5c

The calculator already had a method to retrieve the inverse (Gauss-Jordan through absolute pivot), but it was extremely slow.

I can't see any algorithm that finds the inverse of a 5x5 matrix (assuming it is non-singular) taking as long as one minute. Using Gauss-Jordan on a 5x5 should be almost instantaneous. For larger matrices the usual technique to reduce time (at least the one I was taught) is to use full or partial virtual pivoting.

I don't recall anything from university about determinants (too ancient).

Edited 2 Years Ago by Reverend Jim

You're right again. I've been very imprecise, mostly trying to use few words. When previous version, 8.3.35, tries to calculate A^-1 through Gauss-Jordan times out at 1 minute, the longer configurable time. I suspect this is due to an error, because the polynomials degree involved increase too much. Calculations not only involve numbers but polynomials, although this is not a reason to take so long.
I'll try to fix the method and come back with timing as, in theory, Gauss-Jordan should be quicker.

Really, there was no error but because of the way I programmed Gauss-Jordan method the polynomial degrees of the matrix's entries increased in excess. By finding the G.C.D., if any, the degree may be reduced. Anyway, as can see in the image, Gauss-Jordan (8.3.85) takes 234ms and through the determinant (8.3.86) 94ms.
For the determinant, in 8.3.86, it follows Leiniz formula. (See Leibinz formula at http://en.wikipedia.org/wiki/Determinant)

047ae745985e14fdf723cfc222464eed

Both verifications are:
93bc51ccc2fe7d06bea3676acd1377d7

and:
baa498f38c212ea2a53cb27b8c6a6fc6

This snippet code is the responsible to obtain the determinant, cofactor, adjoint and inverse matrix. A expression class instance may hold a number, polynomial or any expression like, for example, 'sin(x)' in an AST tree. A exprMatrix class is nothing else than a matrix where the entries are expression class instances. The complete code is downloadable here

    Public Function opDeterminant( _
                Optional bGetInverseOrCofactorMtx As Boolean = False, _
                Optional ByRef InvMtx As ExprMatrix = Nothing, _
                Optional ByRef cofactorMtx As ExprMatrix = Nothing, _
                Optional ByRef adjoint As ExprMatrix = Nothing) As Expression
        Dim vTable()() As Int32 = { _
            New Int32() {1, 2, 3}, _
            New Int32() {3, 1, 3, 1}, _
            New Int32() {3, 4, 3, 2, 3}, _
            New Int32() {5, 3, 1, 5, 3, 1}, _
            New Int32() {5, 2, 7, 2, 1, 2, 3}, _
            New Int32() {7, 1, 5, 5, 3, 3, 7, 1}, _
            New Int32() {7, 8, 1, 6, 5, 4, 9, 2, 3}, _
            New Int32() {9, 7, 5, 3, 1, 9, 7, 5, 3, 1}, _
            New Int32() {9, 6, 3, 10, 9, 4, 3, 8, 9, 2, 3}}
        Dim retDetermExpr As Expression = Nothing
        Dim i, j As Int32
        Dim N As Int32 = Me.Rows
        Dim vIndex(N) As Int32
        Dim c As Int32 = 0
        Dim aux, iAux, viAux(N) As Int32
        Dim fact As Int32 = N
        Dim cI(N - 1), cInv(N - 1) As Int32
        Try
            ' See http://en.wikipedia.org/wiki/Heap's_algorithm
            ' and also http://comjnl.oxfordjournals.org/content/6/3/293.full.pdf

            ' see Leibniz formula at:
            ' http://en.wikipedia.org/wiki/Determinant

            Dim oVar As New VarsAndFns(Me.cfg)

            For i = N - 1 To 2 Step -1
                fact *= i
            Next
            For i = 0 To N - 1
                cI(i) = i
            Next

            retDetermExpr = New Expression(0.0)
            If True Then
                Do
                    Dim sgn As Int32 = 1
                    Dim curr As New Expression(1.0)
                    For i = 0 To cI.Length - 1
                        For j = 0 To i - 1
                            If cI(j) > cI(i) Then sgn *= -1
                        Next
                        curr *= getExpr(i, cI(i))
                    Next
                    If sgn = 1 Then
                        retDetermExpr += curr
                    Else
                        retDetermExpr -= curr
                    End If
                    c += 1
                    i = 1
                    If N = 1 Then Exit Try
                    Do
                        vIndex(i) = (vIndex(i) + 1) Mod (i + 1)
                        If vIndex(i) Then Exit Do
                        i += 1
                        If i >= N Then Exit Try
                    Loop
                    If i > 2 Then
                        iAux = viAux(i)
                        iAux = vTable(i - 3)(iAux) - 1
                        aux = cI(iAux)
                        cI(iAux) = cI(i) : cI(i) = aux
                        viAux(i) = (viAux(i) + 1) Mod i
                    Else
                        aux = cI(0)
                        cI(0) = cI(i) : cI(i) = aux
                    End If
                Loop
            End If
        Catch ex As Exception
            Throw ex
        End Try
        Try
            If bGetInverseOrCofactorMtx Then
                Dim Adj(N - 1, N - 1) As ExprMatrix
                If N = 1 Then
                    Adj(0, 0) = New ExprMatrix(Me.getExpr(0, 0))
                Else
                    For iRow As Int32 = 0 To N - 1
                        For iCol As Int32 = 0 To N - 1
                            ' Get exprMatrix Adj(,) out from
                            ' 'me' where row 'iRow' and column 'iCol'
                            ' are being excluded:
                            Dim i1 As Int32 = 0
                            For i = 0 To N - 1
                                If i <> iRow Then ' exclude row = iRow
                                    Dim j1 As Int32 = 0
                                    For j = 0 To N - 1
                                        If j <> iCol Then ' exclude col = iCol
                                            If Adj(iRow, iCol) Is Nothing Then
                                                Adj(iRow, iCol) = New ExprMatrix( _
                                                    cfg, N - 1, N - 1)
                                            End If
                                            If (i + j) Mod 2 Then
                                                Adj(iRow, iCol).getExpr(i1, j1) = _
                                                    New Expression(-Me.getExpr(i, j))
                                            Else
                                                Adj(iRow, iCol).getExpr(i1, j1) = _
                                                    New Expression(Me.getExpr(i, j))
                                            End If
                                            j1 += 1
                                        End If
                                    Next
                                    i1 += 1
                                End If
                            Next
                        Next
                    Next
                End If
                ' Obtain all the determinants of Adj(,)
                ' and assing to InvMtx:
                InvMtx = New ExprMatrix(cfg, N, N)
                For i = 0 To N - 1
                    For j = 0 To N - 1
                        InvMtx.getExpr(i, j) = Adj(i, j).opDeterminant
                    Next
                Next
                cofactorMtx = InvMtx
                adjoint = ExprMatrix.opTranspose(InvMtx)
                ' Finally, transpose InvMtx and divide each
                ' entry by the determinant (retExpr):
                InvMtx = ExprMatrix.opTranspose(InvMtx, retDetermExpr)
            End If

        Catch ex As Exception
            Throw ex
        End Try

        Return retDetermExpr
    End Function

Edited 2 Years Ago by xrj: missing link