Friendsnippets

1. Various algorithms for computing the sign of a permutation.

   view plain - save this

   This is a Python version of the code (and associated commentary) from http://perplexus.info/show.php?pid=4895&op=sol You can run the source file through PyLit (http://pylit.berlios.de/) to get a nicely-formatted reStructuredText file.

     1 #! /usr/bin/env python
     2 #
     3 """Various algorithms for computing the sign of a permutation.
     4 
     5 Ported from http://perplexus.info/show.php?pid=4895&op=sol
     6 """
     7 __docformat__ = 'reStructuredText'
     8 
     9 
    10 ## perm_sign by John Burkardt is an example that can be found online
    11 ## (see it among the collection at
    12 ## http://orion.math.iastate.edu/burkardt/f_src/subset/subset.f90 or
    13 ## see http://www.scs.fsu.edu/~burkardt/math2071/perm_sign.m). A
    14 ## similar method, that, however, avoids explicit searching, is
    15 ## illustrated by the following function::
    16 
    17 def jbsign(p):
    18     """Compute sign of permutation `p` by counting the number of
    19     interchanges required to change the given permutation into the
    20     identity one.
    21 
    22     Examples::
    23       >>> jbsign([0,1,2])
    24       1
    25       >>> jbsign([0,2,1])
    26       -1
    27     """
    28     n = len(p)
    29     s = +1
    30     # get elements back in their home positions
    31     for j in xrange(0,n):
    32         q=p[j]
    33         if q !=j:
    34             p[j],p[q] = p[q],q # interchange p[j] and p[p[j]]
    35             s = -s             # and account for the interchange
    36     # note that q is now in its home position
    37     # whether or not an interchange was required
    38     return s
    39 
    40 
    41 ## A somewhat more sophisticated method follows the orbits of elements
    42 ## as they are permuted through cycles::
    43 
    44 def sign_by_orbits(p):
    45     """Compute sign of permutation `p` by following orbit cycles.
    46 
    47     Examples::
    48       >>> sign_by_orbits([0,1,2])
    49       1
    50       >>> sign_by_orbits([0,2,1])
    51       -1
    52     """
    53     n = len(p)
    54     s = +1
    55     # follow orbit cycles and tote up signs
    56     for j in xrange(0,n):
    57         if (p[j] >= 0):
    58             # new cycle, follow it, marking each place as "visited"
    59             q = j
    60             t = +1
    61             while (q >= 0):
    62                 r = p[q]
    63                 if (r >= 0):
    64                     p[q] = -1
    65                 t = -t
    66                 q = r
    67             # factor in contribution of each cycle
    68             s *= t
    69     return s
    70 
    71 
    72 ## Lang's Algebra, First Ed. 1965, p. 53, gives an equivalent
    73 ## definition of the sign as simply `(-1)^m` where m=n-number of
    74 ## orbits. Thus the above can be modified to just count the number of
    75 ## orbits and do a little arithmetic at the end, as the following
    76 ## function illustrates::
    77 
    78 def sign_lang_formula(p):
    79     """Compute sign of permutation `p` as `(-1)^{length of permutation
    80     - number of orbits}`.
    81 
    82     Examples::
    83       >>> sign_lang_formula([0,1,2])
    84       1
    85       >>> sign_lang_formula([0,2,1])
    86       -1
    87     """
    88     n = len(p)
    89     c = 0
    90     # count the orbit cycles
    91     for j in xrange(0,n):
    92         if (p[j] >= 0):
    93             # new cycle, follow it, marking each place
    94             q = j
    95             while (q >= 0):
    96                 r = p[q]
    97                 if (r >= 0):
    98                     p[q] = -1
    99                 q = r
   100             # increment cycle count
   101             c += 1
   102     if (n-c) % 2 == 0:
   103         return +1
   104     else:
   105         return -1
   106 
   107 
   108 ## Another method counts inversions and leads to the following short
   109 ## program::
   110 
   111 def sign_counting_inversions(p):
   112     """Compute sign of permutation `p` by counting the number of
   113     interchanges required to change the given permutation into the
   114     identity one.
   115 
   116     Examples::
   117       >>> sign_counting_inversions([0,1,2])
   118       1
   119       >>> sign_counting_inversions([0,2,1])
   120       -1
   121     """
   122     n = len(p)
   123     v = 0
   124     # count inversions
   125     for i in xrange(1,n):
   126         for j in xrange(i+1, n):
   127             if (p[i] > p[j]):
   128                 v += 1
   129     if (v % 2 == 0):
   130         return +1
   131     else:
   132         return -1
   133 
   134 
   135 ## And one more: Those with a sufficient background in algebra will
   136 ## realize that a simple algorithm exists based on the Laplace
   137 ## expansion of the determinant of a permutation matrix. However,
   138 ## while rather similar to the algorithm above that counts inversions,
   139 ## it turned out to be a bit more complicated to program::
   140 
   141 def lpsign(p):
   142     """Compute sign of permutation `p` via Laplace expansion of the
   143     determinant of the permutation matrix.
   144 
   145     Examples::
   146       >>> lpsign([0,1,2])
   147       1
   148       >>> lpsign([0,2,1])
   149       -1
   150     """
   151     # Laplace expansion of determinant of permutation 
   152     # matrix with a 1 in in row p[j] of column j
   153     n = len(p)
   154     s = 1
   155     for j in xrange(1, n):
   156         q=p[j]
   157         # apply checkerboard sign 
   158         if (q % 2 == 0):
   159             s = -s
   160         for k in xrange(j+1, n):
   161             r = p[k]
   162             # adjust numbers of remaining rows beyond row p[j]
   163             # to account for having deleted row p[j]
   164             if (r > q):
   165                 p[k] = r-1
   166     return s
   167 
   168 
   169 def sign_recursive(p):
   170     """Recursively compute the sign of permutation `p`.
   171 
   172     A permutation is represented as a linear list: `p` maps
   173     `i` to `p[i]`.
   174 
   175     Examples:
   176     >>> sign_recursive([0,1,2])
   177     1
   178     >>> sign_recursive([0,2,1])
   179     -1
   180     """
   181     n = len(p)
   182     if 1 >= n:
   183         return 1
   184     # find highest-numbered element
   185     k = p.index(n-1)
   186     # remove highest-numbered element for recursion
   187     del p[k] 
   188     # recursively compute
   189     if 0 == ((n+k) % 2):
   190         s = -1
   191     else:
   192         s = 1
   193     return s * sign_recursive(p)
   194 
   195 
   196 ## So, based on simplicity, counting inversions appears to have the
   197 ## edge over all the other algorithms discussed above.  However, based
   198 ## on speed, efficiently returning elements to their homes, and
   199 ## counting orbits both seem to be quite fast.
   200 ##
   201 ## Indeed, running 100000 iterations of each function to compute the
   202 ## sign of a 5-element permutation, gives::
   203 ##   
   204 ## |  $ python2.4 profile.py ./permsign.py 
   205 ## |  ncalls  tottime  percall  cumtime  percall  function
   206 ## |  [...]                                       
   207 ## |  100002    1.950    0.000    2.480    0.000  jbsign
   208 ## |  100002    2.550    0.000    3.140    0.000  sign_counting_inversions
   209 ## |  100002    2.810    0.000    3.350    0.000  sign_lang_formula
   210 ## |  100002    2.820    0.000    3.510    0.000  sign_by_orbits
   211 ## |  100002    3.470    0.000    4.190    0.000  lpsign
   212 ## |  100002   11.450    0.000   16.890    0.000  sign_recursive
   213 ## |  [...]
   214 ##
   215 
   216 
   217 
   218 ## main: run tests
   219 
   220 if "__main__" == __name__:
   221     import doctest
   222     doctest.testmod(optionflags=doctest.NORMALIZE_WHITESPACE)
   

   Posted by rmurri to python python algorithms mathematics permutations ... saved by 1 person ... 0 comments ... 1 year, 1 month

Pages : 1