Too many remainders

$2016 \equiv 1×(1009 + n) + (1007 - n)$ ,

therefore, if we are applying it for n = 0, 1, 2, ..., 1006, 1007 to

divide 2016 by 1009, 1010, 1011, ..., 2016, we get a quotient of 1 in all cases, and a remainder of 1007, 1006, ... , 2, 1, 0 (all integers between 0 and 1007; 1008 distinct values).

Since a division by a positive integer m < 1009 cannot be greater, than (m - 1), the maximum value of the remainder cannot be bigger, than 1008 - 1 = 1007 (and the minimum being 0), so these possible remainders for the [1,1008] interval are all represented already.

It is easy to see, that we get a quotient of 0 and a remainder of 2016 (our 1009th value) if we divide it by any integer m > 2016:

$2016 \equiv 0×(2016 + k) + 2016 \ , k \in \mathbb {Z} , k > 0$

Now, the sum of all possible remainders:

$(0 + 1 + 2 + ...+ 1006 + 1007) + 2016 = \frac {(0 + 1007) × 1008}{2} + 2016 = \boxed {509544}$

We iterate over all the possible divisors and keep adding the remainders to a set. Since all integers above 2017 will yield 2016 as a remainder, it is sufficient to check upto 2017.

1
2
3
4

remainders = set()
for divisor in xrange(1,2018):
    remainders.add(2016%divisor)
print sum(remainders)