; fermat's method

(define (isqrt n)
  (let loop ((x n) (y (quotient (+ n 1) 2)))
    (if (<= 0 (- y x) 1) x
      (loop y (quotient (+ y (quotient n y)) 2)))))

(define (factors n)
  (if (even? n) (cons 2 (factors (/ n 2)))
    (let ((s (+ (isqrt n) 1)))
      (let loop ((u (+ s s 1)) (v 1) (r (- (* s s) n)))
        (cond ((positive? r) (loop u (+ v 2) (- r v)))
              ((negative? r) (loop (+ u 2) v (+ r u)))
              ((= (- u v) 2) (list (/ (+ u v -2) 2)))
              (else (append (factors (/ (+ u v -2) 2))
                            (factors (/ (- u v) 2)))))))))

(display (factors 600851475143))
