snippets / mathematics
All snippets tagged mathematics (3)
-
Generate all partitions of an integer
Iterative algorithms to generate all partitions of a given integer. Only tested with Python 2.5
1 import itertools
2
3 class FixedLengthPartitionIterator:
4 """Iterate over partitions of integer `N` into *exactly* `K`
5 positive integers.
6
7 Each returned partition is a list of positive integers in
8 descending order, such that their sum is `N`.
9
10 Arguments `min_` and `max_` bound the values in each partition.
11
12 Examples::
13 >>> list(FixedLengthPartitionIterator(3,1))
14 [(3,)]
15 >>> list(FixedLengthPartitionIterator(3,2))
16 [(2, 1)]
17 >>> list(FixedLengthPartitionIterator(3,3))
18 [(1, 1, 1)]
19 >>> list(FixedLengthPartitionIterator(6,2))
20 [(5, 1), (4, 2), (3, 3)]
21 >>> list(FixedLengthPartitionIterator(8,3))
22 [(6, 1, 1), (5, 2, 1), (4, 3, 1), (4, 2, 2), (3, 3, 2)]
23 >>> list(FixedLengthPartitionIterator(8,4,2))
24 [(2, 2, 2, 2)]
25 >>> list(FixedLengthPartitionIterator(8,3,2))
26 [(4, 2, 2), (3, 3, 2)]
27 >>> list(FixedLengthPartitionIterator(8,3,2,3))
28 [(3, 3, 2)]
29 """
30 def __init__(self, N, K, min_=1, max_=None):
31 # `max_` really defaults to N
32 if max_ is None:
33 max_ = N
34
35 # rule out trivial cases
36 if (K*min_ > N) or (K*max_ < N):
37 self.done = True
38 return
39
40 self._N = N #: integer to be partitioned
41 self._K = K #: maximum number of nonzero parts
42 self._k = 0 #: current number of nonzero parts
43 self._min = min_ #: minimum value of each part
44 self._max = max_ #: maximum value of each part
45 self._p = [0] * K #: current partition
46 self.done = False #: when `True`, enumeration is over
47
48 def __iter__(self):
49 return self
50
51 def next(self):
52 if self.done:
53 raise StopIteration
54 if self._N == self._K * self._min:
55 self.done = True
56 return tuple([self._min]*self._K)
57 else:
58 while self._k <= self._K:
59 i = self._k - 1
60 while i > 0:
61 if (self._p[i]+1 <= self._p[i-1]-1) \
62 and (self._p[i-1]-1 >= self._min) \
63 and (self._p[i]+1 <= self._max):
64 self._p[i-1] -= 1
65 self._p[i] += 1
66 return tuple(self._p)
67 else:
68 i -= 1
69 # only change the first `k` parts
70 self._k += 1
71 head = (self._N - self._k - self._min*self._K + self._min + 1)
72 if (head < self._min+1) or (head > self._max) \
73 or (self._k > self._K):
74 continue
75 # advance to next partition:
76 # [N-2*(k-1)-(K-k), 2, ..., 2, (k-1 times) 1, ..., 1 (K-k times)]
77 self._p = [head] \
78 + ([self._min + 1] * (self._k - 1)) \
79 + ([self._min] * (self._K - self._k))
80 return tuple(self._p)
81 raise StopIteration
82
83
84 def PartitionIterator(N, K, min_=1, max_=None):
85 """Iterate over partitions of integer `N` into *at most* `K`
86 positive integers.
87
88 Each returned partition is a list of positive integers in
89 descending order, such that their sum is `N`.
90
91 Optional arguments `min_` and `max_` bound the values in each
92 partition.
93
94 Examples::
95 >>> list(PartitionIterator(2,1))
96 [(2,)]
97 >>> list(PartitionIterator(3,3))
98 [(3,), (2, 1), (1, 1, 1)]
99 >>> list(PartitionIterator(8,3,2))
100 [(8,), (6, 2), (5, 3), (4, 4), (4, 2, 2), (3, 3, 2)]
101 >>> list(PartitionIterator(8,3,2,3))
102 [(3, 3, 2)]
103 """
104 return itertools.chain(*[FixedLengthPartitionIterator(N,k,min_,max_)
105 for k in xrange(1,K+1)]) -
A function to enumerate the items of a set product.
A function to enumerate the items of a set product. A recursive version is also provided for completeness, although it is much slower.
1 """A function to enumerate the items of a set product.
2 """
3 __docformat__ = 'reStructuredText'
4
5
6 def enumerate_set_product(*args):
7 """Iterate over all elements in the cartesian product of its arguments.
8 Each argument to this variadic function *must* be a sequence.
9
10 Examples::
11 >>> list(enumerate_set_product())
12 [[]]
13 >>> list(enumerate_set_product([1]))
14 [[1]]
15 >>> list(enumerate_set_product([1],[1]))
16 [[1, 1]]
17 >>> list(enumerate_set_product([1,2],[1]))
18 [[1, 1], [2, 1]]
19 >>> list(enumerate_set_product([1,2],[1,2]))
20 [[1, 1], [2, 1], [1, 2], [2, 2]]
21 """
22
23 # `assert` guard against invocation with wrong arguments.
24 def __all_items_are_sequences(s):
25 from operator import and_, isSequenceType
26 return reduce(and_,
27 [ isSequenceType(item) for item in s],
28 True)
29 assert __all_items_are_sequences(args), \
30 "All arguments to `enumerate_set_products` must be sequences."
31
32 if len(args) == 0:
33 yield []
34 else:
35 L = len(args)
36 M = [ len(s)-1 for s in args ]
37 m = [0] * L
38 i = 0
39 while i < L:
40 # return element corresponding to current multi-index
41 yield [ s[m[i]] for (i,s) in enumerate(args) ]
42 # advance multi-index
43 i = 0
44 while (i < L):
45 if m[i] == M[i]:
46 m[i] = 0
47 i += 1
48 else:
49 m[i] += 1
50 break
51
52
53 def recursively_enumerate_set_product(*args):
54 """Iterate over all elements in the cartesian product of its arguments.
55 Each argument to this variadic function *must* be a sequence.
56
57 You would probably prefer `enumerate_set_product`, which is faster.
58
59 Examples::
60 >>> list(recursively_enumerate_set_product())
61 [[]]
62 >>> list(recursively_enumerate_set_product([1]))
63 [[1]]
64 >>> list(recursively_enumerate_set_product([1],[1]))
65 [[1, 1]]
66 >>> list(recursively_enumerate_set_product([1,2],[1]))
67 [[1, 1], [2, 1]]
68 >>> list(recursively_enumerate_set_product([1,2],[1,2]))
69 [[1, 1], [2, 1], [1, 2], [2, 2]]
70 """
71
72 # `assert` guard against invocation with wrong arguments.
73 def __all_items_are_sequences(s):
74 from operator import and_, isSequenceType
75 return reduce(and_,
76 [ isSequenceType(item) for item in s],
77 True)
78 assert __all_items_are_sequences(args), \
79 "All arguments to `enumerate_set_products` must be sequences."
80
81 if len(args) == 0:
82 yield []
83 else:
84 for i in args[-1]:
85 for js in recursively_enumerate_set_product(*args[:-1]):
86 yield js+[i]
87
88
89 ## main: run tests
90
91 if "__main__" == __name__:
92 import doctest
93 doctest.testmod(optionflags=doctest.NORMALIZE_WHITESPACE) -
Various algorithms for computing the sign of a permutation.
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)
Pages : 1
Flux RSS friendsnippetLatest snippets
- c++ QT saatin yaninda... | 6 months, 2 weeks
- src/css/styles.css | 6 months, 2 weeks
- src/charts/Node.as | 6 months, 2 weeks
- src/charts/NodeSkin.mxml | 6 months, 2 weeks
- src/programmaticSkinClass... | 6 months, 2 weeks
