Lotto Algorithms - Permutations, Combinations

1. Introduction

This thread aims to consolidate previous contributions and recent developments, in order to provide efficient and useful functions for analysing Permutations and Combinations.

We understand these tools are used by Lotto-style game analysers, and while we personally feel this is an exercise whose profitability (in fiscal terms) is extremely dubious from a purely mathematical point of view .... :huh:

....... nevertheless!! :rolleyes: these people have every right to decent algorithms - and they might still make lots of lotto money - (remember us if you do, ok? ;) )

Those seeking simple "generators" for complete or partial sets of combinations or permutations will find simple functions below ("Odometer-style Generators").

Those needing to do more than that will want to look at the slightly more complex "Index Mapping" functions, which can be used to generate sets but also provide additional capabilities.


In order to understand the difference between the two methods, we need to understand the concept of how all possible combinations (or perms) can be listed in a natural sort order called the lexicographical order.

Here are the combinations of 3 from 5 in lexicographical order:

#1 = 1, 2, 3
#2 = 1, 2, 4
#3 = 1, 2, 5
#4 = 1, 3, 4
#5 = 1, 3, 5
#6 = 1, 4, 5
#7 = 2, 3, 4
#8 = 2, 3, 5
#9 = 2, 4, 5
#10 = 3, 4, 5


As you can see, it's really just the natural numerical sort order, based on the leftmost fields having precedence.

Permutations also have a natural lexicographical (sequence) order. Here is the complete listing for K = 4:

#1 = 1, 2, 3, 4
#2 = 1, 2, 4, 3
#3 = 1, 3, 2, 4
#4 = 1, 3, 4, 2
#5 = 1, 4, 2, 3
#6 = 1, 4, 3, 2
#7 = 2, 1, 3, 4
#8 = 2, 1, 4, 3
#9 = 2, 3, 1, 4
#10 = 2, 3, 4, 1
#11 = 2, 4, 1, 3
#12 = 2, 4, 3, 1
#13 = 3, 1, 2, 4
#14 = 3, 1, 4, 2
#15 = 3, 2, 1, 4
#16 = 3, 2, 4, 1
#17 = 3, 4, 1, 2
#18 = 3, 4, 2, 1
#19 = 4, 1, 2, 3
#20 = 4, 1, 3, 2
#21 = 4, 2, 1, 3
#22 = 4, 2, 3, 1
#23 = 4, 3, 1, 2
#24 = 4, 3, 2, 1


Index Mapping

By creating a MAPPING function (CombinationToIndex) that assigns an INDEX value to each combination or permutation according to its position in the full list, we can represent individual combinations by simple numbers.

With an inverse function (IndexToCombination) that returns the combination for a given index value, we can then generate any subset of all possible combinations from a given starting point like this:

StartIndex = CombinationToIndex(50, 6, "1, 7, 23, 35, 47, 49")
For CombNo = StartIndex To StartIndex + 99
Debug.Print IndexToCombination(50, 6, CombNo)
Next



1, 7, 23, 35, 47, 49
1, 7, 23, 35, 47, 50
1, 7, 23, 35, 48, 49
1, 7, 23, 35, 48, 50
1, 7, 23, 35, 49, 50
1, 7, 23, 36, 37, 38
1, 7, 23, 36, 37, 39
....
....



Odometer-style Generators

These are like odometers, you keep a table with a "wheel" for each element in the set. You start it at "1 2 3 4 5 ...", and then you "rotate" it with successive calls to a NextCombination function which advances the wheels to the next position.

SetCombination 50, 6, "1, 7, 23, 35, 47, 49" ' set the wheels
For i = 1 To 100
Debug.Print ThisCombination ' returns current reading
NextCombination 50, 6 ' advance to next position
Next


The rules for rotating the wheels are more complex than for an odometer (which allows duplicate values), but there are wonderfully simple combinatorial algorithms for doing both Perms and Combs.


All VB functions described will be presented below, but first we'll look at the pros and cons of the 2 approaches.

2. Comparison of the 2 Algorithms

The odometer approach is both concise and efficient, so it wins the "elegance" award, and it also can handle HUGE games (like a million balls from a billion!).

