/* 

  Euler #21 in Picat.

  Problem 21
  """
  Let d(n) be defined as the sum of proper divisors of n (numbers less 
  than n which divide evenly into n).
  If d(a) = b and d(b) = a, where a /= b, then a and b are an amicable 
  pair and each of a and b are called amicable numbers.
  
  For example, the proper divisors of 220 are 
  1, 2, 4, 5, 10, 11, 20, 22, 44, 55 and 110; therefore d(220) = 284. 
  The proper divisors of 284 are 1, 2, 4, 71 and 142; so d(284) = 220.
  
  Evaluate the sum of all the amicable numbers under 10000.
  """

  This Picat model was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Picat page: http://www.hakank.org/picat/

*/

import cp.

main => go.


go => time(euler21c).

euler21 => 
  S = new_map(),
  foreach(A in 1..9999) 
     B = sum_divisors3(A),
     C = sum_divisors3(B),
     if A != B, A == C then
        S.put(A, 1),
        S.put(B, 1)
     end
  end,
  println(sum(S.keys())).

% list comprehension
euler21b => 
  Sum = [A : A in 1..9999, 
       B = sum_divisors3(A), 
       C = sum_divisors3(B),
       A != B, A == C].remove_dups().sum(),
  println(Sum).

euler21c =>
  N = 9999,
  S = new_array(N),
  foreach(I in 1..N)
    S[I] := sum_divisors3(I)
  end,
  A = 0,
  foreach(I in 1..N)
    if S[I]<=N, I==S[S[I]], I!=S[I] then
      A := A + I
    end
  end,
  println(A),
  nl.
  

table
sum_divisors2(N) = Sum =>
    D = floor(sqrt(N)),
    Sum1 = 1,
    foreach(I in 2..D, N mod I == 0) 
        Sum1 := Sum1+I,
        if I != N div I then
            Sum1 := Sum1 + N div I
        end
    end,
    Sum = Sum1.

% This recursive version is slightly faster than sum_divisors2/1.
table
sum_divisors3(N) = Sum =>
  sum_divisors3(2,N,1,Sum).

sum_divisors3(I,N,Sum0,Sum), I > floor(sqrt(N)) =>
  Sum = Sum0.

% I is a divisor of N
sum_divisors3(I,N,Sum0,Sum), N mod I == 0 =>
  Sum1 = Sum0 + I,
  (I != N div I -> 
    Sum2 = Sum1 + N div I 
    ; 
    Sum2 = Sum1
  ),
  sum_divisors3(I+1,N,Sum2,Sum).

% I is no divisor of N.
sum_divisors3(I,N,Sum0,Sum) =>
  sum_divisors3(I+1,N,Sum0,Sum).