; modular multiplication without overflow

(define (mod-mul a b m) ; assumes a,b < 2^64 and m < 2^63
  (define (shift x k)
    (let loop ((k k) (x x))
      (if (zero? k) x
        (let* ((x (* x 2))
               (x (if (< x m) x (- x m))))
          (loop (- k 1) x)))))
  (let* ((two32 (expt 2 32))
         (a1 (remainder a two32))
         (a2 (quotient a two32))
         (b1 (remainder b two32))
         (b2 (quotient b two32))
         (p11 (modulo (* a1 b1) m))
         (p12 (shift (modulo (* a1 b2) m) 32))
         (p21 (shift (modulo (* a2 b1) m) 32))
         (p22 (shift (modulo (* a2 b2) m) 64))
         (s (+ p11 p12))
         (s (if (< s m) s (- s m)))
         (s (+ s p21))
         (s (if (< s m) s (- s m)))
         (s (+ s p22))
         (s (if (< s m) s (- s m))))
    s))

(define a 17259738289493410580109721600)
(define b 10327523882682224844906430464)
(define m 8631375519702822467321987072)
(display (mod-mul a b m))

7869126927168251407166865408