(define exit (let* ((s (cons 0 0))
   (z (call-with-current-continuation (lambda (e) (cons s e)))))
      (if (and (pair? z) (eq? (car z) s)) (cdr z) exit)))

(define mer (lambda (x) (let ((r (call-with-current-continuation
  (lambda(d) (cons #f (lambda(y) (d (cons #t y))))))))
   (if (car r) (begin (write (list "Exception:" x (cdr r))) (newline) (exit 0))
   (cdr r)))))

(define (merge-assim A B comp assim) (let mg ((A A)(B B)) 
(if (null? A) B (if (null? B) A (let ((t (comp (car A) (car B))))
    (cond ((negative? t) (cons (car A) (mg (cdr A) B)))
      ((positive? t) (cons (car B) (mg A (cdr B))))
      (else (cons (assim (car A)(car B)) (mg (cdr A) (cdr B))))))))))
(define (ms comp assim) (lambda (A) (let nm ((ll (map list A)))
  (if (pair? ll) (if (pair? (cdr ll))
  (nm (let srt ((a ll)) (if (null? a) '() (if (null? (cdr a)) a
     (cons (merge-assim (car a) (cadr a) comp assim) (srt (cddr a)))))))
  (car ll)) '()))))

(define (Do k proc)(let loop ((j 0)) (if (< j k)
  (begin (proc j) (loop (+ j 1))))))
(define (vg m n) (let* (
   (L (* m n))
   (s (make-vector L '()))
   (si (ms - (mer "Twins")))
   (f 1/4)
   (mf (- f))
   (term (list (cons L 0)))
   (ndx (lambda (i j)(+ (* n i) j)))
   (d (lambda n (cons (car n)
      (let sd ((n (si (cdr n)))) (if (null? n) term
      (cons (cons (car n) f) (sd (cdr n)))))))))
  (vector-set! s (ndx 0 0)
      (d 0 (ndx 0 1) (ndx 1 0)))
  (vector-set! s (ndx 0 (- n 1))
      (d 0 (ndx 1 (- n 1)) (ndx 0 (- n 2))))
  (vector-set! s (ndx (- m 1) 0)
      (d 0 (ndx (- m 2) 0) (ndx (- m 1) 1)))
  (vector-set! s (ndx (- m 1)(- n 1))
      (d 0 (ndx (- m 2)(- n 1))(ndx (- m 1)(- n 2))))
  (Do (- n 2) (lambda (jm) (let ((j (+ jm 1)))
     (vector-set! s (ndx 0 j) (d mf (ndx 1 j) (ndx 0 (+ j 1)) (ndx 0 (- j 1))))
     (vector-set! s (ndx (- m 1) j)
         (d mf (ndx (- m 2) j) (ndx (- m 1) (+ j 1)) (ndx (- m 1) (- j 1)))))))
  (Do (- m 2)(lambda (jm) (let ((j (+ jm 1)))
     (vector-set! s (ndx j 0) (d f (ndx j 1) (ndx (+ j 1) 0) (ndx (- j 1) 0)))
     (vector-set! s (ndx j (- n 1))
         (d f (ndx j (- n 2)) (ndx (+ j 1) (- n 1)) (ndx (- j 1) (- n 1)))))))
  (Do (- m 2) (lambda (im) (let ((i (+ im 1)))
  (Do (- n 2) (lambda (jm) (let ((j (+ jm 1)))
     (vector-set! s (ndx i j)
        (d 0 (ndx (- i 1) j) (ndx (+ i 1) j) (ndx i (- j 1)) (ndx i (+ j 1))))))))))
  (Do L (lambda (x) (if (null? (vector-ref s x)) (write (list x)))))
  (cons (list m n) s)
  ))
 
(define (sl a) ((ms (lambda (x y) (let ((s (- (car x) (car y))))
  (if (zero? s) (- (cdr x) (cdr y)) s))) (mer "bump")) a))
(define (retire s n)(let* ((s (cdr s))(lc (vector-ref s n))(n2 (vector-length s)))
  (if (not (number? (car lc))) ((mer "Dead already") "Pooh"))
  (let qx ((Lc (cdr lc))) (let* ((N (car Lc))(id (car N)))
   (if (not (= id n2)) (begin (let* (
     (nlc (vector-ref s id))
     (f (cdr (assq n (cdr nlc))))
     (g (cdr (assq id (cdr lc))))
     (r (/ (- 1 (* f g)))))
   (vector-set! s id (cons (* r (+ (car nlc) (* f (car lc))))
   (let xl ((lC (cdr lc))(nlC (cdr nlc)))
     (if (= (caar lC)(caar nlC)) (if (= (caar lC) n2) lC (cons (cons (caar lC)
         (* r (+ (* f (cdar lC)) (cdar nlC)))) (xl (cdr lC) (cdr nlC))))
         (if (< (caar lC)(caar nlC))
           (if (= (caar lC) id) (xl (cdr lC) nlC)
             (cons (cons (caar lC) (* r f (cdar lC))) (xl (cdr lC) nlC)))
           (if (= (caar nlC) n) (xl lC (cdr nlC))
             (cons (cons (caar nlC) (* r (cdar nlC))) (xl lC (cdr nlC))))))))))
       (qx (cdr Lc))))))
 (vector-set! s n (cons "d" lc))))

(define (ln s er) (let ((c '())(n2 (vector-length (cdr s))))
  (if (not (= n2 (apply * (car s)))) (er (list "Wrong Length" n2)))
  (Do n2 (lambda (j) (let ((x (vector-ref (cdr s) j)))
    (if (number? (car x))      
    (let zl ((el (cdr x))(p -1)) (let ((e (caar el))) 
      (if (= e j) (er (cons "Reflexive!" j)))
      (if (not (= e n2))
      (begin (if (>= p e) ((mer "sort") (cons e p)))
        (set! c (cons (cons j e) c)) (zl (cdr el) e)))))
      ))))
   (let ((d (sl c))(e (sl (map (lambda (x) (cons (cdr x)(car x))) c))))
 (if (not (equal? d e)) (er (list "non-symmetric" d e))))))

 (define (bs s n) (let ((lc (vector-ref (cdr s) n))(n2 (vector-length (cdr s))))
  (if (string=? (car lc) "d") (begin
  (vector-set! (cdr s) n (cons "k" (let sum ((t (cddr lc))) (if (= (caar t) n2)
    (cadr lc) (+ (* (cdar t) (let* ((id (caar t))(olc (vector-ref (cdr s) id)))
      (if (string=? "k" (car olc)) (cdr olc)
         ((mer "Order") olc)))) (sum (cdr t))))))) lc)
  ((mer "order") lc))))

(define (vector-copy v) (let* ((l (vector-length v))(x (make-vector l)))
  (Do l (lambda (j) (vector-set! x j (vector-ref v j)))) x))
(define s (vg 10 10))
(define ge (mer "Problem!"))
(define (solve s er) (let ((n (vector-length (cdr s)))
     (b (cons (car s)(vector-copy (cdr s)))))
  (ln s er)
  (Do n (lambda (j) (retire b j)))
  (Do n (lambda (j) (bs b (- (- n 1) j)))) b))

(define r (solve s ge))
; and to test for symetry of solution;
(define (vector-copy v) (let* ((l (vector-length v))(x (make-vector l)))
  (Do l (lambda (j) (vector-set! x j (vector-ref v j)))) x))
(define (state->float s)(cons (car s)(let* ((s (cdr s))
  (nv (vector-copy s)) (l (vector-length s)))
  (Do l (lambda (j) (vector-set! nv j (let ((lc (vector-ref s j)))
    (cons (exact->inexact (car lc)) (let cl ((x (cdr lc)))
      (if (= (caar x) l) x (cons (cons (caar x)
      (exact->inexact (cdar x))) (cl (cdr x)))))))))) nv)))
(define (flipH a)(let* ((ds (car a))(m (car ds))(n (cadr ds))
   (ndx (lambda (i j)(+ (* n i) j)))(b (make-vector (* m n))))
  (Do m (lambda (j) (Do n (lambda (k)
     (vector-set! b (ndx j k) (vector-ref (cdr a) (ndx j (- (- n 1) k))))))))
 (cons ds b)))
(define (flipV a)(let* ((ds (car a))(m (car ds))(n (cadr ds))
   (ndx (lambda (i j)(+ (* n i) j)))(b (make-vector (* m n))))
  (Do m (lambda (j) (Do n (lambda (k)
     (vector-set! b (ndx j k) (vector-ref (cdr a) (ndx (- (- m 1) j) k)))))))
 (cons ds b)))
(define (flipD a)(let* ((ds (car a))(m (car ds))(n (cadr ds))(b (make-vector (* m n))))
  (Do m (lambda (j) (Do n (lambda (k)
     (vector-set! b (+ (* m k) j) (vector-ref (cdr a) (+ k (* n j))))))))
      (cons (cons n (cons m (cddr ds))) b)))
(define (cmpv a b)(let ((L (vector-length (cdr a)))(c 0))
  (Do L (lambda (j) (set! c (max c (abs
     (- (cdr (vector-ref (cdr a) j)) (cdr (vector-ref (cdr b) j)))))))) c))
(define (negv x)(let* ((a (cdr x))(L (vector-length a))(b (make-vector L)))
  (Do L (lambda (j) (vector-set! b j (let ((x (vector-ref a j)))
      (cons (car x)(- (cdr x))))))) (cons (car x) b)))
; r still holds the exact rational solution.
(cmpv r (flipH r))
(cmpv r (flipV r))
(cmpv r (negv (flipD r)))
(cmpv r (solve (state->float s) ge)) ; < (expt 2 -50)
(cmpv r (negv (flipD (solve (vg (cadar s) (caar s)) ge))))

(define s (vg 3 14))
(define r (solve s ge))
(cmpv r (flipH r))
(cmpv r (flipV r))
(cmpv r (solve (state->float s) ge)) ; < (expt 2 -50)
(cmpv r (negv (flipD (solve (vg (cadar s) (caar s)) ge))))

(solve (vg 4 25) ge) ; 160 sec
(solve (vg 25 4) ge) ; 6 sec

(solve (state->float (vg 4 25)) ge)
(solve (state->float (vg 25 4)) ge)
(cmpv (negv (solve (state->float (vg 4 25)) ge)) (flipD (solve (state->float (vg 25 4)) ge)))