chiark / gitweb /
base/permute.h, utils/permute.lisp, symm/...: Formalize bit permutations.
[catacomb] / utils / permute.lisp
diff --git a/utils/permute.lisp b/utils/permute.lisp
new file mode 100644 (file)
index 0000000..372f5cf
--- /dev/null
@@ -0,0 +1,280 @@
+;;; -*-lisp-*-
+
+;;; This file isn't a program as such: rather, it's a collection of handy
+;;; functions which can be used in an interactive session.
+
+;;;--------------------------------------------------------------------------
+;;; General permutation utilities.
+
+(defun shuffle (v)
+  "Randomly permute the elements of the vector V.  Return V."
+  (let ((n (length v)))
+    (do ((k n (1- k)))
+       ((<= k 1) v)
+      (let ((i (random k)))
+       (unless (= i (1- k))
+         (rotatef (aref v i) (aref v (1- k))))))))
+
+(defun identity-permutation (n)
+  "Return the do-nothing permutation on N elements."
+  (let ((v (make-array n :element-type 'fixnum)))
+    (dotimes (i n v) (setf (aref v i) i))))
+
+(defun invert-permutation (p)
+  "Given a permutation P, return its inverse."
+  (let* ((n (length p)) (p-inv (make-array n :element-type 'fixnum)))
+    (dotimes (i n) (setf (aref p-inv (aref p i)) i))
+    p-inv))
+
+(defun next-permutation (v)
+  "Adjust V so that it reflects the next permutation in ascending order.
+
+   V should be a vector of real numbers.  Returns V if successful, or nil if
+   there are no more permutations."
+
+  ;; The tail of the vector consists of a sequence ... A, Z, Y, X, ..., where
+  ;; Z > Y > X ... is in reverse order, and A < Z.  The next permutation is
+  ;; then the smallest out of Z, Y, X, ... which is larger than A, followed
+  ;; by the remaining elements in ascending order.
+  ;;
+  ;; Equivalently, reverse the tail Z, Y, X, ... so we have A, ... X, Y, Z,
+  ;; and swap A with the next larger element.
+
+  (let ((n (length v)))
+    (cond ((< n 2) nil)
+         (t (let* ((k (1- n))
+                   (x (aref v k)))
+              (loop (when (zerop k) (return-from next-permutation nil))
+                    (decf k)
+                    (let ((y (aref v k)))
+                      (when (prog1 (< y x)
+                              (setf x y))
+                        (return))))
+              (do ((i (1+ k) (1+ i))
+                   (j (1- n) (1- j)))
+                  ((> i j))
+                (rotatef (aref v i) (aref v j)))
+              (do ((i (- n 2) (1- i)))
+                  ((or (<= i k) (< (aref v i) x))
+                   (rotatef (aref v k) (aref v (1+ i)))))
+              v)))))
+
+(defun make-index-mask (w mask-expr)
+  "Construct a bitmask based on bitwise properties of the bit indices.
+
+   The function returns a W-bit mask in which each bit is set if MASK-EXPR
+   of true of the bit's index.  MASK-EXPR may be one of the following:
+
+     * I -- an integer I is true if bit I of the bit index is set;
+     * (not EXPR) -- is true if EXPR is false;
+     * (and EXPR EXPR ...) -- is true if all of the EXPRs are true; and
+     * (or EXPR EXPR ...) -- is true if any of the EXPRs is true."
+
+  (let ((max-bit (1- (integer-length (1- w))))
+       (mask 0))
+    (dotimes (i w mask)
+      (labels ((interpret (expr)
+                (cond ((and (integerp expr) (<= 0 expr max-bit))
+                       (logbitp expr i))
+                      ((and (consp expr) (eq (car expr) 'not)
+                            (null (cddr expr)))
+                       (not (interpret (cadr expr))))
+                      ((and (consp expr) (eq (car expr) 'and))
+                       (every #'interpret (cdr expr)))
+                      ((and (consp expr) (eq (car expr) 'or))
+                       (some #'interpret (cdr expr)))
+                      (t
+                       (error "unknown mask expression ~S" expr)))))
+       (when (interpret mask-expr)
+         (setf (ldb (byte 1 i) mask) 1))))))
+
+(defun make-permutation-network (w steps)
+  "Construct a permutation network.
+
+   The integer W gives the number of bits to be acted upon.  The STEPS are a
+   list of instructions of the following forms:
+
+     * (SHIFT . MASK) -- a pair of integers is treated literally;
+
+     * (SHIFT MASK-EXPR) -- the SHIFT is literal, but the MASK-EXPR is
+       processed by `make-index-mask' to calculate the mask;
+
+     * (:invert I) -- make an instruction which inverts the sense of the
+       index bit I;
+
+     * (:exchange I J) -- make an instruction which exchanges index bits I
+       and J; or
+
+     * (:exchange-invert I J) -- make an instruction which exchanges and
+       inverts index bits I and J.
+
+   The output is a list of primitive (SHIFT . MASK) steps, indicating that
+   the bits of the input selected by MASK are to be swapped with the bits
+   selected by (ash MASK SHIFT)."
+
+  (let ((max-mask (1- (ash 1 w)))
+       (max-shift (1- w))
+       (max-bit (1- (integer-length (1- w))))
+       (list nil))
+    (dolist (step steps)
+      (cond ((and (consp step)
+                 (integerp (car step)) (<= 0 (car step) max-shift)
+                 (integerp (cdr step)) (<= 0 (cdr step) max-mask))
+            (push step list))
+           ((and (consp step)
+                 (integerp (car step)) (<= 0 (car step) max-shift)
+                 (null (cddr step)))
+            (push (cons (car step) (make-index-mask w (cadr step))) list))
+           ((and (consp step)
+                 (eq (car step) :invert)
+                 (integerp (cadr step)) (<= 0 (cadr step) max-bit)
+                 (null (cddr step)))
+            (let ((i (cadr step)))
+              (push (cons (ash 1 i) (make-index-mask w `(not ,i))) list)))
+           ((and (consp step)
+                 (eq (car step) :exchange)
+                 (integerp (cadr step)) (integerp (caddr step))
+                 (<= 0 (cadr step) (caddr step) max-bit)
+                 (null (cdddr step)))
+            (let ((i (cadr step)) (j (caddr step)))
+              (push (cons (- (ash 1 j) (ash 1 i))
+                          (make-index-mask w `(and ,i (not ,j))))
+                    list)))
+           ((and (consp step)
+                 (eq (car step) :exchange-invert)
+                 (integerp (cadr step)) (integerp (caddr step))
+                 (<= 0 (cadr step) (caddr step) max-bit)
+                 (null (cdddr step)))
+            (let ((i (cadr step)) (j (caddr step)))
+              (push (cons (+ (ash 1 i) (ash 1 j))
+                          (make-index-mask w `(and (not ,i) (not ,j))))
+                    list)))
+           (t
+            (error "unknown permutation step ~S" step))))
+    (nreverse list)))
+
+;;;--------------------------------------------------------------------------
+;;; Permutation network diagnostics.
+
+(defun print-permutation-network (steps &optional (stream *standard-output*))
+  "Print a description of the permutation network STEPS to STREAM.
+
+   A permutation network consists of a list of pairs
+
+       (SHIFT . MASK)
+
+   indicating that the bits selected by MASK, and those SHIFT bits to the
+   left, should be exchanged.
+
+   The output is intended to be human-readable and is subject to change."
+
+  (let ((shiftwd 1) (maskwd 2))
+
+    ;; Determine suitable print widths for shifts and masks.
+    (dolist (step steps)
+      (let ((shift (car step)) (mask (cdr step)))
+       (let ((swd (1+ (floor (log shift 10))))
+             (mwd (ash 1 (- (integer-length (1- (integer-length mask)))
+                            2))))
+         (when (> swd shiftwd) (setf shiftwd swd))
+         (when (> mwd maskwd) (setf maskwd mwd)))))
+
+    ;; Print the display.
+    (pprint-logical-block (stream steps :prefix "(" :suffix ")")
+      (let ((first t))
+       (dolist (step steps)
+         (let ((shift (car step)) (mask (cdr step)))
+
+           ;; Separate entries with newlines.
+           (cond (first (setf first nil))
+                 (t (pprint-newline :mandatory stream)))
+
+           (let ((swaps nil))
+
+             ;; Determine the list of exchanges implied by the mask.
+             (dotimes (i (integer-length mask))
+               (when (logbitp i mask)
+                 (push (cons i (+ i shift)) swaps)))
+             (setf swaps (nreverse swaps))
+
+             ;; Print the entry.
+             (format stream "~@<(~;~vD #x~(~v,'0X~) ~8I~:@_~W~;)~:>"
+                     shiftwd shift maskwd mask swaps))))))
+
+    ;; Print a final newline following the close parenthesis.
+    (terpri stream)))
+
+(defun demonstrate-permutation-network
+    (n steps &optional reference (stream *standard-output*))
+  "Print, on STREAM, a demonstration of the permutation STEPS.
+
+   Begin, on the left, with the integers from 0 up to N - 1.  For each
+   (SHIFT . MASK) element in STEPS, print an additional column showing the
+   effect of that step on the vector.  If REFERENCE is not nil, then it
+   should be a vector of length at least N: on the right, print the REFERENCE
+   vector, showing where the result of the permutation STEPS differs from the
+   REFERENCE.  Return non-nil if the output matches the reference; return nil
+   if the output doesn't match, or no reference was supplied."
+
+  (let ((v (make-array n)))
+
+    ;; Initialize a vector of lists which will record, for each step in the
+    ;; permutation network, which value is in that position.  The lists are
+    ;; reversed, so the `current' value is at the front.
+    (dotimes (i n) (setf (aref v i) (cons i nil)))
+
+    ;; Work through the permutation steps updating the vector.
+    (dolist (step steps)
+      (let ((shift (car step)) (mask (cdr step)))
+
+       (dotimes (i n) (push (car (aref v i)) (aref v i)))
+
+       (dotimes (i n)
+         (when (logbitp i mask)
+           (rotatef (car (aref v i))
+                    (car (aref v (+ i shift))))))))
+
+    ;; Print the result.
+    (let ((ok (not (null reference))))
+      (dotimes (i n)
+       (let* ((entry (aref v i))
+              (final (car entry)))
+         (format stream "~{ ~7D~}" (reverse entry))
+         (when reference
+           (let* ((want (aref reference i))
+                  (match (eql final want)))
+             (format stream " ~:[/=~;==~] ~7D" match want)
+             (unless match (setf ok nil))))
+         (terpri stream)))
+      (when reference
+       (format stream "~:[FAIL~;pass~]~%" ok))
+      ok)))
+
+;;;--------------------------------------------------------------------------
+;;; Examples and useful runes.
+
+#+example
+(let* ((ip #(58 50 42 34 26 18 10  2
+            60 52 44 36 28 20 12  4
+            62 54 46 38 30 22 14  6
+            64 56 48 40 32 24 16  8
+            57 49 41 33 25 17  9  1
+            59 51 43 35 27 19 11  3
+            61 53 45 37 29 21 13  5
+            63 55 47 39 31 23 15  7))
+       (fixed-ip (map '(vector fixnum)
+                     (lambda (x) (- 64 x))
+                     (reverse ip)))
+       (trad-network
+       (make-permutation-network
+        64                             ;  5  4  3  2  1  0
+        '((:exchange-invert 2 5)       ; ~2  4  3 ~5  1  0
+          (:exchange-invert 1 4)       ; ~2 ~1  3 ~5 ~4  0
+          (:exchange-invert 0 3)       ; ~2 ~1 ~0 ~5 ~4 ~3
+          (:exchange-invert 3 4)       ; ~2  0  1 ~5 ~4 ~3
+          (:exchange-invert 4 5)))))   ; ~0  2  1 ~5 ~4 ~3
+
+  (fresh-line)
+
+  (print-permutation-network trad-network)
+  (demonstrate-permutation-network 64 trad-network fixed-ip))