; wheel factorization

(define (td-factors n)
  (let loop ((n n) (x 2) (fs '()))
    (cond ((< n (* x x)) (reverse (cons n fs)))
          ((zero? (modulo n x)) (loop (/ n x) x (cons x fs)))
          (else (loop n (+ x 1) fs)))))

(define (primes n)
  (let* ((max-index (quotient (- n 3) 2))
         (v (make-vector (+ 1 max-index) #t)))
    (let loop ((i 0) (primes '(2)))
      (cond ((< max-index i) (reverse primes))
            ((vector-ref v i)
              (let ((prime (+ i i 3)))
                (do ((j (+ 3 (* 3 i)) (+ j prime)))
                    ((< max-index j))
                  (vector-set! v j #f))
              (loop (+ 1 i) (cons prime primes))))
            (else (loop (+ 1 i) primes))))))

(define (last-pair xs)
  (if (null? (cdr xs)) xs
    (last-pair (cdr xs))))

(define (cycle xs)
  (set-cdr! (last-pair xs) xs) xs)

(define (totatives n)
  (let loop ((i n) (ts '()))
    (cond ((zero? i) ts)
          ((= (gcd i n) 1)
            (loop (- i 1) (cons i ts)))
          (else (loop (- i 1) ts)))))

(define (diffs xs)
  (let loop ((x (car xs)) (xs (cdr xs)) (ds '()))
    (if (null? xs) (reverse ds)
      (loop (car xs) (cdr xs) (cons (- (car xs) x) ds)))))

(define (wheel n)
  (let* ((ps (primes n))
         (t (apply * (cdr (reverse ps))))
         (ts (totatives t))
         (ds (diffs ts)))
    (append (diffs ps)
            (cycle (append (cdr ds)
                           (list 2)
                           (list (car ds)))))))

(define wheel-factors
  (let ((w (wheel 11)))
    (lambda (n)
      (let loop ((n n) (i 2) (fs '()) (w w))
        (cond ((< n (* i i)) (reverse (cons n fs)))
              ((zero? (modulo n i))
                (loop (quotient n i) i (cons i fs) w))
              (else (loop n (+ i (car w)) fs (cdr w))))))))

(display (td-factors 600851475143))
(newline)
(display (wheel-factors 600851475143))

(71 839 1471 6857)
(71 839 1471 6857)