talk about sieving sigma

notebooks/problem0021.ipynb

@@ -67,6 +67,61 @@
     "\n",
     "Therefore, if you have the number's factorization (see [problem 3](https://projecteuler.net/problem=3)), you can use it to compute the sum of its divisors.\n",
     "\n",
+    "## Sieving the sums of divisors\n",
+    "Since we need all the sums of divisors up to 10000, instead of factoring each number individually, we could sieve the values of $\\sigma$. This is possible because the sum of divisors function is multiplicative."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 2,
+   "id": "4cbc68db",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def sum_of_divisors_sieve(limit):\n",
+    "    dsum = [1 for _ in range(0, limit)]\n",
+    "    \n",
+    "    for n in range(0, limit):\n",
+    "        if n == 0 or n == 1 or dsum[n] != 1:\n",
+    "            yield dsum[n]\n",
+    "            continue\n",
+    "            \n",
+    "        m = 1\n",
+    "        while m * n < limit:\n",
+    "            dsum[m*n] = dsum[m] + m*n\n",
+    "            for k in range(2 * m * n, limit, m * n):\n",
+    "                if k % (m*n^2) != 0:\n",
+    "                    dsum[k] *= dsum[m*n]\n",
+    "                    \n",
+    "            m *= n\n",
+    "        \n",
+    "        yield dsum[n]"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "fdbabc18",
+   "metadata": {},
+   "source": [
+    "Here's a check to make sure the sieve matches the values from SageMath."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 3,
+   "id": "0539d1fc",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "s = list(sum_of_divisors_sieve(10000))\n",
+    "assert all(s[n] == sigma(n) for n in range(1, 10000))"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "a847f754",
+   "metadata": {},
+   "source": [
     "## Relevent sequences\n",
     "* Amicable numbers: [A063990](https://oeis.org/A063990)"
    ]