# Speeding the search for amicable numbers

The following text is from Problem 21 - Project Euler :

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$.

The fastest Mathematica snippet I came up with uses the characterization of amicable numbers as proper involutive points for the divisor sum function:

DSum[n_]:=DivisorSigma[1,n]-n
AmicableCandidates=Select[Range[10000], DSum[DSum[#]]==# &];
AmicableNumbers=Select[AmicableCandidates, DSum[#]!=# &];
Sum[AmicableNumbers[[i]], {i,1,Length[AmicableNumbers]}]


Timing this yields 0.09375s on my machine.

1. Can this code be improved?
2. What are some other efficient ways to compute the sum of all amicable numbers under $10000$?

Since DivisorSigma is Listable, one can use a vectorized approach:

amicable[list_] := Module[{f1, f2},
f1 = DivisorSigma[1, list] - list;
f2 = DivisorSigma[1, f1] - f1;
Pick[list, Unitize[f1-list]BitXor[Unitize[f2-list], 1], 1]
]
amicable[Range[2,10^4]] //Total //RepeatedTiming


{0.045, 31626}

I don't include 1 because DivisorSigma[1, 1] - 1 is 0, and DivisorSigma[1, 0] is not defined. Hence the output of f2 becomes unpacked when 1 is included in the list.

Instant speedup: Change the last line to:

Total[AmicableNumbers]


Also try this:

DSum[n_]:=DivisorSigma[1,n]-n;
Total@Select[Range[10000], (DSum[DSum[#]]==# && DSum[#]!=#) &]


Just for another way:

f[n_] := DivisorSigma[1, n] - n;
g[n_] := NestList[f, n, 2]
Total[Cases[g /@ Range[2, 10000], {i_, j_, i_} /; i!= j]][[1]] //
RepeatedTiming


yields:

{0.11, 31626}