```;;; -*- mode: lisp; syntax: common-lisp -*-

;;; Meiji Cube Puzzle Solver
;;; by Peter Salvi <vukung@yahoo.com> 2007

;;; Usage example:
;; (meiji-cube:solve 15)
;; (meiji-cube:vtk-cube (first *) #p"cube.vtk")

;;; Or, for a complete solution set:
;; (dolist (n '(12 13 14 15))
;;   (let ((solutions (meiji-cube:solve n)))
;;     (dotimes (i (length solutions))
;;       (meiji-cube:vtk-cube (nth i solutions)
;;                         (format nil "cube~d-~2,'0d.vtk" n i)))))

;;; Package

(in-package :cl-user)

(defpackage :meiji-cube
(:use :common-lisp)
(:export :solve :vtk-cube))

(in-package :meiji-cube)

;;; Representation

;;; Piece representation is (X Y Z SHORT-AXIS), so e.g. (1 1 0 1) means
;;; the piece with unit cubes at (1 1 0), (1 1 1), (2 1 0) and (2 1 1).

(declaim (inline make-piece get-pos short-axis))

(defun make-piece (x y z short)
"Make a piece with corner position (X, Y, Z) and short-axis SHORT."
(list x y z short))

(defun get-pos (piece)
"Get the corner position of PIECE."
(subseq piece 0 3))

(defun short-axis (piece)
"Returns the short axis of PIECE."
(nth 3 piece))

;;; Utilities

(defun in-piece-p (piece pos)
"POS is inside PIECE."
(every #'(lambda (n)
(let ((root (nth n piece)))
(<= root
(nth n pos)
(+ root (if (= (short-axis piece) n) 0 1)))))
'(0 1 2)))

(defun emptyp (pieces pos)
"The block at POS is empty, according to PIECES."
(when (every #'(lambda (x) (<= 0 x 3)) pos)
(notany #'(lambda (piece) (in-piece-p piece pos)) pieces)))

(defun pos+ (pos n i)
"Add I to the Nth coordinate of POS."
(mapcar #'(lambda (k) (if (= k n) (+ (nth k pos) i) (nth k pos))) '(0 1 2)))

(defun free-axes (used short)
"Return the two axes that are not USED, the first one not being SHORT."
(let ((axes (remove used '(0 1 2))))
(if (= (first axes) short) (reverse axes) axes)))

(defun empty-box-p (pieces pos axes lengths)
"The box with one corner at POS extending LENGTHS-long in AXES is empty."
(every #'(lambda (x)
(every #'(lambda (y)
(emptyp pieces (pos+ (pos+ pos (first axes) x)
(second axes) y)))
(loop for i from 0 below (second lengths) collect i)))
(loop for i from 0 below (first lengths) collect i)))

(defun movablep (pieces piece &optional (axes '(0 1 2)))
"PIECE can be moved in a configuration described by PIECES.
Movement is only checked along AXES."
(let ((pos (get-pos piece)))
(some #'(lambda (n)
(if (= (short-axis piece) n)
(let ((axes (remove n '(0 1 2))))
(or (empty-box-p pieces (pos+ pos n -1) axes '(2 2))
(empty-box-p pieces (pos+ pos n 1) axes '(2 2))))
(let ((axes (free-axes n (short-axis piece))))
(or (empty-box-p pieces (pos+ pos n -1) axes '(2 1))
(empty-box-p pieces (pos+ pos n 2) axes '(2 1))))))
axes)))

(defun solidp (pieces)
"None of the pieces can be moved."
(notany #'(lambda (piece) (movablep pieces piece)) pieces))

;;; Redundancy check

(defun swap-axes (cube axes)
"Swap the AXES of CUBE."
(mapcar #'(lambda (piece)
(let ((pos (get-pos piece)))
(make-piece (nth (nth 0 axes) pos)
(nth (nth 1 axes) pos)
(nth (nth 2 axes) pos)
(position (short-axis piece) axes))))
cube))

(defun negate-axes (cube x y z)
"Negate those axes that evaluate to T."
(mapcar #'(lambda (piece)
(let ((pos (get-pos piece))
(short (short-axis piece)))
(make-piece
(if x (- 3 (nth 0 pos) (if (= short 0) 0 1)) (nth 0 pos))
(if y (- 3 (nth 1 pos) (if (= short 1) 0 1)) (nth 1 pos))
(if z (- 3 (nth 2 pos) (if (= short 2) 0 1)) (nth 2 pos))
short)))
cube))

(defun transformed-cubes (cube)
"Generate a list of cubes that are considered equivalent to CUBE."
(mapcan #'(lambda (axes)
(let ((swapped (swap-axes cube axes)))
(mapcan #'(lambda (x)
(mapcan #'(lambda (y)
(mapcar #'(lambda (z)
(negate-axes swapped
x y z))
'(nil t)))
'(nil t)))
'(nil t))))
'((0 1 2) (0 2 1) (1 0 2) (1 2 0) (2 0 1) (2 1 0))))

(defun equal-cube-p (cube1 cube2)
"CUBE1 and CUBE2 have the same set of pieces."
(and (= (length cube1) (length cube2))
(every #'(lambda (piece) (member piece cube2 :test #'equal)) cube1)))

(defun equivalentp (cube1 cube2)
"CUBE1 and CUBE2 are equivalent - there exists some rotational
and/or mirroring transformation that transforms one into the other."
(member cube1 (transformed-cubes cube2) :test #'equal-cube-p))

(defun remove-redundant (solutions)
"Solve the puzzle with N pieces, return a list of all distinct solutions."
(do ((lst solutions (rest lst))
result)
((null lst) result)
(when (notany #'(lambda (cube) (equivalentp (first lst) cube)) result)
(push (first lst) result))))

;;; Solution

(defparameter *positions*
(loop with result
for z from 0 to 3 do
(loop for y from 0 to 3 do
(loop for x from 0 to 3 do
(loop for short from 2 downto 0
for axes = (remove short '(0 1 2))
for piece = (make-piece x y z short)
unless (or (= (nth (first axes) (get-pos piece)) 3)
(= (nth (second axes) (get-pos piece)) 3))
do (push piece result))))
finally (return (nreverse result)))
"All possible positions.")

(defun placeablep (pieces piece)
"PIECE is placable, that is, it fits in the cube, there is nothing
in the positions it will take, and it wouldn't fall down in the Z-direction."
(let ((pos (get-pos piece)))
(and (empty-box-p pieces pos (remove (short-axis piece) '(0 1 2)) '(2 2))
(not (if (= (short-axis piece) 2)
(empty-box-p pieces (pos+ pos 2 -1) '(0 1) '(2 2))
(empty-box-p pieces (pos+ pos 2 -1)
(if (= (short-axis piece) 0) '(1 0) '(0 1))
'(2 1)))))))

(defun solve (n)
"Solve the puzzle with N pieces, return a list of all distinct solutions."
(let (result)
(labels ((generate (n lst &optional pieces)
(if (= n 0)
(when (solidp pieces)
(format t "Solution:~%~a~%" pieces)
(push pieces result))
(do ((lst lst (rest lst)))
((null lst))
(when (placeablep pieces (first lst))
(generate (1- n) (rest lst)
(cons (first lst) pieces)))))))
(generate n *positions*))
(remove-redundant result)))

;;; VTK Output

(defun other-corner (piece)
"The far corner of PIECE.
Note that this position is _outside_ the PIECE."
(let ((pos (get-pos piece)))
(list (+ (nth 0 pos) (if (= (short-axis piece) 0) 1 2))
(+ (nth 1 pos) (if (= (short-axis piece) 1) 1 2))
(+ (nth 2 pos) (if (= (short-axis piece) 2) 1 2)))))

(defun vtk-box (stream pos1 pos2 &optional (x 0.1))
"Writes the vertex coordinates of a box given by its two opposing
corners (POS1 and POS2) to STREAM. The box is made smaller by X on every side."
(format stream "~f ~f ~f~%"
(+ (first pos1) x) (+ (second pos1) x) (+ (third pos1) x))
(format stream "~f ~f ~f~%"
(- (first pos2) x) (+ (second pos1) x) (+ (third pos1) x))
(format stream "~f ~f ~f~%"
(- (first pos2) x) (+ (second pos1) x) (- (third pos2) x))
(format stream "~f ~f ~f~%"
(+ (first pos1) x) (+ (second pos1) x) (- (third pos2) x))
(format stream "~f ~f ~f~%"
(+ (first pos1) x) (- (second pos2) x) (+ (third pos1) x))
(format stream "~f ~f ~f~%"
(- (first pos2) x) (- (second pos2) x) (+ (third pos1) x))
(format stream "~f ~f ~f~%"
(- (first pos2) x) (- (second pos2) x) (- (third pos2) x))
(format stream "~f ~f ~f~%"
(+ (first pos1) x) (- (second pos2) x) (- (third pos2) x)))

(defun vtk-cube (pieces filename)
"Write a 3D representation of PIECES in FILENAME (VTK format)."
(with-open-file (s filename :direction :output :if-exists :supersede)
(format s "# vtk DataFile Version 1.0~%~
Meiji Cube~%~
ASCII~%~
DATASET POLYDATA~%~%~
POINTS ~d float~%"
(* (length pieces) 8))
(dolist (piece pieces)
(vtk-box s (get-pos piece) (other-corner piece)))
(format s "~%POLYGONS ~d ~d~%"
(* (length pieces) 6) (* (length pieces) 30))
(dotimes (i (length pieces))
(format s "4 ~d ~d ~d ~d~%"
(+ (* i 8) 0) (+ (* i 8) 1) (+ (* i 8) 2) (+ (* i 8) 3))
(format s "4 ~d ~d ~d ~d~%"
(+ (* i 8) 4) (+ (* i 8) 5) (+ (* i 8) 6) (+ (* i 8) 7))
(format s "4 ~d ~d ~d ~d~%"
(+ (* i 8) 0) (+ (* i 8) 1) (+ (* i 8) 5) (+ (* i 8) 4))
(format s "4 ~d ~d ~d ~d~%"
(+ (* i 8) 3) (+ (* i 8) 2) (+ (* i 8) 6) (+ (* i 8) 7))
(format s "4 ~d ~d ~d ~d~%"
(+ (* i 8) 1) (+ (* i 8) 2) (+ (* i 8) 6) (+ (* i 8) 5))
(format s "4 ~d ~d ~d ~d~%"
(+ (* i 8) 0) (+ (* i 8) 3) (+ (* i 8) 7) (+ (* i 8) 4)))))
```