However, for analytical purposes, there is a fatal flaw - you can't tell "distance" on this odometer - you have no idea where in the list you are from the current settings.

So the odometer can't tell you the answer to questions like this:

How many combinations are there between "1, 7, 23, 35, 47, 49" and "11, 14, 23, 27, 47, 49" in a 6-from-50 game?


What's the first combination in the second half of the full (50, 6) list?


To answer these requires the index mapping functions. The first answer is about 11 million, as it turns out:

?combination2index(50,6,"11, 14, 23, 27, 47, 49")
12,231,760
?combination2index(50,6,"1, 7, 23, 35, 47, 49")
926,277


And the second question, the halfway mark?

?combinations(50,6)
15,890,700
? 15890700 / 2
7,945,350
?index2combination(50, 6, 7945351)
6, 8, 19, 20, 36, 37


The mapping functions allow us to divide the full virtual set into a nice regular grid, and section it accordingly. Generating sequential list members is slightly less efficient than the Odometer, but still a tiny fraction of the cost of either printing or otherwise displaying or (worse) storing them all.


Scaling - Maximum Numbers Supported

Because the Odometer method doesn't require any measurement of the size of the full possible combination set, it can work with games of ridiculous size.

For example, let's play "10 balls from 50,000!!":

setcombination 5000,10,"12, 34, 56, 123, 456, 789, 1234, 2345, 3456, 4567"
?thiscombination
12, 34, 56, 123, 456, 789, 1234, 2345, 3456, 4567
nextcombination 5000, 10
?thiscombination
12, 34, 56, 123, 456, 789, 1234, 2345, 3456, 4568


Now you can't do that with index mapping functions because you simply can't even store the numbers required (in the case of a 10-from-50000 game the number of combinations needs a million or so digits to represent!).

But that's only an artificial restriction, as the standard data types in VB give us plenty of scope to model all current games of interest, plus room to expand seamlessly to bigger games like "40 from 100", in fact any game for which the total set size (number of possible combinations) is less than 75 x 10 ^ 27.

This is quite good - you can choose any K (number of balls to pick) for all values of N up to 99. At N = 100 you can only go to K = 46 (just shy of the finishing line at K=50) . But you'd normally use smaller K anyway, so that lets you increase N even further - e.g. you can model 10-from-1000

You'll note that 75 x 10^27 is rather larger than the maximum LONG value in VB. The index mapping functions use Variant (Decimal) variables where required, which gives us access to what are effectively triple-length Longs' (96-bit integers).

We used this technique to extend the range of the original function that we posted to calculate the number of Combinations(N, K) 129902

And that's while you'll see some Variant variables in the Index Mapping functions below, but not the Odometer ones.

3. Odometer Functions for Combinations

If you've followed everything above, you should find this module self-contained and the function usage self-explanatory! ;)


Option Explicit
'
' Odometer-style Combination generating functions
' by MathImagics (Dr Memory) 2004
'
' SetCombination
' NextCombination
' ThisCombination
'
Dim Cwheel() As Long

Public Sub SetCombination(ByVal N As Long, ByVal K As Long, ByVal Combn As String)

ReDim token(K) As String
token = Split(Combn, ",")

If UBound(token) <> K - 1 Then Exit Sub

Dim W As Long

ReDim Cwheel(0 To K)

For W = 1 To K
Cwheel(W) = Val(token(W - 1))
If Cwheel(W) <= Cwheel(W - 1) Then Exit Sub ' invalid combn
If Cwheel(W) > N Then Exit Sub ' ditto
Next
End Sub

Public Function ThisCombination() As String
'
' Current Combination readout
'
Dim i As Long, Comb As String
Comb = Cwheel(1)
For i = 2 To UBound(Cwheel)
Comb = Comb & ", " & Cwheel(i)
Next
ThisCombination = Comb
End Function


Public Sub NextCombination(ByVal N As Long, ByVal K As Long)
' "Combination Odometer"
'
' By MathImagics: the array Cwheel contains the current
' K items combined, in increasing order.
' Each call to this sub will adjust Cwheel
' so it contains
' the NEXT combination in lex order
'

