Project Euler 21

http://projecteuler.net/index.php?section=problems&id=21


Sequenceは遅延評価らしいです。大きなSequenceもメモリを食うことはありません。

let N = int64 1e8
printfn "%d" (Seq.sum (seq { 1L..N }))

Seq.toListもListを中に持って一つずつ出していくようです。

let N = int64 3e6
printfn "%d" (Seq.sum (Seq.toList [1L..N]))

ちなみに、これでも100MB以上メモリを食うようです。

さて、問題のdをエラトステネスのふるい的に求めます。ほぼここに書いた素因数分解の方法と同様です。このときに100までの素数で割っていくのですが、100を求めるのにsqrtを使うのがなんとなく嫌なので、takeWhileを使っています。無限列を使えればいいのですが、それを気軽に書けるシンタックスが無い様なので、十分に長い有限の列を代替にしています。

Seq.initInfinite (fun n -> n + 2)   // 2から始まる無限列
seq { 2..max_n }     // 有限列だがこのほうが簡単に書ける
let rec pow n e = if e = 0 then 1 else n * (pow n (e - 1))

let rec calc_exp n p =
    if n % p <> 0 then
        (0, n)
    else
        let e, m = calc_exp (n / p) p
        (e + 1, m)

let sieve max_n =
   let a = [| 0..max_n |]
   let d = Array.create (max_n + 1) 1
   for p in Seq.takeWhile (fun n -> n * n <= max_n)
               (Seq.filter (fun n -> d.[n] = 1) (seq { 2..max_n })) do
      for n in seq { p..p..max_n } do
        let e, m = calc_exp a.[n] p
        a.[n] <- m
         d.[n] <- d.[n] * (((pow p (e + 1)) - 1) / (p - 1))
   
   for n in Seq.filter (fun n -> a.[n] <> 1) (seq { 2..max_n }) do
    d.[n] <- d.[n] * (a.[n] + 1)
   
   for n in seq { 1..max_n } do
    d.[n] <- d.[n] - n
   
   d

let N = int 1e4
let d = sieve N
let is_amicable n = d.[n] < n && d.[d.[n]] = n
printfn "%d" (Seq.sum (Seq.map (fun n -> n + d.[n])
                    (Seq.filter is_amicable (seq { 2..N }))))