Commit | Line | Data |
---|---|---|
d3f33b9a MW |
1 | ;;; -*-lisp-*- |
2 | ||
3 | ;;; This file isn't a program as such: rather, it's a collection of handy | |
4 | ;;; functions which can be used in an interactive session. | |
5 | ||
6 | ;;;-------------------------------------------------------------------------- | |
7 | ;;; General permutation utilities. | |
8 | ||
9 | (defun shuffle (v) | |
10 | "Randomly permute the elements of the vector V. Return V." | |
11 | (let ((n (length v))) | |
12 | (do ((k n (1- k))) | |
13 | ((<= k 1) v) | |
14 | (let ((i (random k))) | |
15 | (unless (= i (1- k)) | |
16 | (rotatef (aref v i) (aref v (1- k)))))))) | |
17 | ||
18 | (defun identity-permutation (n) | |
19 | "Return the do-nothing permutation on N elements." | |
20 | (let ((v (make-array n :element-type 'fixnum))) | |
21 | (dotimes (i n v) (setf (aref v i) i)))) | |
22 | ||
23 | (defun invert-permutation (p) | |
24 | "Given a permutation P, return its inverse." | |
25 | (let* ((n (length p)) (p-inv (make-array n :element-type 'fixnum))) | |
26 | (dotimes (i n) (setf (aref p-inv (aref p i)) i)) | |
27 | p-inv)) | |
28 | ||
29 | (defun next-permutation (v) | |
30 | "Adjust V so that it reflects the next permutation in ascending order. | |
31 | ||
32 | V should be a vector of real numbers. Returns V if successful, or nil if | |
33 | there are no more permutations." | |
34 | ||
35 | ;; The tail of the vector consists of a sequence ... A, Z, Y, X, ..., where | |
36 | ;; Z > Y > X ... is in reverse order, and A < Z. The next permutation is | |
37 | ;; then the smallest out of Z, Y, X, ... which is larger than A, followed | |
38 | ;; by the remaining elements in ascending order. | |
39 | ;; | |
40 | ;; Equivalently, reverse the tail Z, Y, X, ... so we have A, ... X, Y, Z, | |
41 | ;; and swap A with the next larger element. | |
42 | ||
43 | (let ((n (length v))) | |
44 | (cond ((< n 2) nil) | |
45 | (t (let* ((k (1- n)) | |
46 | (x (aref v k))) | |
47 | (loop (when (zerop k) (return-from next-permutation nil)) | |
48 | (decf k) | |
49 | (let ((y (aref v k))) | |
50 | (when (prog1 (< y x) | |
51 | (setf x y)) | |
52 | (return)))) | |
53 | (do ((i (1+ k) (1+ i)) | |
54 | (j (1- n) (1- j))) | |
55 | ((> i j)) | |
56 | (rotatef (aref v i) (aref v j))) | |
57 | (do ((i (- n 2) (1- i))) | |
58 | ((or (<= i k) (< (aref v i) x)) | |
59 | (rotatef (aref v k) (aref v (1+ i))))) | |
60 | v))))) | |
61 | ||
62 | (defun make-index-mask (w mask-expr) | |
63 | "Construct a bitmask based on bitwise properties of the bit indices. | |
64 | ||
65 | The function returns a W-bit mask in which each bit is set if MASK-EXPR | |
66 | of true of the bit's index. MASK-EXPR may be one of the following: | |
67 | ||
68 | * I -- an integer I is true if bit I of the bit index is set; | |
69 | * (not EXPR) -- is true if EXPR is false; | |
70 | * (and EXPR EXPR ...) -- is true if all of the EXPRs are true; and | |
71 | * (or EXPR EXPR ...) -- is true if any of the EXPRs is true." | |
72 | ||
73 | (let ((max-bit (1- (integer-length (1- w)))) | |
74 | (mask 0)) | |
75 | (dotimes (i w mask) | |
76 | (labels ((interpret (expr) | |
77 | (cond ((and (integerp expr) (<= 0 expr max-bit)) | |
78 | (logbitp expr i)) | |
79 | ((and (consp expr) (eq (car expr) 'not) | |
80 | (null (cddr expr))) | |
81 | (not (interpret (cadr expr)))) | |
82 | ((and (consp expr) (eq (car expr) 'and)) | |
83 | (every #'interpret (cdr expr))) | |
84 | ((and (consp expr) (eq (car expr) 'or)) | |
85 | (some #'interpret (cdr expr))) | |
86 | (t | |
87 | (error "unknown mask expression ~S" expr))))) | |
88 | (when (interpret mask-expr) | |
89 | (setf (ldb (byte 1 i) mask) 1)))))) | |
90 | ||
91 | (defun make-permutation-network (w steps) | |
92 | "Construct a permutation network. | |
93 | ||
94 | The integer W gives the number of bits to be acted upon. The STEPS are a | |
95 | list of instructions of the following forms: | |
96 | ||
97 | * (SHIFT . MASK) -- a pair of integers is treated literally; | |
98 | ||
99 | * (SHIFT MASK-EXPR) -- the SHIFT is literal, but the MASK-EXPR is | |
100 | processed by `make-index-mask' to calculate the mask; | |
101 | ||
102 | * (:invert I) -- make an instruction which inverts the sense of the | |
103 | index bit I; | |
104 | ||
105 | * (:exchange I J) -- make an instruction which exchanges index bits I | |
106 | and J; or | |
107 | ||
108 | * (:exchange-invert I J) -- make an instruction which exchanges and | |
109 | inverts index bits I and J. | |
110 | ||
111 | The output is a list of primitive (SHIFT . MASK) steps, indicating that | |
112 | the bits of the input selected by MASK are to be swapped with the bits | |
113 | selected by (ash MASK SHIFT)." | |
114 | ||
115 | (let ((max-mask (1- (ash 1 w))) | |
116 | (max-shift (1- w)) | |
117 | (max-bit (1- (integer-length (1- w)))) | |
118 | (list nil)) | |
119 | (dolist (step steps) | |
120 | (cond ((and (consp step) | |
121 | (integerp (car step)) (<= 0 (car step) max-shift) | |
122 | (integerp (cdr step)) (<= 0 (cdr step) max-mask)) | |
123 | (push step list)) | |
124 | ((and (consp step) | |
125 | (integerp (car step)) (<= 0 (car step) max-shift) | |
126 | (null (cddr step))) | |
127 | (push (cons (car step) (make-index-mask w (cadr step))) list)) | |
128 | ((and (consp step) | |
129 | (eq (car step) :invert) | |
130 | (integerp (cadr step)) (<= 0 (cadr step) max-bit) | |
131 | (null (cddr step))) | |
132 | (let ((i (cadr step))) | |
133 | (push (cons (ash 1 i) (make-index-mask w `(not ,i))) list))) | |
134 | ((and (consp step) | |
135 | (eq (car step) :exchange) | |
136 | (integerp (cadr step)) (integerp (caddr step)) | |
137 | (<= 0 (cadr step) (caddr step) max-bit) | |
138 | (null (cdddr step))) | |
139 | (let ((i (cadr step)) (j (caddr step))) | |
140 | (push (cons (- (ash 1 j) (ash 1 i)) | |
141 | (make-index-mask w `(and ,i (not ,j)))) | |
142 | list))) | |
143 | ((and (consp step) | |
144 | (eq (car step) :exchange-invert) | |
145 | (integerp (cadr step)) (integerp (caddr step)) | |
146 | (<= 0 (cadr step) (caddr step) max-bit) | |
147 | (null (cdddr step))) | |
148 | (let ((i (cadr step)) (j (caddr step))) | |
149 | (push (cons (+ (ash 1 i) (ash 1 j)) | |
150 | (make-index-mask w `(and (not ,i) (not ,j)))) | |
151 | list))) | |
152 | (t | |
153 | (error "unknown permutation step ~S" step)))) | |
154 | (nreverse list))) | |
155 | ||
156 | ;;;-------------------------------------------------------------------------- | |
157 | ;;; Permutation network diagnostics. | |
158 | ||
159 | (defun print-permutation-network (steps &optional (stream *standard-output*)) | |
160 | "Print a description of the permutation network STEPS to STREAM. | |
161 | ||
162 | A permutation network consists of a list of pairs | |
163 | ||
164 | (SHIFT . MASK) | |
165 | ||
166 | indicating that the bits selected by MASK, and those SHIFT bits to the | |
167 | left, should be exchanged. | |
168 | ||
169 | The output is intended to be human-readable and is subject to change." | |
170 | ||
171 | (let ((shiftwd 1) (maskwd 2)) | |
172 | ||
173 | ;; Determine suitable print widths for shifts and masks. | |
174 | (dolist (step steps) | |
175 | (let ((shift (car step)) (mask (cdr step))) | |
176 | (let ((swd (1+ (floor (log shift 10)))) | |
177 | (mwd (ash 1 (- (integer-length (1- (integer-length mask))) | |
178 | 2)))) | |
179 | (when (> swd shiftwd) (setf shiftwd swd)) | |
180 | (when (> mwd maskwd) (setf maskwd mwd))))) | |
181 | ||
182 | ;; Print the display. | |
183 | (pprint-logical-block (stream steps :prefix "(" :suffix ")") | |
184 | (let ((first t)) | |
185 | (dolist (step steps) | |
186 | (let ((shift (car step)) (mask (cdr step))) | |
187 | ||
188 | ;; Separate entries with newlines. | |
189 | (cond (first (setf first nil)) | |
190 | (t (pprint-newline :mandatory stream))) | |
191 | ||
192 | (let ((swaps nil)) | |
193 | ||
194 | ;; Determine the list of exchanges implied by the mask. | |
195 | (dotimes (i (integer-length mask)) | |
196 | (when (logbitp i mask) | |
197 | (push (cons i (+ i shift)) swaps))) | |
198 | (setf swaps (nreverse swaps)) | |
199 | ||
200 | ;; Print the entry. | |
201 | (format stream "~@<(~;~vD #x~(~v,'0X~) ~8I~:@_~W~;)~:>" | |
202 | shiftwd shift maskwd mask swaps)))))) | |
203 | ||
204 | ;; Print a final newline following the close parenthesis. | |
205 | (terpri stream))) | |
206 | ||
207 | (defun demonstrate-permutation-network | |
7306ec27 MW |
208 | (n steps |
209 | &key reference | |
210 | (stream *standard-output*)) | |
d3f33b9a MW |
211 | "Print, on STREAM, a demonstration of the permutation STEPS. |
212 | ||
45be3aa8 MW |
213 | The output is intended to be useful to human readers and is subject to |
214 | change. Currently, it prints a sequence of diagrams on STREAM. The left | |
215 | hand side of each row shows a map of which bits are affected: `-' means | |
216 | that the bit remains in the same position, `*' means that it moves | |
217 | forward, and `#' means that it moves back; each `*' pairs with the | |
218 | earliest unpaired `#' marker. The right hand side shows the arrangement | |
219 | of the original input bits. | |
220 | ||
221 | If REFERENCE is not nil, then print a final pair of diagrams. This time, | |
222 | the map shows `-' for correct bits and `x' for incorrect ones, with the | |
223 | right hand side showing the expected arrangement of input bits. | |
224 | ||
225 | The function returns non-nil if the STEPS resulted in the REFERENCE | |
226 | permutation, and nil if either the STEPS are incorrect or no REFERENCE was | |
227 | provided." | |
228 | ||
229 | (flet ((apply-step (shift mask v) | |
230 | (dotimes (k n) | |
231 | (when (logbitp k mask) | |
232 | (rotatef (aref v k) (aref v (+ k shift))))))) | |
233 | ||
234 | (let* ((v (identity-permutation n)) | |
235 | (end (or reference | |
236 | (let ((e (identity-permutation n))) | |
237 | (dolist (step steps e) | |
238 | (let ((shift (car step)) (mask (cdr step))) | |
239 | (apply-step shift mask e)))))) | |
240 | (end-inv (invert-permutation end)) | |
241 | (mapwd (ceiling (sqrt n))) | |
242 | (ixwd (length (format nil "~D" (1- n))))) | |
243 | ||
244 | (flet ((show-stage (shift mask v) | |
245 | (do ((i 0 (+ i mapwd))) | |
246 | ((>= i n)) | |
247 | (write-string ";; " stream) | |
248 | (dotimes (j mapwd) | |
249 | (let ((k (+ i j))) | |
250 | (when (plusp j) (write-char #\space stream)) | |
251 | (write-char (cond ((>= k n) | |
252 | #\space) | |
253 | ((logbitp k mask) | |
254 | #\*) | |
255 | ((and (>= k shift) | |
256 | (logbitp (- k shift) mask)) | |
257 | #\#) | |
258 | (t | |
259 | #\-)) | |
260 | stream))) | |
261 | (write-string " | " stream) | |
262 | ||
263 | (fresh-line) | |
264 | (show-stage 0 0 v) | |
d3f33b9a | 265 | |
45be3aa8 MW |
266 | (dolist (step steps) |
267 | (let ((shift (car step)) (mask (cdr step))) | |
268 | (apply-step shift mask v) | |
269 | (format stream ";;~%") | |
270 | (show-stage shift mask v))) | |
d3f33b9a | 271 | |
45be3aa8 | 272 | (let ((ok (not (null reference)))) |
d3f33b9a | 273 | (when reference |
45be3aa8 MW |
274 | (format stream ";;~%") |
275 | (do ((i 0 (+ i mapwd))) | |
276 | ((>= i n)) | |
277 | (write-string ";; " stream) | |
278 | (dotimes (j mapwd) | |
279 | (let ((k (+ i j))) | |
280 | (when (plusp j) (write-char #\space stream)) | |
281 | (write-char (cond ((>= k n) #\space) | |
282 | ((/= (aref v k) (aref reference k)) #\x) | |
283 | (t #\-)) | |
284 | stream))) | |
285 | (write-string " | " stream) | |
286 | (dotimes (j (min mapwd (- n i))) | |
287 | (let ((k (+ i j))) | |
288 | (when (plusp j) (write-char #\space stream)) | |
289 | (format stream "~vD" ixwd (aref reference k)))) | |
290 | (terpri)) | |
291 | (unless (every #'= v reference) | |
292 | (setf ok nil)) | |
293 | (format stream "~:[FAIL~;pass~]~%" ok)) | |
294 | ok))))) | |
d3f33b9a MW |
295 | |
296 | ;;;-------------------------------------------------------------------------- | |
c7c44436 MW |
297 | ;;; Beneš networks. |
298 | ||
299 | (defun compute-benes-step (n p p-inv bit clear-input) | |
300 | "Compute a single layer of a Beneš network. | |
301 | ||
302 | N is a fixnum. P is a vector of fixnums defining a permutation: for each | |
303 | output bit position i (numbering the least significant bit 0), element i | |
304 | gives the number of the input which should end up in that position; P-INV | |
305 | gives the inverse permutation in the same form. BIT is a power of 2 which | |
306 | gives the distance between bits we should consider. CLEAR-INPUT is | |
307 | a (generalized) boolean: if true, we attempt to do no work on the input | |
308 | side; if false, we try to do no work on the output side. The length of P | |
309 | must be at least (logior N BIT). | |
310 | ||
311 | The output consists of a pair of masks M0 and M1, to be used on the input | |
312 | and output sides respectively. The permutation stage, applied to an input | |
313 | X, is as follows: | |
314 | ||
315 | (let ((tmp (logand (logxor x (ash x (- bit))) mask))) | |
316 | (logxor x tmp (ash tmp bit))) | |
317 | ||
318 | The critical property of the masks is that it's possible to compute an | |
319 | inner permutation, mapping the output of the M0 stage to the input of the | |
320 | M1 stage, such that (a) the overall composition of the three permutations | |
321 | is the given permutation P, and (b) the distances that the bits need to | |
322 | be moved by the inner permutation all have BIT clear. | |
323 | ||
324 | The resulting permutation will attempt to avoid touching elements at | |
325 | indices greater than N. This attempt will succeed if all such elements | |
326 | contain values equal to their indices. | |
327 | ||
328 | By appropriately composing layers computed by this function, then, it's | |
329 | possible to perform an arbitrary permutation of 2^n bits in 2 n - 1 simple | |
330 | steps like the one above." | |
331 | ||
332 | ;; Consider the following problem. You're given two equally-sized | |
333 | ;; collections of things, called `left' and `right'. Each left thing is | |
334 | ;; attached to exactly one right thing with a string, and /vice versa/. | |
335 | ;; Furthermore, the left things, and the right things, are each linked | |
336 | ;; together in pairs, so each pair has two strings coming out of it. Our | |
337 | ;; job is to paint the strings so that each linked pair of things has one | |
338 | ;; red string and one blue string. | |
339 | ;; | |
340 | ;; This is quite straightforward. Pick a pair whose strings aren't yet | |
341 | ;; coloured, and colour one of its strings, chosen arbitrarily, red. Find | |
342 | ;; the pair at the other end of this red string. If the two other things | |
343 | ;; in these two pairs are connected, then paint that string blue and move | |
344 | ;; on. Otherwise, both things have an uncoloured string, so paint both of | |
345 | ;; them blue and trace along these now blue strings to find two more thing | |
346 | ;; pairs. Again, the other thing in each pair has an uncoloured string. | |
347 | ;; If these things share the /same/ string, paint it red and move on. | |
348 | ;; Otherwise, paint both strings red, trace, and repeat. Eventually, we'll | |
349 | ;; find two things joined by the same string, each paired with another | |
350 | ;; thing whose strings we've just painted the same colour. Once we're | |
351 | ;; done, we'll have painted a cycle's worth of strings, and each pair of | |
352 | ;; things will have either both of its strings painted different colours, | |
353 | ;; or neither of them. | |
354 | ;; | |
355 | ;; The right things are the integers 0, 1, ..., n - 1, where n is some | |
356 | ;; power of 2, paired according to whether they differ only in BIT. The | |
357 | ;; left things are the same integers, connected to the right things | |
358 | ;; according to the permutation P: the right thing labelled i is connected | |
359 | ;; to the left thing labelled P(i). Similarly, two left things are paired | |
360 | ;; if their labels P(i) and P(j) differ only in BIT. We're going to paint | |
361 | ;; a string red if we're going to arrange to clear BIT in the labels at | |
362 | ;; both ends, possibly by swapping the two labels, and paint it red if | |
363 | ;; we're going to arrange to set BIT. Once we've done this, later stages | |
364 | ;; of the filter will permute the red- and blue-painted things | |
365 | ;; independently. | |
366 | ||
367 | (let ((m0 0) (m1 0) (done 0)) | |
368 | ||
369 | ;; Now work through the permutation cycles. | |
370 | (do ((i (1- n) (1- i))) | |
371 | ((minusp i)) | |
372 | ||
373 | ;; Skip if we've done this one already. | |
374 | (unless (or (plusp (logand i bit)) | |
375 | (logbitp i done)) | |
376 | ||
377 | ;; Find the other associated values. | |
378 | (let* ((i0 i) (i1 (aref p-inv i)) | |
379 | (sense (cond ((>= (logior i0 bit) n) 0) | |
380 | (clear-input (logand i0 bit)) | |
381 | (t (logand i1 bit))))) | |
382 | ||
383 | #+noise | |
384 | (format t ";; new cycle: i0 = ~D, j0 = ~D; i1 = ~D, j1 = ~D~%" | |
385 | i0 (logxor i0 bit) | |
386 | i1 (logxor i1 bit)) | |
387 | ||
388 | ;; Mark this index as done. | |
389 | (setf (ldb (byte 1 i0) done) 1) | |
390 | #+noise (format t ";; done = #x~2,'0X~%" done) | |
391 | ||
392 | ;; Now trace round the cycle. | |
393 | (loop | |
394 | ||
395 | ;; Mark this index as done. | |
396 | (setf (ldb (byte 1 (logandc2 i0 bit)) done) 1) | |
397 | #+noise (format t ";; done = #x~2,'0X~%" done) | |
398 | ||
399 | ;; Swap the input and output pairs if necessary. | |
400 | (unless (= (logand i0 bit) sense) | |
401 | #+noise | |
402 | (format t ";; swap input: ~D <-> ~D~%" | |
403 | (logandc2 i0 bit) (logior i0 bit)) | |
404 | (setf (ldb (byte 1 (logandc2 i0 bit)) m0) 1)) | |
405 | (unless (= (logand i1 bit) sense) | |
406 | #+noise | |
407 | (format t ";; swap output: ~D <-> ~D~%" | |
408 | (logandc2 i1 bit) (logior i1 bit)) | |
409 | (setf (ldb (byte 1 (logandc2 i1 bit)) m1) 1)) | |
410 | ||
411 | ;; Advance around the cycle. | |
412 | (let* ((j0 (logxor i0 bit)) | |
413 | (j1 (logxor i1 bit)) | |
414 | (next-i1 (aref p-inv j0)) | |
415 | (next-i0 (aref p j1))) | |
416 | (when (= next-i0 j0) (return)) | |
417 | (setf i0 next-i0 | |
418 | i1 next-i1 | |
419 | sense (logxor sense bit))) | |
420 | ||
421 | #+noise | |
422 | (format t ";; advance: i0 = ~D, j0 = ~D; i1 = ~D, j1 = ~D~%" | |
423 | i0 (logxor i0 bit) | |
424 | i1 (logxor i1 bit)))))) | |
425 | ||
426 | (values m0 m1))) | |
427 | ||
428 | (defun compute-final-benes-step (n p p-inv bit) | |
429 | "Determine the innermost stage of a Beneš network. | |
430 | ||
431 | N is a fixnum. P is a vector of fixnums defining a permutation: for each | |
432 | output bit position i (numbering the least significant bit 0), element i | |
433 | gives the number of the input which should end up in that position; P-INV | |
434 | gives the inverse permutation in the same form. BIT is a power of 2 which | |
435 | gives the distance between bits we should consider. The length of P must | |
436 | be at least (logior N BIT). | |
437 | ||
438 | Furthermore, the ith element of P must be equal either to i or to | |
439 | (logxor i BIT); and therefore P-INV must be equal to P. | |
440 | ||
441 | Return the mask such that | |
442 | ||
443 | (let ((tmp (logand (logxor x (ash x (- bit))) mask))) | |
444 | (logxor x tmp (ash tmp bit))) | |
445 | ||
446 | applies the permutation P to the bits of x." | |
447 | ||
448 | (declare (ignorable p-inv)) | |
449 | ||
450 | (let ((m 0)) | |
451 | (dotimes (i n) | |
452 | (unless (plusp (logand i bit)) | |
453 | (let ((x (aref p i))) | |
454 | #+paranoid | |
455 | (assert (= (logandc2 x bit) i)) | |
456 | #+paranoid | |
457 | (assert (= x (aref p-inv i))) | |
458 | ||
459 | (unless (= x i) | |
460 | (setf (ldb (byte 1 i) m) 1))))) | |
461 | m)) | |
462 | ||
463 | (defun apply-benes-step (p p-inv bit m0 m1) | |
464 | "Apply input and output steps for a Beneš network to a permutation. | |
465 | ||
466 | Given the permutation P and its inverse, and the distance BIT, as passed | |
467 | to `compute-benes-step', and the masks M0 and M1 returned, determine and | |
468 | return the necessary `inner' permutation to be applied between these | |
469 | steps, and its inverse. | |
470 | ||
471 | A permutation-network step, and, in particular, a Beneš step, is an | |
472 | involution, so the change to the vectors P and P-INV can be undone by | |
473 | calling the function again with the same arguments." | |
474 | ||
475 | (flet ((swaps (p p-inv mask) | |
476 | (dotimes (i0 (length p)) | |
477 | (when (logbitp i0 mask) | |
478 | (let* ((j0 (logior i0 bit)) | |
479 | (i1 (aref p-inv i0)) | |
480 | (j1 (aref p-inv j0))) | |
481 | (setf (aref p i1) j0 | |
482 | (aref p j1) i0) | |
483 | (rotatef (aref p-inv i0) (aref p-inv j0))))))) | |
484 | (swaps p p-inv m0) | |
485 | (swaps p-inv p m1) | |
486 | ||
487 | #+paranoid | |
488 | (let* ((n (length p))) | |
489 | (dotimes (i n) | |
490 | (assert (= (aref p (aref p-inv i)) i)) | |
491 | (assert (= (aref p-inv (aref p i)) i)))))) | |
492 | ||
493 | (defun benes-search (p) | |
494 | "Given a bit permutation P, describe a Beneš network implementing P. | |
495 | ||
496 | P is a sequence of fixnums defining a permutation: for each output bit | |
497 | position i (numbering the least significant bit 0), element i gives the | |
498 | number of the input which should end up in that position. | |
499 | ||
500 | The return value is a list of steps of the form | |
501 | ||
502 | (BIT MASK (X . Y) (X' . Y') ...) | |
503 | ||
504 | To implement this permutation step: | |
505 | ||
506 | * given an input X, compute | |
507 | ||
508 | (let ((tmp (logand (logxor x (ash x (- bit))) mask))) | |
509 | (logxor x tmp (ash tmp bit))) | |
510 | ||
511 | or, equivalently, | |
512 | ||
513 | * exchange the bits in the positions given in each of the pairs X, Y, | |
514 | ..., where each Y = X + BIT." | |
515 | ||
516 | (let* ((n (length p)) | |
517 | (w (ash 1 (integer-length (1- n)))) | |
518 | (p (let ((new (make-array w :element-type 'fixnum))) | |
519 | (replace new p) | |
520 | (do ((i n (1+ i))) | |
521 | ((>= i w)) | |
522 | (setf (aref new i) i)) | |
523 | new)) | |
524 | (p-inv (invert-permutation p))) | |
525 | ||
526 | (labels ((recurse (todo) | |
527 | ;; Main recursive search. DONE is a mask of the bits which | |
528 | ;; have been searched. Return the number of skipped stages | |
529 | ;; and a list of steps (BIT M0 M1), indicating that (BIT M0) | |
530 | ;; should be performed before the following stages, and | |
531 | ;; (BIT M1) should be performed afterwards. | |
532 | ;; | |
533 | ;; The permutation `p' and its inverse `p-inv' will be | |
534 | ;; modified and restored. | |
535 | ||
536 | (cond ((zerop (logand todo (1- todo))) | |
537 | ;; Only one more bit left. Use the more efficient | |
538 | ;; final-step computation. | |
539 | ||
540 | (let ((m (compute-final-benes-step n p p-inv todo))) | |
541 | (values (if m 0 1) (list (list todo m 0))))) | |
542 | ||
543 | (t | |
544 | ;; More searching to go. We'll keep the result which | |
545 | ;; maximizes the number of skipped stages. | |
546 | (let ((best-list nil) | |
547 | (best-skips -1)) | |
548 | ||
549 | (flet ((try (bit clear-input) | |
550 | ;; Try a permutation with the given BIT and | |
551 | ;; CLEAR-INPUT arguments to | |
552 | ;; `compute-benes-step'. | |
553 | ||
554 | ;; Compute the next step. | |
555 | (multiple-value-bind (m0 m1) | |
556 | (compute-benes-step n p p-inv | |
557 | bit clear-input) | |
558 | ||
559 | ;; Apply the step and recursively | |
560 | ;; determine the inner permutation. | |
561 | (apply-benes-step p p-inv bit m0 m1) | |
562 | (multiple-value-bind (nskip tail) | |
563 | (recurse (logandc2 todo bit)) | |
564 | (apply-benes-step p p-inv bit m0 m1) | |
565 | ||
566 | ;; Work out how good this network is. | |
567 | ;; Keep it if it improves over the | |
568 | ;; previous attempt. | |
569 | (when (zerop m0) (incf nskip)) | |
570 | (when (zerop m1) (incf nskip)) | |
571 | (when (> nskip best-skips) | |
572 | (setf best-list | |
573 | (cons (list bit m0 m1) | |
574 | tail) | |
575 | best-skips | |
576 | nskip)))))) | |
577 | ||
578 | ;; Work through each bit that we haven't done | |
579 | ;; already, and try skipping both the start and end | |
580 | ;; steps. | |
581 | (do ((bit 1 (ash bit 1))) | |
582 | ((>= bit w)) | |
583 | (when (plusp (logand bit todo)) | |
584 | (try bit nil) | |
585 | (try bit t)))) | |
586 | (values best-skips best-list)))))) | |
587 | ||
588 | ;; Find the best permutation network. | |
589 | (multiple-value-bind (nskip list) (recurse (1- w)) | |
590 | (declare (ignore nskip)) | |
591 | ||
592 | ;; Turn the list returned by `recurse' into a list of (SHIFT MASK) | |
593 | ;; entries as expected by everything else. | |
594 | (let ((head nil) (tail nil)) | |
595 | (dolist (step list (nconc (nreverse head) tail)) | |
596 | (destructuring-bind (bit m0 m1) step | |
597 | (when (plusp m0) (push (cons bit m0) head)) | |
598 | (when (plusp m1) (push (cons bit m1) tail))))))))) | |
599 | ||
600 | ;;;-------------------------------------------------------------------------- | |
601 | ;;; Special functions for DES permutations. | |
602 | ||
603 | (defun benes-search-des (p &optional attempts) | |
604 | "Search for a Beneš network for a DES 64-bit permutation. | |
605 | ||
606 | P must be a sequence of 64 fixnums, each of which is between 0 and 64 | |
607 | inclusive. In the DES convention, bits are numbered with the most- | |
608 | significant bit being bit 1, and increasing towards the least-significant | |
609 | bit, which is bit 64. Each nonzero number must appear at most once, and | |
610 | specifies which input bit must appear in that output position. There may | |
611 | also be any number of zero entries, which mean `don't care'. | |
612 | ||
613 | This function searches for and returns a Beneš network which implements a | |
614 | satisfactory permutation. If ATTEMPTS is nil or omitted, then search | |
615 | exhaustively, returning the shortest network. Otherwise, return the | |
616 | shortest network found after considering ATTEMPTS randomly chosen | |
617 | matching permutations." | |
618 | ||
619 | (let* ((n (length p)) | |
620 | (p (map '(vector fixnum) | |
621 | (lambda (x) | |
622 | (if (zerop x) -1 | |
623 | (- 64 x))) | |
624 | (reverse p))) | |
625 | (seen (make-hash-table)) | |
626 | (nmissing 0) (missing nil) (indices nil)) | |
627 | ||
628 | ;; Find all of the `don't care' slots, and keep track of the bits which | |
629 | ;; have homes to go to. | |
630 | (dotimes (i n) | |
631 | (let ((x (aref p i))) | |
632 | (cond ((minusp x) | |
633 | (push i indices) | |
634 | (incf nmissing)) | |
635 | (t (setf (gethash x seen) t))))) | |
636 | ||
637 | ;; Fill in numbers of the input bits which don't have fixed places to go. | |
638 | (setf missing (make-array nmissing :element-type 'fixnum)) | |
639 | (let ((j 0)) | |
640 | (dotimes (i n) | |
641 | (unless (gethash i seen) | |
642 | (setf (aref missing j) i) | |
643 | (incf j))) | |
644 | (assert (= j nmissing))) | |
645 | ||
646 | ;; Run the search, printing successes as we find them to keep the user | |
647 | ;; amused. | |
648 | (let ((best nil) (best-length nil)) | |
649 | (loop | |
650 | (cond ((eql attempts 0) (return best)) | |
651 | (attempts (shuffle missing) (decf attempts)) | |
652 | ((null (next-permutation missing)) (return best))) | |
653 | (do ((idx indices (cdr idx)) | |
654 | (i 0 (1+ i))) | |
655 | ((endp idx)) | |
656 | (setf (aref p (car idx)) (aref missing i))) | |
657 | (let* ((benes (benes-search p)) (len (length benes))) | |
658 | (when (or (null best-length) | |
659 | (< len best-length)) | |
660 | (setf best-length len | |
661 | best benes) | |
662 | (print-permutation-network benes) | |
663 | (force-output))))))) | |
664 | ||
665 | ;;;-------------------------------------------------------------------------- | |
d3f33b9a MW |
666 | ;;; Examples and useful runes. |
667 | ||
668 | #+example | |
669 | (let* ((ip #(58 50 42 34 26 18 10 2 | |
670 | 60 52 44 36 28 20 12 4 | |
671 | 62 54 46 38 30 22 14 6 | |
672 | 64 56 48 40 32 24 16 8 | |
673 | 57 49 41 33 25 17 9 1 | |
674 | 59 51 43 35 27 19 11 3 | |
675 | 61 53 45 37 29 21 13 5 | |
676 | 63 55 47 39 31 23 15 7)) | |
677 | (fixed-ip (map '(vector fixnum) | |
678 | (lambda (x) (- 64 x)) | |
679 | (reverse ip))) | |
70f0901a | 680 | |
45be3aa8 | 681 | ;; The traditional network. |
d3f33b9a MW |
682 | (trad-network |
683 | (make-permutation-network | |
684 | 64 ; 5 4 3 2 1 0 | |
685 | '((:exchange-invert 2 5) ; ~2 4 3 ~5 1 0 | |
686 | (:exchange-invert 1 4) ; ~2 ~1 3 ~5 ~4 0 | |
687 | (:exchange-invert 0 3) ; ~2 ~1 ~0 ~5 ~4 ~3 | |
688 | (:exchange-invert 3 4) ; ~2 0 1 ~5 ~4 ~3 | |
48af823d | 689 | (:exchange-invert 4 5)))) ; ~0 2 1 ~5 ~4 ~3 |
70f0901a | 690 | |
45be3aa8 | 691 | ;; The new twizzle-optimized network. |
48af823d MW |
692 | (new-network |
693 | (make-permutation-network | |
694 | 64 ; 5 4 3 2 1 0 | |
695 | '((:exchange-invert 2 5) ; ~2 4 3 ~5 1 0 | |
696 | (:exchange-invert 4 5) ; ~4 2 3 ~5 1 0 | |
697 | (:exchange 1 5) ; 1 2 3 ~5 ~4 0 | |
698 | (:exchange 3 5) ; 3 2 1 ~5 ~4 0 | |
699 | (:exchange-invert 0 5))))) ; ~0 2 1 ~5 ~4 ~3 | |
d3f33b9a MW |
700 | |
701 | (fresh-line) | |
702 | ||
c7c44436 MW |
703 | (let ((benes-network (benes-search fixed-ip))) |
704 | (print-permutation-network benes-network) | |
7306ec27 | 705 | (demonstrate-permutation-network 64 benes-network :reference fixed-ip)) |
c7c44436 | 706 | (terpri) |
d3f33b9a | 707 | (print-permutation-network trad-network) |
7306ec27 | 708 | (demonstrate-permutation-network 64 trad-network :reference fixed-ip) |
48af823d MW |
709 | (terpri) |
710 | (print-permutation-network new-network) | |
7306ec27 | 711 | (demonstrate-permutation-network 64 new-network :reference fixed-ip)) |
c7c44436 MW |
712 | |
713 | #+example | |
714 | (benes-search-des #( 0 0 0 0 | |
715 | 57 49 41 33 25 17 9 1 | |
716 | 58 50 42 34 26 18 10 2 | |
717 | 59 51 43 35 27 19 11 3 | |
718 | 60 52 44 36 | |
719 | 0 0 0 0 | |
720 | 63 55 47 39 31 23 15 7 | |
721 | 62 54 46 38 30 22 14 6 | |
722 | 61 53 45 37 29 21 13 5 | |
723 | 28 20 12 4)) | |
724 | ||
725 | #+example | |
726 | (let ((pc2 (make-array '(8 6) | |
727 | :element-type 'fixnum | |
728 | :initial-contents '((14 17 11 24 1 5) | |
729 | ( 3 28 15 6 21 10) | |
730 | (23 19 12 4 26 8) | |
731 | (16 7 27 20 13 2) | |
732 | (41 52 31 37 47 55) | |
733 | (30 40 51 45 33 48) | |
734 | (44 49 39 56 34 53) | |
735 | (46 42 50 36 29 32))))) | |
736 | (benes-search-des | |
737 | (make-array 64 | |
738 | :element-type 'fixnum | |
739 | :initial-contents | |
740 | (loop for i in '(2 4 6 8 1 3 5 7) | |
741 | nconc (list 0 0) | |
742 | nconc (loop for j below 6 | |
743 | for x = (aref pc2 (1- i) j) | |
744 | collect (if (<= x 32) (+ x 4) (+ x 8))))) | |
745 | 1000)) |