Dim i As Long
Dim j As Long
i = K
While Cwheel(i) >= N - K + i
i = i - 1 ' find rightmost wheel that allows an increment
If i = 0 Then
' wraps around (natch!)
i = 1
Cwheel(1) = 0
End If
Wend
Cwheel(i) = Cwheel(i) + 1
For j = i + 1 To K
Cwheel(j) = Cwheel(i) + j - i
Next
End Sub

4. Odometer for Permutations

Using similar principles, we can enumerate Permutations of K objects, starting with any given point, for any practical value of K.


Option Explicit
'
' Odometer-style Combo/Permn functions
' by MathImagics (Dr Memory) 2004
'
' SetPermutation
' NextPermutation
' ThisPermutation
'
Dim Pwheel() As Long

Public Sub SetPermutation(ByVal K As Long, ByVal Perm As String)
ReDim token(K) As String, flag(1 To K) As Boolean
token = Split(Perm, ",")

If UBound(token) <> K - 1 Then Exit Sub

Dim W As Long, i As Long

For W = 1 To K
i = Val(token(W - 1))
If i < 1 Or i > K Then Exit Sub ' Invalid permn
If flag(i) Then Exit Sub ' ditto
flag(i) = True
Pwheel(W) = i
Next
End Sub

Sub NextPermutation(ByVal K As Long)

' "Permutation Odometer"
'
' By MathImagics: the array Pwheel contains the values
' 1 to K in any order. Each call to this
' sub will adjust Pwheel so it contains
' the NEXT permutation in lex order
'
' Note: like NextCombn, this also cycles from last to
' first

Dim i As Long, j As Long

j = K - 1
Do While Pwheel(j) > Pwheel(j + 1)
j = j - 1
If j = 0 Then Exit Do
Loop

If j Then
i = K
Do While Pwheel(j) > Pwheel(i)
If i = 1 Then Exit Do
i = i - 1
Loop
GoSub SwapIJ
End If

i = K
j = j + 1

While i > j
GoSub SwapIJ
i = i - 1
j = j + 1
Wend
Exit Sub

SwapIJ: ' swap wheel positions i and j
Dim temp As Long
temp = Pwheel(i)
Pwheel(i) = Pwheel(j)
Pwheel(j) = temp
Return
End Sub

Public Function ThisPermutation() As String
'
' Current Permutation readout
'
Dim i As Long, Perm As String
Perm = Pwheel(1)
For i = 2 To UBound(Pwheel)
Perm = Perm & ", " & Pwheel(i)
Next
ThisPermutation = Perm
End Function

5. Index mapping functions for Combinations

These are self-contained, and use no external variables. One optimisation you might consider, however, if heavy use of the functions proves slow, is simply to pre-store the required Combinations(N,K) values (easily predicted) in a lookup table, which should speed it up dramatically.


Option Explicit

Function Index2Combination( _
ByVal N As Long, ByVal K As Long, ByVal CNO As Variant) As String
'
' by Dr Memory (MathImagics) May 2004
'
' N = # of values (balls)
' K = # selected to combine
' CNO = index number of combination to return
' 1 <= CNO <= Combinations(N, K) < 75 x 10^27
'
'=============================================
Dim NC, V, cCount ' Variants (Decimal)
Dim Wheel As Long, W As Long, Combo As String

NC = Combinations(N, K)
If CNO < 1 Or CNO > NC Then Exit Function

Wheel = 1

cCount = CDec(0) ' running count of combinations

For W = 1 To K - 1
Do
V = Combinations(N - Wheel, K - W) ' how many till next value?
If V = 0 Then Exit Do
If cCount + V >= CNO Then Exit Do
cCount = cCount + V
Wheel = Wheel + 1
Loop
If W = 1 Then Combo = Wheel Else Combo = Combo & ", " & Wheel
Wheel = Wheel + 1
Next
Index2Combination = Combo & ", " & Wheel + (CNO - cCount - 1)
End Function


Function Combination2Index(ByVal N As Long, ByVal K As Long, Combn As String) As Variant
'
' by Dr Memory (MathImagics) May 2004
'
' N = # of values (balls)
' K = # selected to combine
' Combn= list of values, e.g. "2, 7, 9, 23, 40, 45"
' must be in ascending order
'=============================================

ReDim token(K) As String
token = Split(Combn, ",")

