fix(maths): add the leftover prime factor in sum_of_divisors

sum_of_divisors only multiplied in the prime factors it found while
trial-dividing, so any prime factor above the trial-division bound was
dropped: 94 of the first 119 inputs returned the wrong answer
(sum_of_divisors(3) -> 1 instead of 4, sum_of_divisors(12) -> 7 instead
of 28). number_of_divisors in the same file already handles this case
with `if n > 1: div *= 2`.

Also compute with // instead of /: the function accumulated a float and
truncated with int(s), so sum_of_divisors(2**62) returned
9223372036854775808 instead of 9223372036854775807 and large inputs were
off by orders of magnitude (3**40 was wrong by ~10^4). math.isqrt
replaces int(math.sqrt(n)) for the same reason.

Verified against an independent factorization of sigma(n) for
n = 1..50000: 0 mismatches.

Author: bougieL

maths/basic_maths.py

@@ -67,6 +67,14 @@ def sum_of_divisors(n: int) -> int:
     """Calculate Sum of Divisors.
     >>> sum_of_divisors(100)
     217
+    >>> sum_of_divisors(12)
+    28
+    >>> sum_of_divisors(7)
+    8
+    >>> sum_of_divisors(2**40)
+    2199023255551
+    >>> sum_of_divisors(2**62)
+    9223372036854775807
     >>> sum_of_divisors(0)
     Traceback (most recent call last):
         ...
@@ -82,17 +90,21 @@ def sum_of_divisors(n: int) -> int:
     temp = 1
     while n % 2 == 0:
         temp += 1
-        n = int(n / 2)
+        n //= 2
     if temp > 1:
-        s *= (2**temp - 1) / (2 - 1)
-    for i in range(3, int(math.sqrt(n)) + 1, 2):
+        s *= (2**temp - 1) // (2 - 1)
+    for i in range(3, math.isqrt(n) + 1, 2):
         temp = 1
         while n % i == 0:
             temp += 1
-            n = int(n / i)
+            n //= i
         if temp > 1:
-            s *= (i**temp - 1) / (i - 1)
-    return int(s)
+            s *= (i**temp - 1) // (i - 1)
+    # Whatever is left is a single prime factor above the trial-division
+    # bound, so it contributes 1 + n to the sum.
+    if n > 1:
+        s *= n + 1
+    return s
 
 
 def euler_phi(n: int) -> int: