Let $r$ be the remainder when $(a−1)^n + (a+1)^n$ is divided by $a^2$.
For example, if $a = 7$ and $n = 3$, then $r = 42$: $6^3 + 8^3 = 728 \equiv 42 \pmod{49}$.
And as $n$ varies, so too will $r$, but for $a = 7$ it turns out that $r_{max} = 42$.
For $3 \le a \le 1000$, find $\sum r_{max}$.
In [1]:
def rmax(a):
t1 = t2 = 1
aa = a*a
Table = {(1,1):1}
rm = 2 % aa
while 1:
t1 = ((a-1)*t1) % aa
t2 = ((a+1)*t2) % aa
rm = max(rm, (t1+t2)%aa)
if (t1, t2) in Table:
return rm
Table[(t1,t2)] = 1
print(sum(rmax(a) for a in range(3,1001)))