If UBound(token) <> K - 1 Then Exit Function

Dim NC, wcount, Wheel() As Long, W As Long
Dim V As Long, msg As String
ReDim Wheel(K)
Combination2Index = CDec(0)
For W = 1 To K
Wheel(W) = Val(token(W - 1))
If Wheel(W) <= Wheel(W - 1) Then Exit Function
If Wheel(W) > N Then Exit Function
Next

For W = 1 To K - 1
For V = Wheel(W - 1) + 1 To Wheel(W) - 1
Combination2Index = Combination2Index + Combinations(N - V, K - W) - 1
Next
Next
Combination2Index = Combination2Index + Wheel(K) - K + 1
End Function

6. Index mapping functions for Permutations

Here's the index functions for Permutations. Note that the same optimisation described for Combinations might be worthwhile (i.e. prestore Factorial() values in a lookup table).


Function Permutation2Index(ByVal K As Long, Perm As String) As Variant
'
' by Dr Memory (MathImagics) May 2004
'
' K = # items to Permute
'
' Permn= list of values, e.g. "2, 7, 9, 23, 40, 45"
' must be in ascending order
'=============================================

ReDim token(K) As String, flag(K) As Boolean
token = Split(Perm, ",")

If UBound(token) <> K - 1 Then Exit Function

Dim V As Long, W As Long, i As Long, iPrev As Long, j As Long

V = CDec(0)

For W = 1 To K
i = Val(token(W - 1))
If i < 1 Or i > K Then Exit Function
If flag(i) Then Exit Function
flag(i) = True
Select Case W
Case 1: If i > 1 Then V = CDec(i - 1) * Factorial(K - 1) ' 1st digit index is absolute
Case K - 1: iPrev = i
Case K: If i > iPrev Then V = V + CDec(1) Else V = V + CDec(2)
Case Else
'
' for intermediate wheel positions, need to know the relative index
' of the wheel value given for this position. This depends on the
' items already used, but they are flag()-ed for easy id
'
Dim ix As Long
ix = 0: j = 1
Do
Do While flag(j)
If j = K Then Exit Do
j = j + 1
Loop
If j >= i Then Exit Do
ix = ix + 1
j = j + 1
Loop
If ix Then V = V + CDec(ix) * Factorial(K - W) ' i.e. Factorial(k-w)
End Select
Next
Permutation2Index = V
End Function

Function Index2Permutation(ByVal K As Long, ByVal CNO) As String
'
' by Dr Memory (MathImagics) May 2004
'
' K = # of items to permute
' CNO = index number of permutation to return
' 1 <= CNO <= Factorial(K-1)
'
'=============================================
Dim NC, V, cCount ' Variants (Decimal)
Dim Wheel As Long, W As Long, Perm As String


NC = Factorial(K)
If CNO < 1 Or CNO > NC Then Exit Function

Static flag(30) As Boolean ' we can't handle K>27 anyway
Erase flag

cCount = CDec(0) ' running count of combinations

For W = 1 To K
'
' set first available value for next wheel
'
Wheel = 1
While flag(Wheel): Wheel = Wheel + 1: Wend
If W < K Then
Do
V = Factorial(K - W)
If cCount + V >= CNO Then Exit Do
cCount = cCount + V
Do
Wheel = Wheel + 1
Loop While flag(Wheel)
Loop
flag(Wheel) = True
End If

If W = 1 Then Perm = Wheel Else Perm = Perm & ", " & Wheel
Next
Index2Permutation = Perm
End Function

Function Factorial(ByVal N As Long)
'
' MathImagics - Good for results to 96-bits (N <= 27)
'
If N <= 2 Then Factorial = CDec(N): Exit Function
If N > 27 Then Exit Function ' no can do!
Factorial = CDec(N) * Factorial(N - 1)
End Function

7. Alternative Algorithm for Generating Permutations

There is a very interesting alternative algorithm for permutation generating - it's called the "Johnson-Trotter" method. I'll refer to it here as JTS (Johnson-Trotter Shuffle), and compare it to the method used in section 4 above (ODO method).

