Project EulerをSchemeで(46-50)

そろそろProject Euler用のモジュールつくったほうがいいかな。

(use util.combinations)
(use gauche.sequence)
(use srfi-1)

(define (integer->list i)
  (letrec ((i->rl (^i (cons (modulo i 10)
			    (if (< i 10)
				'()
				(i->rl (quotient i 10)))))))
    (reverse (i->rl i))))

(define (list->integer l)
  (define (shift-d n)
    (let loop ((d 10))
      (if (< n d) d (loop (* 10 d)))))
  (fold (^(n p) (+ (* p (shift-d n)) n)) 0 l))

(define (primes n)
  (if (<= n 2)
      '()
      (let1 u (truncate (sqrt n))
	(let loop ((ps '(2))
		   (l (unfold (cut < n <>) values (cut + 2 <>) 3)))
	  (let1 m (car l)
	    (if (> m u)
		(append (reverse ps) l)
		(loop (cons m ps)
		      (remove (^x (zero? (modulo x m))) l))))))))

; 試し割り
(define (prime? n)
  (cond ((= n 2) #t)
	((or (< n 2) (zero? (modulo n 2))) #f)
	(else (let1 m (floor->exact (sqrt n))
		(let loop ((i 3))
		  (cond ((< m i) #t)
			((zero? (modulo n i)) #f)
			(else (loop (+ i 2)))))))))

問46.

平方数の2倍と素数の和で表せない最小の奇合成数を求める問題。

; 総当り
(define (e46)
  (let* ((ps (primes 10000))
	 (prime? ((^h (dolist (p ps)
			(hash-table-put! h p #t))
		      (^n (hash-table-get h n #f)))
		  (make-hash-table)))
	 (twice-a-square? ((^h (dolist (n (iota 100 1))
				 (hash-table-put! h (* 2 n n) #t))
			       (^n (hash-table-get h n #f)))
			   (make-hash-table))))
    (let loop ((n 3))
      (if (prime? n)
	  (loop (+ n 2))
	  (let loop2 ((ps ps))
	    (cond ((null? ps) n)
		  ((twice-a-square? (- n (car ps)))
		   (loop (+ n 2)))
		  (else (loop2 (cdr ps)))))))))

問47.

連続する4つの数がそれぞれ4つの異なる素因数を持つ場合を考え, 連続する数の中で最小のものを求める問題。

; 連続する4数で素因数の個数が全て4となっているものを探す
; 遅い
(define (e47)
  (define (factor-count n ps)
    (let loop ((ans 0) (n n) (ps ps))
      (cond ((or (< n (car ps)) (null? ps)) ans)
	    ((zero? (modulo n (car ps)))
	     (loop (+ ans 1) (/ n (car ps)) (cdr ps)))
	    (else (loop ans n (cdr ps))))))
  (let1 ps (primes 200000)
    (let loop ((i 210) (fc '()))
      (cond ((< i 214) (loop (+ i 1)
			     (cons (factor-count i ps) fc)))
	    ((= 4 (car fc) (cadr fc) (caddr fc) (cadddr fc)) (- i 4))
	    (else (loop (+ i 1)
			(cons (factor-count i ps) fc)))))))

問48.

1^1 + 2^2 + 3^3 + ... + 1000^1000 の最後の10桁を求める問題。

; 全て足して下10桁のみ取り出す
(define (e48)
  (mod (apply + (map (^n (expt n n))
		     (iota 1000 1)))
       (expt 10 10)))

問49.

それぞれ素数で各項は他の項の置換で表せ、等差数列となるような3数を求める問題。

; 1000<p<10000 なる素数pと同じ数字の組み合わせで作られる数iが
; 素数かつ p<i で
; さらにp i j が等差数列となるようにjをとったとき
; jが素数かつpと同じ数の組み合わせで作られていれば
; そのp i jが求める答えになる
(define (e49)
  (call/cc (^(return)
	     (for-each
	      (^p (for-each
		   (^i (when (and (< p i) (prime? i))
			 (let1 j (+ i (- i p))
			   (when (and (< j 10000)
				      (prime? j)
				      (= (list->integer (sort (integer->list p)))
					 (list->integer (sort (integer->list j)))))
			     (return (list->integer (list p i j)))))))
		   (map list->integer (permutations* (integer->list p)))))
	      ($ delete 1487 $ filter (cut < 1000 <>) $ primes 10000)))))

問50.

連続する素数の和で表したときに最長になる100万未満の素数を求める問題。

; 最長は2以降の素数を100万を超えないように順番に足していったときの長さなので
; そこから素数が見つかるまで長さを縮めていく
(define (e50)
  (let* ((ps (primes 5000))
	 (v (make-vector (+ (length ps) 1) 0)))
    (for-each-with-index
     (^(i p) (vector-set! v (+ i 1) (+ p (vector-ref v i))))
     ps)
    (let1 max-len (let loop ((i 0))
		    (if (< 1000000 (vector-ref v (+ i 1)))
			i (loop (+ i 1))))
      (call/cc (^(return)
		 (for-each
		  (^l (call/cc (^(break)
				 (for-each
				  (^i (let1 s (- (vector-ref v (+ i l))
						 (vector-ref v i))
					(cond ((< 1000000 s) (break))
					      ((prime? s) (return s)))))
				  (iota (- (+ max-len 1) l) 1)))))
		  (iota max-len max-len -1)))))))