Although the Leibniz formula is inefficient for calculation in general I decided to use it, which meant I needed to create permutations of n items, σ, and the sign of each permutation
Swapping signs
The sign of a permutation is defined as the parity of the number of swaps between any two items of a permutation needed to generate that permutation from the initial permutation.So for example, given a starting perm of [0,1,2] the permutation [2,0,1] can be achieved by first swapping the 2 and the 0 to form [2,1,0], then swapping the 1 and 0 to form [2,0,1. That is two swaps in total. The sign of the permutation is given as -1 if an odd number of swaps are used in generating it otherwise as +1 if an even (or no) swaps are needed. (The sign of the permutation will not change even if other swaps, or redundant swaps are used to generate it from the same starting perm).
Pythons itertools.permutations function was of unknown sign, but I picked up on the Steinhaus–Johnson–Trotter algorithm that generates permutations in such a way that it guarantees that successive permutations only differ by the swapping of two items - this meant that the sign would start as +1 and alternate in sign between successive items.
The Steinhaus–Johnson–Trotter Python program
I decided to create a new Rosetta code task for generating both permutations and their sign with a view for their subsequent use in calculating determinants:from operator import itemgetter DEBUG = False # like the built-in __debug__ def spermutations(n): """permutations by swapping. Yields: perm, sign""" sign = 1 p = [[i, 0 if i == 0 else -1] # [num, direction] for i in range(n)] if DEBUG: print ' #', p yield tuple(pp[0] for pp in p), sign while any(pp[1] for pp in p): # moving i1, (n1, d1) = max(((i, pp) for i, pp in enumerate(p) if pp[1]), key=itemgetter(1)) sign *= -1 if d1 == -1: # Swap down i2 = i1 - 1 p[i1], p[i2] = p[i2], p[i1] # If this causes the chosen element to reach the First or last # position within the permutation, or if the next element in the # same direction is larger than the chosen element: if i2 == 0 or p[i2 - 1][0] > n1: # The direction of the chosen element is set to zero p[i2][1] = 0 elif d1 == 1: # Swap up i2 = i1 + 1 p[i1], p[i2] = p[i2], p[i1] # If this causes the chosen element to reach the first or Last # position within the permutation, or if the next element in the # same direction is larger than the chosen element: if i2 == n - 1 or p[i2 + 1][0] > n1: # The direction of the chosen element is set to zero p[i2][1] = 0 if DEBUG: print ' #', p yield tuple(pp[0] for pp in p), sign for i3, pp in enumerate(p): n3, d3 = pp if n3 > n1: pp[1] = 1 if i3 < i2 else -1 if DEBUG: print ' # Set Moving' if __name__ == '__main__': from itertools import permutations for n in (3, 4): print '\nPermutations and sign of %i items' % n sp = set() for i in spermutations(n): sp.add(i[0]) print('Perm: %r Sign: %2i' % i) #if DEBUG: raw_input('?') # Test p = set(permutations(range(n))) assert sp == p, 'Two methods of generating permutations do not agree'
Determinants using Leibniz's formula
I then had all I needed to create the determinant (and permanent) of a matrix code:
from itertools import permutations from operator import mul from math import fsum from spermutations import spermutations def prod(lst): return reduce(mul, lst, 1) def perm(a): n = len(a); r = range(n); s = list(permutations(r)) return fsum(prod(a[i][sigma[i]] for i in r) for sigma in s) def det(a): n = len(a); r = range(n); s = list(spermutations(n)) return fsum(sign*prod(a[i][sigma[i]] for i in r) for sigma, sign in s) if __name__ == '__main__': from pprint import pprint as pp for a in ( [ [1, 2], [3, 4]], [ [1, 2, 3, 4], [4, 5, 6, 7], [7, 8, 9, 10], [10, 11, 12, 13]], [ [ 0, 1, 2, 3, 4], [ 5, 6, 7, 8, 9], [10, 11, 12, 13, 14], [15, 16, 17, 18, 19], [20, 21, 22, 23, 24]], ): print('') pp(a) print('Perm: %s Det: %s' % (perm(a), det(a)))
Sign() of the times
I was still thinking about how to generate the sign of Pythons in-built permutation generator, and a little thought made me think that somewhere there was a sort routine that sorted by swapping two elements. If you count the swaps then you can calculate the sign if the number of swaps is even or odd.
This lead me to Selection Sort.
def selectionSort(lst): for i in range(0,len(lst)-1): mn = min(range(i,len(lst)), key=lst.__getitem__) lst[i],lst[mn] = lst[mn],lst[i] return lst
Unfortunately the sorter does unnecessary swaps so I made a slight modification to create a routine that generated a sign of any permutation of the integers 0..n:
def perm_parity(lst): '''\ Given a permutation of the digits 0..N in order as a list, returns its parity (or sign): +1 for even parity; -1 for odd. ''' parity = 1 for i in range(0,len(lst)-1): if lst[i] != i: parity *= -1 mn = min(range(i,len(lst)), key=lst.__getitem__) lst[i],lst[mn] = lst[mn],lst[i] return parity if __name__ == '__main__': from itertools import permutations for p in permutations(range(3)): l = list(p) print "%2i %r" % (perm_parity(l), p)
I published the sign generator on ActiveState as I didn't think it was enough for a Rosetta Code task.
(If you run the code you can see that their is in fact a pattern to the signs generated by successive members Pythons permutations).
Patterns
If you look at a few of the series generated by the Steinhaus–Johnson–Trotter Python program you can see patterns that point to a recursive method of generating the same sequences. I coded that up too.
Serendipity doo dah
A few days later I came across a Stack-exchange question that needed to generate all permutations when there were duplicate items in the initial list. (It goes on to need the permutations in a specific order, but I wasn't thinking of that bit at first).
Of coure I started by generating possible perm then using a set to remove duplicates, but it gets very (very) wasteful, so I started looking around for an algorithm to generate what are sometimes called lexicographically correct permutations. for example the l-perms of [1,1,1,2,2,2] are the twenty terms:
1 1 1 2 2 2
1 1 2 1 2 2
1 2 1 1 2 2
2 1 1 1 2 2
2 1 1 2 1 2
1 2 1 2 1 2
1 1 2 2 1 2
1 2 2 1 1 2
2 1 2 1 1 2
2 2 1 1 1 2
2 2 1 1 2 1
2 1 2 1 2 1
1 2 2 1 2 1
1 1 2 2 2 1
1 2 1 2 2 1
2 1 1 2 2 1
2 1 2 2 1 1
1 2 2 2 1 1
2 2 1 2 1 1
2 2 2 1 1 1
I spent an hour playing with my recursive permutation generator until I hit on the idea that if you had to do permutations of 123XX then you could think of each (identical) X as being X0 and X1, and that valid l-perms are the perms when X0 and X1 appear in only one order.
Now the recursive algorithm for perms introduces each item of the initial perm one-by-one. If you restrict each X as it is introduced to never be introduced 'above' and X that is already in the partial perm item that is being expanded, then, in effect, X0 and X1 can never 'cross' and so will stay in the same 'order'.
I made the change and it seems to work for me.
Lexicographic permutation generator in Python
def l_perm(items): if not items: return [[]] else: dir = 1 new_items = [] this = [items.pop()] for item in l_perm(items): lenitem = len(item) try: # Never insert 'this' above any other 'this' in the item maxinsert = item.index(this[0]) except ValueError: maxinsert = lenitem if dir == 1: # step down new_items += [item[:i] + this + item[i:] for i in range(lenitem, -1, -1) if i <= maxinsert] else: # step up new_items += [item[:i] + this + item[i:] for i in range(lenitem + 1) if i <= maxinsert] dir *= -1 return new_items if __name__ == '__main__': n = [1, 1, 1, 2, 2, 2] print '\nLexicograpic Permutations of %i items: %r' % (len(n), n) for i, x in enumerate(l_perm(n)): print('%3i %r' % (i, x))
End bit.
I'm tired and I want to go to bed so turning l-perm into a RC task or sticking it on Activestate can wait. I am wondering if their is anyone else interested in this as I just pottered around the less deep bits.- Paddy.