JST is very fast, much faster than ODO, but it doesn't produce the permutations in lexical order like ODO does, instead it produces a sequence called the "Transpositional Order" (see here (http://www.theory.csc.uvic.ca/~cos/inf/perm/PermInfo.html) for a useful picture of this).

The reason for its speed advantage is that it always produces the next permutation by making just one positional swap, and it has a very clever method of deciding which pair should be swapped at each iteration.


Example
The permutations of order 4 in LEX order are:
1234,1243,1324,1342,1423,1432,
2134,2143,2314,2341,2413,2431,
3124,3142,3214,3241,3412,3421,
4123,4132,4213,4231,4312,4321
The permutations of order 4 in JTS order are:
1234,1243,1423,4123,4132,1432,
1342,1324,3124,3142,3412,4312,
4321,3421,3241,3214,2314,2341,
2431,4231,4213,2413,2143,2134

The Johnson-Trotter Algorithm

The N items are represented by the numbers 1 to N.

We also use a LINK array - at any stage, each item is "linked" with either its left-hand or right-hand neighbour. LINK(j)=0 indicates ITEM(j) is "linked" with ITEM(j-1), LINK(j)=1 indicates it's "linked" with ITEM(j+1).

At each iteration, we look for what are caled "mobile" items (items we are allowed to move, ie. swap). An item in position J is regarded as mobile if and only if its linked item exists, and is smaller in value. (The "exists" condition means that an item in position 1 can only be mobile if its LINK value is 1, and an item in position N can only be mobile if its LINK value is 0).

We begin with the items placed in ascending order in the ITEM array, and all LINK values set to 0, and compute the next permutation as follows:

Find the mobile item with the highest item value.
Let K = this value, and J be its item position. Swap item J with its neighbour, either J-1 or J+1, depending on LINK(J). Note that we swap both the ITEM and LINK values.
Toggle the LINK indicator for any position T where ITEM(T) > K.

We simply repeat this sequence until there are no mobile items at all, which means we have reached the end of the sequence.


A VB6 encoding is given below.

8. The Johnson-Trotter algorithm in VB6

This is a self-contained version, set to list the permutations in the debug window (cf: Step 1).
Sub GenPermutations(ByVal N As Long)
'
' "Johnson-Trotter" VB6 implementation by MathImagics (Dec 2004)
' Each permutation is obtained from the previous by
' swapping just ONE pair of adjacent items.
'

Dim Item() As Long ' items to permute
Dim Link() As Long ' 0 = link left, 1 = right
Dim j As Long
Dim K As Long, kSpot As Long ' largest mobile K and its position
Dim P As Long, pSpot As Long ' iterator value P, its position
Dim mobile As Boolean ' "mobility" test flag
Dim kLink As Long
'
' 0. Setup initial state
'
Redim Item(N), Link(N)
For j = 1 To N
Item(j) = j
Next

Do
'
' 1. report current permutation
'
Debug.Print Item(1);
For j = 2 To N: Debug.Print ","; Item(j);: Next
Debug.Print
'
' 2. select "mobile" position with highest value
'
K = 0
pSpot = 0

Do While pSpot < N
pSpot = pSpot + 1
P = Item(pSpot)

mobile = False

If Link(pSpot) = 0 Then
If pSpot > 1 Then
If Item(pSpot - 1) < P Then mobile = True
End If
ElseIf pSpot < N Then
If Item(pSpot + 1) < P Then mobile = True
End If

If mobile Then
If P > K Then
K = P
kSpot = pSpot
If K = N Then Exit Do ' look no further
End If
End If
Loop

If K = 0 Then Exit Do ' all done!

'
' 3. Swap item kSpot with "neighbour"
'
kLink = Link(kSpot)
If kLink Then
Item(kSpot) = Item(kSpot + 1): Link(kSpot) = Link(kSpot + 1)
Item(kSpot + 1) = K: Link(kSpot + 1) = 1
Else
Item(kSpot) = Item(kSpot - 1): Link(kSpot) = Link(kSpot - 1)
Item(kSpot - 1) = K: Link(kSpot - 1) = 0
End If
'
' 4. Toggle Links for any items > K
'
For pSpot = 1 To N
If Item(pSpot) > K Then Link(pSpot) = 1 - Link(pSpot)
Next
Loop
End Sub