filifa/eulerbooks
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases 3 Wiki Activity

move notebooks to subdir

Browse Source
This commit is contained in:
filifa
2025-04-15 19:26:01 -04:00
parent c0dfb64979
commit 6a8215888b
29 changed files with 3 additions and 3 deletions
Download Patch File Download Diff File Expand all files Collapse all files

61
notebooks/problem0001.ipynb Normal file
View File

@@ -0,0 +1,61 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "fdd3f4eb",
"metadata": {},
"source": [
"# [Multiples of 3 or 5](https://projecteuler.net/problem=1)\n",
"\n",
"Lots of different approaches to this problem! Writing a program may be the easiest, but I always find it interesting when these problems can be solved without a computer doing the heavy lifting, and with a little bit of [summation](https://en.wikipedia.org/wiki/Summation) knowledge and the [inclusion-exclusion principle](https://en.wikipedia.org/wiki/Inclusion%E2%80%93exclusion_principle), we can calculate this answer by hand.\n",
"\n",
"First off, the sum of the multiples of 3 below 1000 is\n",
"$$3 + 6 + 9 + \\cdots + 999$$\n",
"Likewise, the sum of the multiples of 5 below 1000 is\n",
"$$5 + 10 + 15 + \\cdots + 995$$\n",
"Although these sums may seem tedious to calculate by hand, we can apply the distributive property and factor out the 3 and 5, respectively, from each summation.\n",
"$$\\begin{align}\n",
"3 + 6 + 9 + \\cdots + 999 &= 3(1 + 2 + 3 + \\cdots + 333) \\\\\n",
"5 + 10 + 15 + \\cdots + 995 &= 5(1 + 2 + 3 + \\cdots + 199)\n",
"\\end{align}$$\n",
"We are then left with a summation of the form\n",
"$$1 + 2 + 3 + \\cdots + n$$\n",
"in both instances. If we successively evaluate this summation at $n = 1,2,3,\\ldots$ we get the [triangular numbers](https://en.wikipedia.org/wiki/Triangular_number), which have a well-known closed formula:\n",
"$$1 + 2 + 3 + \\cdots + n = \\frac{n(n+1)}{2}$$\n",
"(This formula can be proven valid with [mathematical induction](https://en.wikipedia.org/wiki/Mathematical_induction), but you don't need to formally prove it to gain an intuition for it. Click the link on triangular numbers to find other explanations.)\n",
"\n",
"Therefore, the first summation equals $3(333(334)/2) = 166833$, and the second equals $5(199(200)/2) = 99500$.\n",
"\n",
"However, we can't just sum these two values. That's because any number that's a multiple of both 3 and 5 is included in *both* of these summations, so if we add the results up, we'll count those values twice! To adjust for this, we'll need to subtract the sum of those numbers from our total so they are only counted once. Luckily, a number is a multiple of both 3 and 5 if and only if it is a multiple of 15, so we can apply our triangular number formula from above.\n",
"$$15 + 30 + 45 + \\cdots + 990 = 15(1 + 2 + 3 + \\cdots + 66) = 15(66(67)/2) = 33165$$\n",
"\n",
"Our answer is then\n",
"$$(3 + 6 + 9 + \\cdots + 999) + (5 + 10 + 15 + \\cdots + 995) - (15 + 30 + 45 + \\cdots + 990) = 233168$$\n",
"\n",
"## Relevant sequences\n",
"* Triangular numbers: [A000217](https://oeis.org/A000217)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

91
notebooks/problem0002.ipynb Normal file
View File

@@ -0,0 +1,91 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "1899e933",
"metadata": {},
"source": [
"# [Even Fibonacci Numbers](https://projecteuler.net/problem=2)\n",
"\n",
"There's a lot to learn about the [Fibonacci numbers](https://en.wikipedia.org/wiki/Fibonacci_sequence), and this will not be the last time we see them in a Project Euler problem, but for now this problem keeps it relatively simple.\n",
"\n",
"Since this is a SageMath notebook, we'll use its functionality so we can keep our code simple like the problem."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "57b68465",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"4613732"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"sum(n for n in fibonacci_xrange(1, 4000000) if n % 2 == 0)"
]
},
{
"cell_type": "markdown",
"id": "eef788da",
"metadata": {},
"source": [
"Even without SageMath, it's easy to implement our own generator of the sequence up to $n$."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "24a7c257",
"metadata": {},
"outputs": [],
"source": [
"def fib(n):\n",
" a = 0\n",
" b = 1\n",
" while a <= n:\n",
" yield a\n",
" a, b = b, a + b"
]
},
{
"cell_type": "markdown",
"id": "b3293a7d",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Fibonacci numbers: [A000045](https://oeis.org/A000045)\n",
"* Partial sums of even Fibonacci numbers: [A099919](https://oeis.org/A099919)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

129
notebooks/problem0003.ipynb Normal file
View File

@@ -0,0 +1,129 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "1a8d8609",
"metadata": {},
"source": [
"# [Largest Prime Factor](https://projecteuler.net/problem=3)\n",
"\n",
"The [fundamental theorem of arithmetic](https://en.wikipedia.org/wiki/Fundamental_theorem_of_arithmetic) says that every integer greater than 1 has a unique prime factorization.\n",
"\n",
"There are numerous factorization algorithms. Some of the easiest to implement are [trial division](https://en.wikipedia.org/wiki/Trial_division), [wheel factorization](https://en.wikipedia.org/wiki/Wheel_factorization) and [Pollard's rho algorithm](https://en.wikipedia.org/wiki/Pollard%27s_rho_algorithm).\n",
"\n",
"In general, as an integer gets larger, [factoring it gets harder](https://en.wikipedia.org/wiki/Integer_factorization). We don't have any efficient algorithm for factoring in general, but somewhat strangely, we don't know for certain that no such algorithm exists.\n",
"\n",
"Fortunately for us, this number is not too large for any of these algorithms, or for SageMath."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "d2f8884c",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"[71, 839, 1471, 6857]"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"prime_factors(600851475143)"
]
},
{
"cell_type": "markdown",
"id": "d4343a42",
"metadata": {},
"source": [
"If we wanted to implement our own factorization algorithm, though, Pollard's rho is relatively easy to implement and more efficient (on average) than trial division. It's also got some interesting math at its foundation, so let's take a look at it.\n",
"\n",
"## Pollard's rho algorithm\n",
"Let $n=pq$ be a positive integer where $p$ and $q$ are non-trivial factors of $n$ (i.e. not 1 or $n$), and let $\\{x_i\\}$ be a sequence such that $x_0$ is some random initial value and\n",
"$$x_{i+1} = (x_i^2 - 1) \\bmod{n}$$\n",
"Since the number of values this sequence can take on is finite, and each term is determined entirely by the previous term, it must ultimately become periodic, although there may be some initial terms before we reach a term that will repeat (refer to [this image](https://en.wikipedia.org/wiki/File:Pollard_rho_cycle.svg) for a visual aid - the fact this diagram looks like the Greek letter $\\rho$ is where the algorithm gets its name).\n",
"\n",
"How many terms will we see before the sequence repeats? If we assume the sequence is random, the [birthday paradox](https://en.wikipedia.org/wiki/Birthday_problem) implies that it will take an average of $O(\\sqrt{n})$ iterations before it repeats (this is using [Big O notation](https://en.wikipedia.org/wiki/Big_O_notation)). Of course, the sequence is not truly random, but it *is* [pseudorandom](https://en.wikipedia.org/wiki/Pseudorandomness), and experimental results seem to corroborate this analysis despite the somewhat faulty assumption.\n",
"\n",
"Now consider the related sequence $\\{x_i \\bmod p\\}$. Obviously, we can't calculate this sequence directly since we don't know the value of $p$, but similarly to $\\{x_i\\}$, this sequence is eventually periodic, with an average of $O(\\sqrt{p})$ iterations before a repetition. In other words, this sequence is expected to repeat earlier than $\\{x_i\\}$.\n",
"\n",
"When $\\{x_i \\bmod p\\}$ cycles, it implies that $x_j \\equiv x_k \\pmod{p}$ for some $j \\neq k$. This implies that $p | (x_j - x_k)$. Since $p$ also divides $n$, this means that $\\gcd(x_j - x_k, n)$ will be some multiple of $p$.\n",
"\n",
"This all means that we simply need to iterate through terms in $\\{x_i\\}$ until we find $x_j$ and $x_k$ such that $\\gcd(x_j - x_k, n)$ equals something other than 1. When it does, by definition this will be a factor of $n$. However, this factor may be $n$ itself, which is trivial. In this case the algorithm fails, but it can be attempted again with a different value for $x_0$.\n",
"\n",
"To facilitate finding $x_j$ and $x_k$, [Floyd's cycle-finding algorithm](https://en.wikipedia.org/wiki/Cycle_detection) is used. Briefly, two variables are used to iterate through $\\{x_i\\}$ at different rates. Because of this, the two variables will never have the same index, but eventually $x_i \\equiv x_j \\pmod{p}$, which we will detect with the GCD."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "d12993bf",
"metadata": {},
"outputs": [],
"source": [
"from itertools import count\n",
"\n",
"def factorization(n):\n",
" x = randint(0, n-1)\n",
" y = x\n",
" k = 2\n",
" for i in count(1):\n",
" x = (x ** 2 - 1) % n\n",
" d = gcd(y - x, n)\n",
" \n",
" assert d != n\n",
" \n",
" if d != 1:\n",
" return d\n",
" \n",
" if i == k:\n",
" y = x\n",
" k *= 2"
]
},
{
"cell_type": "markdown",
"id": "569b615c",
"metadata": {},
"source": [
"Note that this algorithm gives us a factor, but not necessarily a *prime* factor. However, we could just repeatedly apply the algorithm until we only have prime factors left. But how do we know when a factor is prime? We'll need a [primality test](https://en.wikipedia.org/wiki/Primality_test) like the Miller-Rabin algorithm - see [problem 7](https://projecteuler.net/problem=7) for a discussion of that algorithm."
]
},
{
"cell_type": "markdown",
"id": "12233560",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Greatest prime factors: [A006530](https://oeis.org/A006530)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

60
notebooks/problem0004.ipynb Normal file
View File

@@ -0,0 +1,60 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "9c7a281d",
"metadata": {},
"source": [
"# [Largest Palindrome Product](https://projecteuler.net/problem=4)\n",
"\n",
"Our search space is only $\\binom{900}{2} = 404550$ 3-digit pairs - a lot to check by hand but peanuts to a modern computer (we're using notation for [combinations](https://en.wikipedia.org/wiki/Combination) if you're unfamiliar)."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "e5fd66f6",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"906609"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"from itertools import combinations\n",
"\n",
"is_palindrome = lambda x: str(x) == str(x)[::-1]\n",
"three_digit_pairs = combinations(range(100, 1000), 2)\n",
"max(x*y for (x,y) in three_digit_pairs if is_palindrome(x*y))"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

74
notebooks/problem0005.ipynb Normal file
View File

@@ -0,0 +1,74 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "088fe896",
"metadata": {},
"source": [
"# [Smallest Multiple](https://projecteuler.net/problem=5)\n",
"\n",
"Another problem we can solve by hand! We're looking for the [least common multiple](https://en.wikipedia.org/wiki/Least_common_multiple), or LCM, of all the numbers from 1 to 20. Obviously,\n",
"$$20! = 20 \\times 19 \\times 18 \\times \\cdots \\times 1 = 2432902008176640000$$\n",
"is *a* common multiple, but the LCM is indeed smaller than this.\n",
"\n",
"One way to find the LCM is to first determine the prime factorizations of these numbers.\n",
"* $1 = 1$\n",
"* $2 = 2$\n",
"* $3 = 3$\n",
"* $4 = 2^2$\n",
"* $5 = 5$\n",
"* $6 = 2 \\times 3$\n",
"* $7 = 7$\n",
"* $8 = 2^3$\n",
"* $9 = 3^2$\n",
"* $10 = 2 \\times 5$\n",
"* $11 = 11$\n",
"* $12 = 2^2 \\times 3$\n",
"* $13 = 13$\n",
"* $14 = 2 \\times 7$\n",
"* $15 = 3 \\times 5$\n",
"* $16 = 2^4$\n",
"* $17 = 17$\n",
"* $18 = 2 \\times 3^2$\n",
"* $19 = 19$\n",
"* $20 = 2^2 \\times 5$\n",
"\n",
"Now, for each prime factor that appears in these factorizations, we'll take note of the largest power they're raised to. For instance, $2^4$ (the prime factorization of 16) is the largest power of 2 appearing.\n",
"\n",
"The LCM is the product of all these \"largest prime powers.\"\n",
"$$2^4 \\times 3^2 \\times 5 \\times 7 \\times 11 \\times 13 \\times 17 \\times 19 = 232792560$$\n",
"How do we know this is the LCM? Well, it's trivial to see it is *a* common multiple - just divide by each of the above numbers to confirm. Furthermore, we know it's the *least* common multiple because if we divided by any of these prime factors to get a smaller result, it would no longer be divisible by one or more of our numbers. For instance, if we divide by 3, our result is no longer divisible by $9 = 3^2$ or $18 = 2 \\times 3^2$.\n",
"\n",
"## Alternate method\n",
"I personally find the prime factorization method above the easiest way to solve this specific problem by hand. However, if there were a lot more numbers or the numbers were larger, I wouldn't want to do all those factorizations. Fortunately, there is another method that doesn't involve factoring at all, and can be implemented efficiently on a computer.\n",
"\n",
"For two numbers $m$ and $n$,\n",
"$$\\mathrm{lcm}(m,n) = \\frac{mn}{\\gcd(m,n)}$$\n",
"where $\\gcd$ is the [greatest common divisor](https://en.wikipedia.org/wiki/Greatest_common_divisor). The [Euclidean algorithm](https://en.wikipedia.org/wiki/Euclidean_algorithm) gives an efficient method for calculating the GCD.\n",
"\n",
"This formula is only defined for two numbers, but if we want to find an LCM of a set of three or more numbers, we can simply find the LCM of any two numbers in the set, then proceed to find the LCM of that value and another number from the set, repeating until we have used each number. For example, to find the LCM of 5, 8, and 14, you can first find the LCM of 5 and 8 using the above formula (40), then find the LCM of 40 and 14 (280)."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

49
notebooks/problem0006.ipynb Normal file
View File

@@ -0,0 +1,49 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "a39ba505",
"metadata": {},
"source": [
"# [Sum Square Difference](https://projecteuler.net/problem=6)\n",
"\n",
"In [problem 1](https://projecteuler.net/problem=1), we applied the following formula for [triangular numbers](https://en.wikipedia.org/wiki/Triangular_number):\n",
"$$1 + 2 + 3 + \\cdots + n = \\frac{n(n+1)}{2}$$\n",
"We can apply it again here and determine that\n",
"$$(1 + 2 + 3 + \\cdots + 100)^2 = \\left(\\frac{100(101)}{2}\\right)^2 = 25502500$$\n",
"\n",
"A similar formula exists for computing sums of squares, also called the [square pyramidal numbers](https://en.wikipedia.org/wiki/Square_pyramidal_number):\n",
"$$1^2 + 2^2 + 3^2 + \\cdots + n^2 = \\frac{n(n+1)(2n+1)}{6}$$\n",
"(In fact, [Faulhaber's formula](https://en.wikipedia.org/wiki/Faulhaber%27s_formula) gives a formula for the sum of $k$th powers, but we obviously only need the cases $k=1$ and $k=2$ for this problem.) Consequently,\n",
"$$1^2 + 2^2 + 3^2 + \\cdots + 100^2 = \\frac{100(101)(201)}{6} = 338350$$\n",
"\n",
"Therefore, the difference is $25502500 - 338350 = 25164150$.\n",
"\n",
"## Relevant sequences\n",
"* Triangular numbers: [A000217](https://oeis.org/A000217)\n",
"* Square pyramidal numbers: [A000330](https://oeis.org/A000330)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

188
notebooks/problem0007.ipynb Normal file
View File

@@ -0,0 +1,188 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "87fa1305",
"metadata": {},
"source": [
"# [100001st Prime](https://projecteuler.net/problem=7)\n",
"\n",
"Unlike integer factorization, we *can* efficiently determine whether a number is prime. One of the simplest and fastest methods is the [Miller-Rabin primality test](https://en.wikipedia.org/wiki/Miller%E2%80%93Rabin_primality_test), which is a probabilistic test.\n",
"\n",
"Of course, SageMath has tools to let us compute the 10001st prime without implementing a primality test ourselves! "
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "3a27f454",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"104743"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Primes.unrank is 0-indexed, so the 10001st prime is at index 10000\n",
"Primes().unrank(10000)"
]
},
{
"cell_type": "markdown",
"id": "6f644292",
"metadata": {},
"source": [
"Let's look into how the Miller-Rabin test works anyway, though.\n",
"\n",
"## Fermat's little theorem\n",
"First, a little bit of elementary number theory. [Fermat's little theorem](https://en.wikipedia.org/wiki/Fermat%27s_little_theorem) says that if $p$ is a prime number, then for any integer $a$,\n",
"$$a^p \\equiv a \\pmod{p}$$\n",
"It follows from the [contrapositive](https://en.wikipedia.org/wiki/Contraposition) that if there exists an integer $a$ such that $a^p \\not\\equiv a \\pmod{p}$, then $p$ is not prime. We then call $a$ a *witness* that proves $p$ is composite.\n",
"\n",
"Consider the fact that $3^{6} \\equiv 3 \\pmod{6}$, but 6 is composite. In cases like this, we say that $p$ is a *base $a$ Fermat pseudoprime*, e.g. 6 is a base 3 Fermat pseudoprime. 6 is also a base 4 Fermat pseudoprime, but 2 and 5 both serve as witnesses to 6 being composite (we need only consider bases modulo $p$).\n",
"\n",
"This may lead us to believe that we can just check every base $a$, and if there are no witnesses, then $p$ is prime. However, (along with this method being inefficient for large $p$), there are numbers that are Fermat pseudoprimes for *every* base, called [Carmichael numbers](https://en.wikipedia.org/wiki/Carmichael_number). Therefore, the [converse](https://en.wikipedia.org/wiki/Converse_(logic)) of Fermat's little theorem does not hold.\n",
"\n",
"Nevertheless, Carmichael numbers are relatively rare, so if you pick a random number to test and a few random bases, and find no witnesses, you can be pretty sure the number is in fact prime.\n",
"\n",
"## Miller-Rabin test\n",
"\n",
"The Miller-Rabin test employs another useful fact about prime numbers: if $p$ is prime, then the only values of $x$ such that \n",
"$$x^2 \\equiv 1 \\pmod{p}$$\n",
"are $x \\equiv \\pm 1 \\pmod{p}$. (*Proof: if $x^2 \\equiv 1 \\pmod{p}$, then $p \\vert (x+1)(x-1)$. Then by [Euclid's lemma](https://en.wikipedia.org/wiki/Euclid%27s_lemma), $p$ must divide one of $x+1$ or $x-1$, implying $x \\equiv \\pm 1 \\pmod{p}$.*)\n",
"\n",
"Similarly to Fermat's little theorem, the converse of this statement doesn't hold. For example, $x^2 \\equiv 1 \\pmod{9}$ only has $x=1$ and $x=-1$ as solutions, but 9 is composite.\n",
"\n",
"We now have two questions to lead an investigation into the primality of a number $n$:\n",
"* Is there an integer $a$ such that $a^n \\not\\equiv a \\pmod{n}$?\n",
"* Does 1 have a non-trivial square root modulo $n$?\n",
"\n",
"If the answer to either of these is yes, then $n$ is *definitely* composite. If both answers are no, this is strong evidence that $n$ is prime, *but not proof*, since we have already demonstrated that both of these conditions can have false positives. A composite number that is erroneously suggested to be prime by both of these conditions is called a *base $a$ strong pseudoprime*. Fortunately, no number (including Carmichael numbers) is a base $a$ strong pseudoprime for all bases $a$. In fact, Rabin (one of the namesakes of the algorithm) proved that if $n$ is composite, at least 3/4 of the numbers 1 to $n-1$ are witnesses.\n",
"\n",
"Let's finally look at the actual test. Suppose we want to test if an odd number $n$ is prime (testing whether an even number is prime is trivial). There must exist $t \\geq 1$ and odd $u$ such that $n - 1 = 2^t u$."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "720f7e8c",
"metadata": {},
"outputs": [],
"source": [
"def factor_twos(n):\n",
" t = 0\n",
" u = n - 1\n",
" while u % 2 == 0:\n",
" u //= 2\n",
" t += 1\n",
" \n",
" return t, u"
]
},
{
"cell_type": "markdown",
"id": "980bd2cf",
"metadata": {},
"source": [
"It follows that\n",
"$$a^{n-1} \\equiv a^{2^t u} \\pmod{n}$$\n",
"for any base $a$, so we'll randomly select a base $1 \\leq a \\leq n-1$.\n",
"\n",
"The test proceeds by repeatedly squaring $a^u$ modulo $n$ until we reach $a^{n-1}$. If at any point in this squaring, we find that $a^{2^{i-1} u}$ does not equal 1 or -1, but the following squaring gives $a^{2^i u} = 1$, then we have found a non-trivial square root of 1 modulo $n$. Therefore, $n$ must be composite.\n",
"\n",
"After $t$ squarings, we will reach $a^{n-1}$. If this equals anything other than 1, then by the contrapositive of Fermat's last theorem, this also tells us that $n$ is composite.\n",
"\n",
"If, after $t$ squarings, we have not found a non-trivial square root of 1 modulo $n$ and $a^{n-1} \\equiv 1 \\pmod{n}$, then either $n$ is prime, or is a base $a$ strong pseudoprime."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "123d1068",
"metadata": {},
"outputs": [],
"source": [
"def witness(a, n):\n",
" t, u = factor_twos(n)\n",
" x = pow(a, u, n)\n",
" for _ in range(0, t):\n",
" x, y = pow(x, 2, n), x\n",
" if x == 1 and y != 1 and y != n - 1:\n",
" return True\n",
" \n",
" return x != 1"
]
},
{
"cell_type": "markdown",
"id": "f8454603",
"metadata": {},
"source": [
"You can reduce the chance of a false report of primality by repeating the procedure with different bases. If we test several bases and do not find any witnesses, it is exceedingly unlikely that $n$ is composite."
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "e71ba6ec",
"metadata": {},
"outputs": [],
"source": [
"from random import randint\n",
"\n",
"def is_prime(n, k=10):\n",
" if n < 2:\n",
" return False\n",
" elif n == 2:\n",
" return True\n",
" elif n % 2 == 0:\n",
" return False\n",
" \n",
" for _ in range(0, k):\n",
" a = randint(1, n - 1)\n",
" if witness(a, n):\n",
" return False\n",
" \n",
" return True"
]
},
{
"cell_type": "markdown",
"id": "ec384441",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Prime numbers: [A000040](https://oeis.org/A000040)\n",
"* Carmichael numbers: [A002997](https://oeis.org/A002997)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

83
notebooks/problem0008.ipynb Normal file
View File

@@ -0,0 +1,83 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "8744af43",
"metadata": {},
"source": [
"# [Largest Product in a Series](https://projecteuler.net/problem=8)\n",
"\n",
"First, let's load the number as a tuple."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "699aa447",
"metadata": {},
"outputs": [],
"source": [
"num = tuple(int(d) for d in \"7316717653133062491922511967442657474235534919493496983520312774506326239578318016984801869478851843858615607891129494954595017379583319528532088055111254069874715852386305071569329096329522744304355766896648950445244523161731856403098711121722383113622298934233803081353362766142828064444866452387493035890729629049156044077239071381051585930796086670172427121883998797908792274921901699720888093776657273330010533678812202354218097512545405947522435258490771167055601360483958644670632441572215539753697817977846174064955149290862569321978468622482839722413756570560574902614079729686524145351004748216637048440319989000889524345065854122758866688116427171479924442928230863465674813919123162824586178664583591245665294765456828489128831426076900422421902267105562632111110937054421750694165896040807198403850962455444362981230987879927244284909188845801561660979191338754992005240636899125607176060588611646710940507754100225698315520005593572972571636269561882670428252483600823257530420752963450\")"
]
},
{
"cell_type": "markdown",
"id": "8b53af56",
"metadata": {},
"source": [
"Now we'll just do a sliding window over the number to find the thirteen adjacent digits with the greatest product."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "6b760673",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(5, 5, 7, 6, 6, 8, 9, 6, 6, 4, 8, 9, 5)"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"subnums = (num[i:i+13] for i in range(0, len(num)-13))\n",
"max(subnums, key=prod)"
]
},
{
"cell_type": "markdown",
"id": "4148d160",
"metadata": {},
"source": [
"Their product is 23514624000."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

112
notebooks/problem0009.ipynb Normal file
View File

@@ -0,0 +1,112 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "df906ad5",
"metadata": {},
"source": [
"# [Special Pythagorean Triplet](https://projecteuler.net/problem=9)\n",
"\n",
"The key mathematical insight for this problem is that Euclid gave a parameterization for generating [Pythagorean triplets](https://en.wikipedia.org/wiki/Pythagorean_triple): given integers $m, n$ such that $m > n > 0$, a triplet is given by\n",
"$$\n",
"\\begin{align}\n",
"a &= m^2 - n^2 \\\\\n",
"b &= 2mn \\\\\n",
"c &= m^2 + n^2\n",
"\\end{align}\n",
"$$\n",
"Furthermore, if $m$ and $n$ are [coprime](https://en.wikipedia.org/wiki/Coprime_integers) (i.e. $\\gcd(m,n)=1$) and exactly one of $m$ and $n$ is even, they will generate a primitive triplet, i.e. a triplet such that $\\gcd(a,b,c)=1$."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "644fb73e",
"metadata": {},
"outputs": [],
"source": [
"from itertools import count\n",
"\n",
"def primitive_pythagorean_triplets():\n",
" for m in count(2):\n",
" for n in range(1, m):\n",
" if not ((m % 2) ^^ (n % 2)) or gcd(m, n) != 1:\n",
" continue\n",
"\n",
" a = m ** 2 - n ** 2\n",
" b = 2 * m * n\n",
" c = m ** 2 + n ** 2\n",
"\n",
" yield (a, b, c)"
]
},
{
"cell_type": "markdown",
"id": "db3c3419",
"metadata": {},
"source": [
"We can simply iterate through these primitive triplets, and for each one, scale them until we either find our special triplet or their sum is greater than 1000, at which point we'll try the next primitive triplet."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "7b2c3efa",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"(375, 200, 425)"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def find_special_triplet():\n",
" for (a, b, c) in primitive_pythagorean_triplets():\n",
" for k in count(1):\n",
" if k * (a + b + c) == 1000:\n",
" return (k * a, k * b, k * c)\n",
" \n",
" if k * (a + b + c) > 1000:\n",
" break\n",
" \n",
" \n",
"find_special_triplet()"
]
},
{
"cell_type": "markdown",
"id": "8809ed3b",
"metadata": {},
"source": [
"So our product is $abc=31875000$."
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

94
notebooks/problem0010.ipynb Normal file
View File

@@ -0,0 +1,94 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "ec6f67ee",
"metadata": {},
"source": [
"# [Summation of Primes](https://projecteuler.net/problem=10)\n",
"\n",
"Once again, SageMath comes with functions that do the hard part for us."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "127873d5",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"142913828922"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"sum(prime_range(2000000))"
]
},
{
"cell_type": "markdown",
"id": "46f100db",
"metadata": {},
"source": [
"If `prime_range` didn't exist, we may have attempted using a primality test (like in [problem 7](https://projecteuler.net/problem=7)) to test each number less than 2,000,000, but for this problem, there's an even simpler way than implementing a test: the [sieve of Eratosthenes](https://en.wikipedia.org/wiki/Sieve_of_Eratosthenes), which can even be carried out on pen and paper (although maybe not for primes up to 2,000,000).\n",
"\n",
"Given an upper limit $n$, the sieve determines all the primes less than or equal to $n$. It does this by first assuming every number from 2 to $n$ is prime. Then, since 2 is assumed to be prime, every multiple of 2 up to $n$ is marked as composite. After this, we skip to the next number still assumed prime (3 for the second iteration) and mark all its multiples as composite. Repeat for each successive \"assumed\" prime. At the end, any number not marked composite is prime. We can then add those primes up to get our answer."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "35ec90ad",
"metadata": {},
"outputs": [],
"source": [
"def sieve(n):\n",
" # leave out multiples of 2 to save memory\n",
" P = set(range(3, n + 1, 2))\n",
" \n",
" for i in range(3, isqrt(n) + 1, 2):\n",
" if i in P:\n",
" P -= set(range(i ** 2, n + 1, 2 * i))\n",
" \n",
" P.add(2)\n",
" return P"
]
},
{
"cell_type": "markdown",
"id": "a4a187cf",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Partial sums of primes: [A007504](https://oeis.org/A007504)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

103
notebooks/problem0011.ipynb Normal file
View File

@@ -0,0 +1,103 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "124fe12d",
"metadata": {},
"source": [
"# [Largest Product in a Grid](https://projecteuler.net/problem=11)\n",
"\n",
"Here's that whole grid in a NumPy array."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "7d93e5b9",
"metadata": {},
"outputs": [],
"source": [
"import numpy as np\n",
"grid = np.array([[08,02,22,97,38,15,00,40,00,75,04,05,07,78,52,12,50,77,91,08],\n",
" [49,49,99,40,17,81,18,57,60,87,17,40,98,43,69,48,04,56,62,00],\n",
" [81,49,31,73,55,79,14,29,93,71,40,67,53,88,30,03,49,13,36,65],\n",
" [52,70,95,23,04,60,11,42,69,24,68,56,01,32,56,71,37,02,36,91],\n",
" [22,31,16,71,51,67,63,89,41,92,36,54,22,40,40,28,66,33,13,80],\n",
" [24,47,32,60,99,03,45,02,44,75,33,53,78,36,84,20,35,17,12,50],\n",
" [32,98,81,28,64,23,67,10,26,38,40,67,59,54,70,66,18,38,64,70],\n",
" [67,26,20,68,02,62,12,20,95,63,94,39,63,08,40,91,66,49,94,21],\n",
" [24,55,58,05,66,73,99,26,97,17,78,78,96,83,14,88,34,89,63,72],\n",
" [21,36,23,09,75,00,76,44,20,45,35,14,00,61,33,97,34,31,33,95],\n",
" [78,17,53,28,22,75,31,67,15,94,03,80,04,62,16,14,09,53,56,92],\n",
" [16,39,05,42,96,35,31,47,55,58,88,24,00,17,54,24,36,29,85,57],\n",
" [86,56,00,48,35,71,89,07,05,44,44,37,44,60,21,58,51,54,17,58],\n",
" [19,80,81,68,05,94,47,69,28,73,92,13,86,52,17,77,04,89,55,40],\n",
" [04,52,08,83,97,35,99,16,07,97,57,32,16,26,26,79,33,27,98,66],\n",
" [88,36,68,87,57,62,20,72,03,46,33,67,46,55,12,32,63,93,53,69],\n",
" [04,42,16,73,38,25,39,11,24,94,72,18,08,46,29,32,40,62,76,36],\n",
" [20,69,36,41,72,30,23,88,34,62,99,69,82,67,59,85,74,04,36,16],\n",
" [20,73,35,29,78,31,90,01,74,31,49,71,48,86,81,16,23,57,05,54],\n",
" [01,70,54,71,83,51,54,69,16,92,33,48,61,43,52,01,89,19,67,48]])"
]
},
{
"cell_type": "markdown",
"id": "d4924c1d",
"metadata": {},
"source": [
"Now we'll just iterate through every 4x4 subgrid of this grid, and for each subgrid, compute the product of the rows, columns, and diagonals, and record the largest product we see."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "7ac0e53c",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"70600674\n"
]
}
],
"source": [
"maximum = 0\n",
"for i in range(0, 16):\n",
" for j in range(0, 16):\n",
" subgrid = grid[i:i+4, j:j+4]\n",
" \n",
" cols = subgrid.prod(axis=0)\n",
" rows = subgrid.prod(axis=1)\n",
" diag = np.diag(subgrid).prod()\n",
" anti_diag = np.diag(np.fliplr(subgrid)).prod()\n",
" \n",
" maximum = max(maximum, *cols, *rows, diag, anti_diag)\n",
"\n",
"print(maximum)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

86
notebooks/problem0012.ipynb Normal file
View File

@@ -0,0 +1,86 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "802ccaa1",
"metadata": {},
"source": [
"# [Highly Divisible Triangular Number](https://projecteuler.net/problem=12)\n",
"\n",
"Yet another problem with a straightforward solution thanks to SageMath's preinstalled functions. We'll also once again apply the closed formula for triangle numbers (see [problem 1](https://projecteuler.net/problem=1))."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "c6249c71",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"76576500\n"
]
}
],
"source": [
"for n in PositiveIntegers():\n",
" t = n * (n+1) / 2\n",
" if number_of_divisors(t) > 500:\n",
" print(t)\n",
" break"
]
},
{
"cell_type": "markdown",
"id": "b504d3c6",
"metadata": {},
"source": [
"If we wanted to implement our own [divisor counting function](https://en.wikipedia.org/wiki/Divisor_function) (denoted $\\sigma_0$), we could apply a factorization algorithm (as discussed for [problem 3](https://projecteuler.net/problem=3)) to compute each triangle number's prime factorization\n",
"$$n = 2^{r_1} 3^{r_2} 5^{r_3} 7^{r_4} \\cdots$$\n",
"(Note that every number can be written in this form, but some exponents will be 0 if the corresponding prime number is not a divisor of $n$.)\n",
"From there, the number of divisors is simply\n",
"$$\\sigma(2^{r_1} 3^{r_2} 5^{r_3} 7^{r_4} \\cdots) = (r_1 + 1)(r_2 + 1)(r_3 + 1)\\cdots$$\n",
"For example, $12 = 2^2 \\times 3$, so $\\sigma_0(12) = (2 + 1)(1 + 1) = 6$ (the factors are 1, 2, 3, 4, 6, and 12).\n",
"\n",
"If this formula doesn't make sense to you, try thinking about it like this: a divisor of $n$ will have the same construction $2^{s_1} 3^{s_2} 5^{s_3} 7^{s_4} \\cdots$, but each exponent must be less than or equal to its corresponding exponent in the factorization of $n$ (if it's larger, it won't divide $n$). This means $s_1$ can be any value from $0,1,2,\\ldots,r_1$. There are therefore $r_1 + 1$ options for the value of $s_1$, $r_2 + 1$ options for the value of $s_2$, and so on. Each choice is independent of the others, so we can multiply the number of options for each exponent together to get the total number of options, which gets us the total number of divisors.\n",
"\n",
"Continuing with our previous example, to construct a divisor of $12=2^2 \\times 3$, our options for the power of 2 are 0, 1, or 2, and our options for the power of 3 are 0 or 1. This gives us 6 options altogether for the two exponents. Sure enough, here are the 6 divisors of 12, each composed of one of $2^0$, $2^1$, or $2^2$, and one of $3^0$ or $3^1$:\n",
"\n",
"* $1=1$\n",
"* $2=2$\n",
"* $3=3$\n",
"* $4=2^2$\n",
"* $6=2 \\times 3$\n",
"* $12=2^2 \\times 3$\n",
"\n",
"## Relevant sequences\n",
"* Number of divisors: [A000005](https://oeis.org/A000005)\n",
"* Triangular numbers: [A000217](https://oeis.org/A000217)\n",
"* Number of divisors of triangular numbers: [A063440](https://oeis.org/A063440)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

173
notebooks/problem0013.ipynb Normal file
View File

@@ -0,0 +1,173 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "6c851d8c",
"metadata": {},
"source": [
"# [Large Sum](https://projecteuler.net/problem=13)\n",
"\n",
"Nothing particularly interesting here. Here's all the numbers in a list."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "0bed92ce",
"metadata": {},
"outputs": [],
"source": [
" nums = [37107287533902102798797998220837590246510135740250,\n",
" 46376937677490009712648124896970078050417018260538,\n",
" 74324986199524741059474233309513058123726617309629,\n",
" 91942213363574161572522430563301811072406154908250,\n",
" 23067588207539346171171980310421047513778063246676,\n",
" 89261670696623633820136378418383684178734361726757,\n",
" 28112879812849979408065481931592621691275889832738,\n",
" 44274228917432520321923589422876796487670272189318,\n",
" 47451445736001306439091167216856844588711603153276,\n",
" 70386486105843025439939619828917593665686757934951,\n",
" 62176457141856560629502157223196586755079324193331,\n",
" 64906352462741904929101432445813822663347944758178,\n",
" 92575867718337217661963751590579239728245598838407,\n",
" 58203565325359399008402633568948830189458628227828,\n",
" 80181199384826282014278194139940567587151170094390,\n",
" 35398664372827112653829987240784473053190104293586,\n",
" 86515506006295864861532075273371959191420517255829,\n",
" 71693888707715466499115593487603532921714970056938,\n",
" 54370070576826684624621495650076471787294438377604,\n",
" 53282654108756828443191190634694037855217779295145,\n",
" 36123272525000296071075082563815656710885258350721,\n",
" 45876576172410976447339110607218265236877223636045,\n",
" 17423706905851860660448207621209813287860733969412,\n",
" 81142660418086830619328460811191061556940512689692,\n",
" 51934325451728388641918047049293215058642563049483,\n",
" 62467221648435076201727918039944693004732956340691,\n",
" 15732444386908125794514089057706229429197107928209,\n",
" 55037687525678773091862540744969844508330393682126,\n",
" 18336384825330154686196124348767681297534375946515,\n",
" 80386287592878490201521685554828717201219257766954,\n",
" 78182833757993103614740356856449095527097864797581,\n",
" 16726320100436897842553539920931837441497806860984,\n",
" 48403098129077791799088218795327364475675590848030,\n",
" 87086987551392711854517078544161852424320693150332,\n",
" 59959406895756536782107074926966537676326235447210,\n",
" 69793950679652694742597709739166693763042633987085,\n",
" 41052684708299085211399427365734116182760315001271,\n",
" 65378607361501080857009149939512557028198746004375,\n",
" 35829035317434717326932123578154982629742552737307,\n",
" 94953759765105305946966067683156574377167401875275,\n",
" 88902802571733229619176668713819931811048770190271,\n",
" 25267680276078003013678680992525463401061632866526,\n",
" 36270218540497705585629946580636237993140746255962,\n",
" 24074486908231174977792365466257246923322810917141,\n",
" 91430288197103288597806669760892938638285025333403,\n",
" 34413065578016127815921815005561868836468420090470,\n",
" 23053081172816430487623791969842487255036638784583,\n",
" 11487696932154902810424020138335124462181441773470,\n",
" 63783299490636259666498587618221225225512486764533,\n",
" 67720186971698544312419572409913959008952310058822,\n",
" 95548255300263520781532296796249481641953868218774,\n",
" 76085327132285723110424803456124867697064507995236,\n",
" 37774242535411291684276865538926205024910326572967,\n",
" 23701913275725675285653248258265463092207058596522,\n",
" 29798860272258331913126375147341994889534765745501,\n",
" 18495701454879288984856827726077713721403798879715,\n",
" 38298203783031473527721580348144513491373226651381,\n",
" 34829543829199918180278916522431027392251122869539,\n",
" 40957953066405232632538044100059654939159879593635,\n",
" 29746152185502371307642255121183693803580388584903,\n",
" 41698116222072977186158236678424689157993532961922,\n",
" 62467957194401269043877107275048102390895523597457,\n",
" 23189706772547915061505504953922979530901129967519,\n",
" 86188088225875314529584099251203829009407770775672,\n",
" 11306739708304724483816533873502340845647058077308,\n",
" 82959174767140363198008187129011875491310547126581,\n",
" 97623331044818386269515456334926366572897563400500,\n",
" 42846280183517070527831839425882145521227251250327,\n",
" 55121603546981200581762165212827652751691296897789,\n",
" 32238195734329339946437501907836945765883352399886,\n",
" 75506164965184775180738168837861091527357929701337,\n",
" 62177842752192623401942399639168044983993173312731,\n",
" 32924185707147349566916674687634660915035914677504,\n",
" 99518671430235219628894890102423325116913619626622,\n",
" 73267460800591547471830798392868535206946944540724,\n",
" 76841822524674417161514036427982273348055556214818,\n",
" 97142617910342598647204516893989422179826088076852,\n",
" 87783646182799346313767754307809363333018982642090,\n",
" 10848802521674670883215120185883543223812876952786,\n",
" 71329612474782464538636993009049310363619763878039,\n",
" 62184073572399794223406235393808339651327408011116,\n",
" 66627891981488087797941876876144230030984490851411,\n",
" 60661826293682836764744779239180335110989069790714,\n",
" 85786944089552990653640447425576083659976645795096,\n",
" 66024396409905389607120198219976047599490197230297,\n",
" 64913982680032973156037120041377903785566085089252,\n",
" 16730939319872750275468906903707539413042652315011,\n",
" 94809377245048795150954100921645863754710598436791,\n",
" 78639167021187492431995700641917969777599028300699,\n",
" 15368713711936614952811305876380278410754449733078,\n",
" 40789923115535562561142322423255033685442488917353,\n",
" 44889911501440648020369068063960672322193204149535,\n",
" 41503128880339536053299340368006977710650566631954,\n",
" 81234880673210146739058568557934581403627822703280,\n",
" 82616570773948327592232845941706525094512325230608,\n",
" 22918802058777319719839450180888072429661980811197,\n",
" 77158542502016545090413245809786882778948721859617,\n",
" 72107838435069186155435662884062257473692284509516,\n",
" 20849603980134001723930671666823555245252804609722,\n",
" 53503534226472524250874054075591789781264330331690]"
]
},
{
"cell_type": "markdown",
"id": "745ec1fc",
"metadata": {},
"source": [
"Now we'll just sum them up and extract the first 10 digits."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "4219a904",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"5537376230"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"int(str(sum(nums))[:10])"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

96
notebooks/problem0014.ipynb Normal file
View File

@@ -0,0 +1,96 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "43f43a42",
"metadata": {},
"source": [
"# [Longest Collatz Sequence](https://projecteuler.net/problem=14)\n",
"\n",
"The [Collatz conjecture](https://en.wikipedia.org/wiki/Collatz_conjecture) is a famous unsolved problem, and a classic example of a seemingly simple question that has proven very difficult, if not impossible, to answer.\n",
"\n",
"It's easy enough to *define* a [recursive](https://en.wikipedia.org/wiki/Recursion_(computer_science%29) function to calculate the chain length for a starting number $n$.\n",
"$$\n",
"f(n) = \\begin{cases}\n",
"1 & n = 1 \\\\\n",
"1+f(n/2) & n \\equiv 0 \\pmod{2} \\\\\n",
"1+f(3n+1) & n \\neq 1\\ \\text{and}\\ n \\equiv 1 \\pmod{2}\n",
"\\end{cases}\n",
"$$\n",
"However, we want its *implementation* to be efficient. We can optimize greatly if we cache the outputs we compute (the computer science term for this is [memoization](https://en.wikipedia.org/wiki/Memoization)). For instance, if we store the fact that $f(4) = 3$ after we've computed it, when we later compute $f(8) = 1 + f(4)$, the program can use the stored value of 3 rather than recomputing $f(4)$. For large inputs, this will save us (or really the computer, I guess) from redoing work.\n",
"\n",
"Python has a nice decorator called [cache](https://docs.python.org/3/library/functools.html) that will automagically memoize our function."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "d1017296",
"metadata": {},
"outputs": [],
"source": [
"from functools import cache\n",
"\n",
"@cache\n",
"def collatz_length(n):\n",
" if n == 1:\n",
" return 1\n",
" elif n % 2 == 0:\n",
" return 1 + collatz_length(n // 2)\n",
" else:\n",
" return 1 + collatz_length(3 * n + 1)"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "fd145a5b",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"837799"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"max(range(1, 1000000), key=collatz_length)"
]
},
{
"cell_type": "markdown",
"id": "bc0666ce",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Collatz chain lengths: [A008908](https://oeis.org/A008908)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

116
notebooks/problem0015.ipynb Normal file
View File

@@ -0,0 +1,116 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "6cc372f3",
"metadata": {},
"source": [
"# [Lattice Paths](https://projecteuler.net/problem=15)\n",
"\n",
"Let $f(m, n)$ equal the number of routes in an $m \\times n$ grid. Intuitively, if you start at the top left corner, you immediately have two choices - you can go down or right. If you go down, there are $f(m-1,n)$ routes in the $(m - 1) \\times n$ subgrid; if you go right, there are $f(m,n-1)$ routes in the $m \\times (n-1)$ subgrid. Therefore\n",
"$$f(m,n) = f(m-1,n) + f(m,n-1)$$\n",
"\n",
"There are two trivial base cases: $f(m,0) = 1$ for any $m$ and $f(0,n) = 1$ for any $n$. From these facts, you can write a simple [memoized](https://en.wikipedia.org/wiki/Memoization) recursive function to solve this problem."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "9878534a",
"metadata": {},
"outputs": [],
"source": [
"from functools import cache\n",
"\n",
"@cache\n",
"def f(m, n):\n",
" if m == 0:\n",
" return 1\n",
" if n == 0:\n",
" return 1\n",
" \n",
" return f(m - 1, n) + f(m, n - 1)"
]
},
{
"cell_type": "markdown",
"id": "5beab0ed",
"metadata": {},
"source": [
"Then the answer is given by"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "f21257bc",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"137846528820"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"f(20,20)"
]
},
{
"cell_type": "markdown",
"id": "3561530b",
"metadata": {},
"source": [
"However, we can do better and give a non-recursive formula for $f$. It turns out, for an $m \\times n$ grid, there are $\\binom{m+n}{m}$ possible routes. Here's an outline of how you can prove this [inductively](https://en.wikipedia.org/wiki/Mathematical_induction):\n",
"\n",
"1. Prove $f(m, 0) = \\binom{m+0}{m} = 1$ for all $m$\n",
" 1. Prove $f(0, 0) = \\binom{0 + 0}{0} = 1$\n",
" 2. Assuming $f(j, 0) = \\binom{j+0}{j} = 1$ for some $j$, prove $f(j+1,0) = \\binom{j+1+0}{j+1} = 1$\n",
"2. Assuming there exists some $k$ such that for all $m$, $f(m, k) = \\binom{m + k}{m}$, prove $f(m, k+1) = \\binom{m + k + 1}{m}$ for all $m$\n",
" 1. Prove $f(0, k+1) = \\binom{0+k+1}{0} = 1$\n",
" 2. Assuming $f(j, k+1) = \\binom{j+k+1}{j}$ for some $j$, prove $f(j+1,k+1) = \\binom{j+k+2}{j+1}$\n",
"\n",
"Note that there are *three* inductive proofs here - an overall induction on $m$ and nested inductions in both the base case and the induction step (yo dawg, I heard you like induction...).\n",
"\n",
"One useful identity for proving this, particularly in step 2, comes from [Pascal's triangle](https://en.wikipedia.org/wiki/Pascal%27s_triangle):\n",
"$$\\binom{n}{k} = \\binom{n-1}{k} + \\binom{n-1}{k-1}$$"
]
},
{
"cell_type": "markdown",
"id": "2c6e7d8c",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Central binomial coefficients: [A000984](https://oeis.org/A000984)\n",
"* General formula: [A046899](https://oeis.org/A046899)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

68
notebooks/problem0016.ipynb Normal file
View File

@@ -0,0 +1,68 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "78b44de2",
"metadata": {},
"source": [
"# [Power Digit Sum](https://projecteuler.net/problem=16)\n",
"\n",
"Since Python/SageMath support integers of arbitrary size, this problem is pretty straightforward. Just compute $2^{1000}$ directly and sum the digits up."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "d3a790fe",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"1366\n"
]
}
],
"source": [
"n = 2 ^ 1000\n",
"total = 0\n",
"while n != 0:\n",
" total += n % 10\n",
" n //= 10\n",
"\n",
"print(total)"
]
},
{
"cell_type": "markdown",
"id": "60ba69ac",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Sums of digits of $2^n$: [A001370](https://oeis.org/A001370)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

147
notebooks/problem0017.ipynb Normal file
View File

@@ -0,0 +1,147 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "fb5ffc43",
"metadata": {},
"source": [
"# [Number Letter Counts](https://projecteuler.net/problem=17)\n",
"\n",
"Time to think about the rules for [English numerals](https://en.wikipedia.org/wiki/English_numerals), which you probably know intuitively but have never had to explain to a computer.\n",
"\n",
"First, we'll hard code some base cases:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "ba4aa4f0",
"metadata": {},
"outputs": [],
"source": [
"numerals = {\n",
" 1: \"one\",\n",
" 2: \"two\",\n",
" 3: \"three\",\n",
" 4: \"four\",\n",
" 5: \"five\",\n",
" 6: \"six\",\n",
" 7: \"seven\",\n",
" 8: \"eight\",\n",
" 9: \"nine\",\n",
" 10: \"ten\",\n",
" 11: \"eleven\",\n",
" 12: \"twelve\",\n",
" 13: \"thirteen\",\n",
" 14: \"fourteen\",\n",
" 15: \"fifteen\",\n",
" 16: \"sixteen\",\n",
" 17: \"seventeen\",\n",
" 18: \"eighteen\",\n",
" 19: \"nineteen\",\n",
" 20: \"twenty\",\n",
" 30: \"thirty\",\n",
" 40: \"forty\",\n",
" 50: \"fifty\",\n",
" 60: \"sixty\",\n",
" 70: \"seventy\",\n",
" 80: \"eighty\",\n",
" 90: \"ninety\",\n",
" }"
]
},
{
"cell_type": "markdown",
"id": "8b11e55b",
"metadata": {},
"source": [
"Obviously, if a number is explicitly given in this dictionary, that's its English numeral. For any other two digit number, the numeral is just the concatenation of the numerals for the tens place and the ones place (e.g. 37 is \"thirty\" and \"seven\").\n",
"\n",
"For three digit numbers, we write the numeral for the hundreds place, followed by \"hundred\" (e.g. 300 is \"three hundred\"). If it's an exact multiple of 100 (e.g. 100, 200, etc.) then that's it. Otherwise, we stick on \"and\" and follow the rules for a two digit number for the remaining digits.\n",
"\n",
"For this problem, all that's left is a special case for 1,000. Here's our function:"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "18b8dacb",
"metadata": {},
"outputs": [],
"source": [
"def numeral(num):\n",
" if num == 1000:\n",
" return \"onethousand\"\n",
" elif num in numerals:\n",
" return numerals[num]\n",
" elif num % 100 == 0:\n",
" return numerals[num // 100] + \"hundred\"\n",
" elif num > 100:\n",
" return numerals[num // 100] + \"hundredand\" + numeral(num % 100)\n",
" \n",
" ones = num % 10\n",
" rest = num - ones\n",
" return numerals[rest] + numerals[ones]"
]
},
{
"cell_type": "markdown",
"id": "1e72e33b",
"metadata": {},
"source": [
"Now we can just use this to find all the numerals we need and sum their lengths."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "73b14cb1",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"21124"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"sum(len(numeral(n)) for n in range(1, 1001))"
]
},
{
"cell_type": "markdown",
"id": "417c208b",
"metadata": {},
"source": [
"## Relevent sequences\n",
"* Number of letters in British numeral: [A362123](https://oeis.org/A362123)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

179
notebooks/problem0018.ipynb Normal file
View File

@@ -0,0 +1,179 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "b0604090",
"metadata": {},
"source": [
"# [Maximum Path Sum I](https://projecteuler.net/problem=18)\n",
"\n",
"As the problem notes, we could brute force every single path through this triangle, since there aren't that many, relatively speaking. But we're going to get a tougher version in [problem 67](https://projecteuler.net/problem=67), so if we make the effort now, we can kill two birds with one stone.\n",
"\n",
"First, let's encode the triangle as a list of lists:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "c3db04f0",
"metadata": {},
"outputs": [],
"source": [
" triangle = [[75],\n",
" [95, 64],\n",
" [17, 47, 82],\n",
" [18, 35, 87, 10],\n",
" [20, 4, 82, 47, 65],\n",
" [19, 1, 23, 75, 3, 34],\n",
" [88, 2, 77, 73, 7, 63, 67],\n",
" [99, 65, 4, 28, 6, 16, 70, 92],\n",
" [41, 41, 26, 56, 83, 40, 80, 70, 33],\n",
" [41, 48, 72, 33, 47, 32, 37, 16, 94, 29],\n",
" [53, 71, 44, 65, 25, 43, 91, 52, 97, 51, 14],\n",
" [70, 11, 33, 28, 77, 73, 17, 78, 39, 68, 17, 57],\n",
" [91, 71, 52, 38, 17, 14, 91, 43, 58, 50, 27, 29, 48],\n",
" [63, 66, 4, 68, 89, 53, 67, 30, 73, 16, 69, 87, 40, 31],\n",
" [4, 62, 98, 27, 23, 9, 70, 98, 73, 93, 38, 53, 60, 4, 23]]"
]
},
{
"cell_type": "markdown",
"id": "51f3aa5a",
"metadata": {},
"source": [
"There are a couple of important things to note about this problem. First, looking at the triangle, try to see how there are two \"sub-triangles\" below the top value, 75 (it's easier to see this on the original page rather than in the list-of-lists). One sub-triangle has a top value of 95, then the next row is 17 and 47, and the next row is 18, 35, and 87, and so on.\n",
"```\n",
" 95\n",
" 17 47\n",
" 18 35 87\n",
" ⋰ ⋮ ⋱\n",
"```\n",
"\n",
"The other sub-triangle has a top value of 64, and the next row is 47 and 82, and so on.\n",
"```\n",
" 64\n",
" 47 82\n",
" 35 87 10\n",
" ⋰ ⋮ ⋱\n",
"```\n",
"\n",
"If we think a little abstractly, the maximum path sum (MPS) through *the whole triangle* must equal the sum of the topmost value (75) with the larger of the MPSs through *each sub-triangle*. This fact holds for the lower layers as well: the MPS for the sub-triangle with top value 95 depends on the MPSs for the sub-triangles with top values 17 and 47, and the MPS for the sub-triangle with top value 64 depends on the MPSs for the sub-triangles with top values 47 and 82, and so on and so on.\n",
"\n",
"Because of this, the problem is said to have [optimal substructure](https://en.wikipedia.org/wiki/Optimal_substructure). In short, we find the MPS through the whole triangle by solving smaller and smaller subproblems, i.e. MPSs through sub-triangles. Eventually we will reach a base case: one of the values in the bottom-most layer of the triangle, where the maximum path sum of that \"sub-triangle\" is just the value itself.\n",
"\n",
"There's another important thing to consider: in trying to find the MPS through the 95-sub-triangle, we have to find the MPS through the 47-sub-triangle. Similarly, in trying to find the MPS through the 64-sub-triangle, we *also* have to find the MPS through the 47-sub-triangle. This reoccurrence of subproblems happens more and more frequently in the lower layers, as well. In other words, there are [overlapping subproblems](https://en.wikipedia.org/wiki/Overlapping_subproblems).\n",
"\n",
"When a problem has optimal substructure and overlapping subproblems, we can use [dynamic programming](https://en.wikipedia.org/wiki/Dynamic_programming). We can employ either a top-down or bottom-up approach.\n",
"\n",
"## Top-down Method\n",
"This approach recusively finds the maximum path sums of each sub-triangle starting from the topmost value, similarly to what's described above. The `cache` decorator is used to memoize and avoid recomputing MPSs for any overlapping sub-triangles. This is crucial to the speed of a top-down approach; without memoization, we will recompute the solutions to overlapping subproblems way too often."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "778ba200",
"metadata": {},
"outputs": [],
"source": [
"from functools import cache\n",
"\n",
"@cache\n",
"def max_path_sum(tri):\n",
" if tri == tuple():\n",
" return 0\n",
"\n",
" left = tuple(row[:-1] for row in tri[1:])\n",
" right = tuple(row[1:] for row in tri[1:])\n",
" return tri[0][0] + max(max_path_sum(left), max_path_sum(right))"
]
},
{
"cell_type": "markdown",
"id": "dec8855a",
"metadata": {},
"source": [
"One slight wrinkle to this approach is that to use the `cache` decorator, the inputs to our function need to be [hashable](https://docs.python.org/3/glossary.html#term-hashable). Lists are not hashable, but tuples are, so we can convert our triangle to a tuple of tuples."
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "789460a8",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"1074"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"max_path_sum(tuple(tuple(row) for row in triangle))"
]
},
{
"cell_type": "markdown",
"id": "4e267302",
"metadata": {},
"source": [
"## Bottom-up Method\n",
"Another approach is to instead find the MPSs of the very smallest sub-triangles first, then use those MPSs to find the MPSs in the next layer up. This doesn't require memoization, since by starting with the smallest subproblems, we ensure they are only computed once."
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "23352679",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"1074"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def max_path_sum(tri):\n",
" for i in reversed(range(0, len(tri) - 1)):\n",
" for j in range(0, len(tri[i])):\n",
" tri[i][j] += max(tri[i+1][j], tri[i+1][j+1])\n",
" \n",
" return tri[0][0]\n",
"\n",
"max_path_sum(triangle)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

109
notebooks/problem0019.ipynb Normal file
View File

@@ -0,0 +1,109 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "819ff42e",
"metadata": {},
"source": [
"# [Counting Sundays](https://projecteuler.net/problem=19)\n",
"\n",
"First, here's a simple function to give us the number of days in any given month and year."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "0f09cadc",
"metadata": {},
"outputs": [],
"source": [
"def days_in_month(month, year):\n",
" key = {1: 31,\n",
" 2: 28,\n",
" 3: 31,\n",
" 4: 30,\n",
" 5: 31,\n",
" 6: 30,\n",
" 7: 31,\n",
" 8: 31,\n",
" 9: 30,\n",
" 10: 31,\n",
" 11: 30,\n",
" 12: 31}\n",
" \n",
" if year % 400 == 0 or (year % 4 == 0 and year % 100 != 0):\n",
" key[2] = 29\n",
" \n",
" return key[month]"
]
},
{
"cell_type": "markdown",
"id": "1537a949",
"metadata": {},
"source": [
"Let's talk briefly about the days of the week. If someone asked you on a Monday what day it will be in eight days, you could naively count up (one day from now is Tuesday, two days from now is Wednesday, etc.), but it's quicker to observe that since the days of the week repeat every seven days, seven days from now will also be a Monday, and therefore eight days later is a Tuesday. Likewise, if someone wanted to know what day is 17 days from now, you know 14 days from now is also a Monday, so 17 (14+3) days later would be Thursday (3 days after Monday).\n",
"\n",
"Now, let's reframe this concept mathematically. Let's assign a number to every day of the week: Sunday is 0, Monday is 1, Tuesday is 2, and so on. Using [modular arithmetic](https://en.wikipedia.org/wiki/Modular_arithmetic), we can add the number of days later to our day-number, then compute the remainder modulo 7 to find the day of our new date. Restating what we found above, if today is a Monday (1), we just add 17 (to get 18) and compute the remainder modulo 7 to see that 17 days from now will be Thursday (4). As a formula: if the day-number is $s$, the day it will be $n$ days later is $(s + n) \\bmod 7$.\n",
"\n",
"Knowing this, to find the day each month started on, we can just iterate through every month from 1900 through 2000, and each time add the length of the month to the day the last month started on, compute the modulus, and check if the result is 0, indicating that month started on a Sunday. (Note that the question gives us the starting day of 1900, but only asked about months from 1901 through 2000, so we'll take care to not add any months from 1900 that fell on a Sunday.)"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "1c8f38fa",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"171\n"
]
}
],
"source": [
"total = 0\n",
"day = 1\n",
"for year in range(1900, 2001):\n",
" for month in range(1, 13):\n",
" day += days_in_month(month, year)\n",
" day %= 7\n",
" if day == 0 and year >= 1901:\n",
" total += 1\n",
"\n",
"print(total)"
]
},
{
"cell_type": "markdown",
"id": "83b6e18e",
"metadata": {},
"source": [
"One last note, just for fun: did you know you can learn how to calculate what day of the week any date falls on [in your head](https://en.wikipedia.org/wiki/Doomsday_rule)?"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

65
notebooks/problem0020.ipynb Normal file
View File

@@ -0,0 +1,65 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "c1d94bd8",
"metadata": {},
"source": [
"# [Factorial Digit Sum](https://projecteuler.net/problem=20)\n",
"\n",
"Piece of cake for SageMath."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "894e7489",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"648"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"sum(map(int, str(factorial(100))))"
]
},
{
"cell_type": "markdown",
"id": "6dab9653",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Sum of digits of $n!$: [A004152](https://oeis.org/A004152)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

96
notebooks/problem0021.ipynb Normal file
View File

@@ -0,0 +1,96 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "dea82444",
"metadata": {},
"source": [
"# [Amicable Numbers](https://projecteuler.net/problem=21)\n",
"\n",
"The sum of the proper divisors of a number is called the [aliquot sum](https://en.wikipedia.org/wiki/Aliquot_sum).\n",
"\n",
"SageMath [provides the `sigma` function](https://doc.sagemath.org/html/en/reference/rings_standard/sage/arith/misc.html#sage.arith.misc.Sigma), which can compute the sum of the divisors of $n$. The aliquot sum is then just `sigma(n) - n`.\n",
"\n",
"If a number is equal to its own aliquot sum, it's called a [perfect number](https://en.wikipedia.org/wiki/Perfect_number), and we exclude those numbers from our total."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "144e1f54",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"31626\n"
]
}
],
"source": [
"aliquot_sum = lambda n: sigma(n) - n\n",
"\n",
"amicables = set()\n",
"for a in range(2, 10000):\n",
" if a in amicables:\n",
" continue\n",
" \n",
" b = aliquot_sum(a)\n",
" if a == b:\n",
" continue\n",
" \n",
" if a == aliquot_sum(b):\n",
" amicables.update({a, b})\n",
"\n",
"print(sum(amicables))"
]
},
{
"cell_type": "markdown",
"id": "ca7d5f36",
"metadata": {},
"source": [
"Funny enough, there's only five pairs of amicable numbers below 10,000. If you looked up [amicable numbers](https://en.wikipedia.org/wiki/Amicable_numbers), you may have stumbled on all the numbers you need to answer the problem!\n",
"\n",
"## Sum of divisors\n",
"Of course, you could implement your own [divisor sum function](https://en.wikipedia.org/wiki/Divisor_function). In [problem 12](https://projecteuler.net/problem=12), we implemented a divisor *counting* function, which is related.\n",
"\n",
"One important property of the divisor sum function $\\sigma(n)$ is that it is [multiplicative](https://en.wikipedia.org/wiki/Multiplicative_function). This means that if $a$ and $b$ are [coprime](https://en.wikipedia.org/wiki/Coprime_integers), $\\sigma(ab) = \\sigma(a)\\sigma(b)$. It follows that if $n = 2^{r_1}3^{r_2}5^{r_3}7^{r_4}\\cdots$, then\n",
"$$\\sigma(n) = \\sigma(2^{r_1})\\sigma(3^{r_2})\\sigma(5^{r_3})\\sigma(7^{r_4})\\cdots$$\n",
"\n",
"Furthermore, for a prime $p$, $\\sigma(p^k) = 1 + p + p^2 + \\cdots + p^k$. This expression is actually a partial sum of a [geometric series](https://en.wikipedia.org/wiki/Geometric_series), which has a closed formula:\n",
"$$1 + p + p^2 + \\cdots + p^k = \\frac{p^{k+1}-1}{p-1}$$\n",
"\n",
"Putting it all together, we can compute $\\sigma(n)$ as\n",
"$$\\sigma(n) = \\frac{2^{r_1+1}-1}{2-1} \\cdot \\frac{3^{r_2+1}-1}{3-1} \\cdot \\frac{5^{r_3+1}-1}{5-1} \\cdot \\frac{7^{r_4+1}-1}{7-1} \\cdot \\cdots$$\n",
"\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",
"## Relevent sequences\n",
"* Amicable numbers: [A063990](https://oeis.org/A063990)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

73
notebooks/problem0023.ipynb Normal file
View File

@@ -0,0 +1,73 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "22f94a3b",
"metadata": {},
"source": [
"# [Non-Abundant Sums](https://projecteuler.net/problem=23)\n",
"\n",
"Just like in [problem 21](https://projecteuler.net/problem=21), we'll define an `aliquot_sum` function and use that to find all the [abundant numbers](https://en.wikipedia.org/wiki/Abundant_number) below 28,124. Then we check every integer less than 28,124 to see if it's the sum of any two abundant numbers, and if it is, remove it from a set containing all those integers. Whatever's left in that set are the non-abundant sums."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "6747d2a3",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"4179871\n"
]
}
],
"source": [
"aliquot_sum = lambda n: sigma(n) - n\n",
"\n",
"abundant_numbers = {k for k in range(1, 28124) if aliquot_sum(k) > k}\n",
"\n",
"non_abundant_sums = set(range(1, 28124))\n",
"for n in range(1, 28124):\n",
" for m in abundant_numbers:\n",
" if n - m in abundant_numbers:\n",
" non_abundant_sums.discard(n)\n",
" break\n",
"\n",
"print(sum(non_abundant_sums))"
]
},
{
"cell_type": "markdown",
"id": "91a43bd1",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Numbers that are not the sum of two abundant numbers: [A048242](https://oeis.org/A048242)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

57
notebooks/problem0025.ipynb Normal file
View File

@@ -0,0 +1,57 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "b5da29ec",
"metadata": {},
"source": [
"# [1000-digit Fibonacci Number](https://projecteuler.net/problem=25)\n",
"\n",
"This problem can be solved with a scientific calculator!\n",
"\n",
"We're looking at [Fibonacci numbers](https://en.wikipedia.org/wiki/Fibonacci_sequence) again, after last seeing them in [problem 2](https://projecteuler.net/problem=2). Our goal is to find the first 1000 digit number in the sequence. Mathematically, we can state this is as the lowest $n$ such that\n",
"$$1+\\log{F_n} \\geq 1000$$\n",
"(Here, $\\log$ is the base 10 logarithm).\n",
"\n",
"To assist in finding $n$, we can employ Binet's formula:\n",
"$$F_n = \\frac{\\phi^n - (-\\phi)^{-n}}{\\sqrt{5}}$$\n",
"Here, $\\phi$ is the [golden ratio](https://en.wikipedia.org/wiki/Golden_ratio).\n",
"\n",
"Binet's formula is exact, but we can make our work a little easier by approximating.\n",
"$$F_n \\approx \\frac{\\phi^n}{\\sqrt{5}}$$\n",
"This approximation works since $(-\\phi)^{-n}$ approaches 0 as $n$ increases. This also means this approximation gets more accurate as $n$ increases, which is especially convenient since we're looking for a large $n$. Substituting this approximation into the above inequality, we have\n",
"$$1 + \\log{\\frac{\\phi^n}{\\sqrt{5}}} \\geq 1000$$\n",
"\n",
"Now with a little bit of algebra, we can solve for $n$.\n",
"$$n \\geq 999\\log_{\\phi}10 + \\log_{\\phi}\\sqrt{5}$$\n",
"If you're trying to solve this with a scientific calculator, you probably don't have a $\\log_{\\phi}$ button (*please let me know if you do*), but we can just use the [logarithmic change of base formula](https://en.wikipedia.org/wiki/List_of_logarithmic_identities). This ends up getting us\n",
"$$n \\geq 4781.85927\\ldots$$\n",
"Therefore, we want $n=4782$.\n",
"\n",
"## Relevant sequences\n",
"* Fibonacci numbers: [A000045](https://oeis.org/A000045)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

129
notebooks/problem0026.ipynb Normal file
View File

@@ -0,0 +1,129 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "b7b7784e",
"metadata": {},
"source": [
"# [Reciprocal Cycles](https://projecteuler.net/problem=26)\n",
"\n",
"Another easy problem with SageMath."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "20a3987a",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"983"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"max(range(2, 1000), key=lambda d: (1/d).period())"
]
},
{
"cell_type": "markdown",
"id": "6dc68724",
"metadata": {},
"source": [
"But let's talk about how we would write our own algorithm for calculating the period.\n",
"\n",
"If $d = 2^a 5^b n$ where $n$ is coprime to 2 and 5, then the [period](https://mathworld.wolfram.com/DecimalPeriod.html) of $\\frac{1}{d}$ is the [multiplicative order](https://en.wikipedia.org/wiki/Multiplicative_order) of 10 modulo $n$. This is the same as finding the smallest positive $k$ such that\n",
"$$10^k \\equiv 1 \\pmod{n}$$\n",
"(Wondering why this is called multiplicative order? It has to do with the mathematical concept of [groups](https://en.wikipedia.org/wiki/Group_(mathematics)), but you don't need to be familiar with them to apply the formula.)\n",
"\n",
"Based on the definition, we can easily write a function that computes multiplicative order that is efficient enough for this problem, but it's not very efficient in general. The computation is a special case of the [discrete logarithm](https://en.wikipedia.org/wiki/Discrete_logarithm), which has no known efficient algorithm in general."
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "d967be38",
"metadata": {},
"outputs": [],
"source": [
"def multiplicative_order(a, n):\n",
" if n == 1:\n",
" return 1\n",
" \n",
" g = 1\n",
" for k in range(1, n):\n",
" g *= a\n",
" g %= n\n",
" if g == 1:\n",
" return k"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "dc51eb0e",
"metadata": {},
"outputs": [],
"source": [
"def period(d):\n",
" while d % 2 == 0:\n",
" d //= 2\n",
" while d % 5 == 0:\n",
" d //= 5\n",
" \n",
" return multiplicative_order(10, d)"
]
},
{
"cell_type": "markdown",
"id": "14f53b42",
"metadata": {},
"source": [
"Personally, I don't feel it's immediately obvious *why* multiplicative order can be used to calculate these periods. Here's a deeper dive into why this works, if you're curious.\n",
"\n",
"## Explanation\n",
"As above, let $d = 2^a 5^b n$, where $n$ is coprime to 2 and 5, and consider the unit fraction $u=\\frac{1}{d}$. The decimal representation of this fraction (as with any fraction) may have a [fractional part](https://en.wikipedia.org/wiki/Fractional_part) with $q$ leading digits, followed by a repetend of $r$ digits.\n",
"\n",
"Let's make note of a trick you probably already know: multiplying a number by 10 gives the same result as if you just moved the original number's decimal point to the right one digit (e.g. $1.25 \\times 10 = 12.5$). This means that $10^q u$ has a decimal representation with only the repetend in its fractional part. If you were to then multiply it by $10^r$, the integer part would increase, but the fractional part would stay the same since it repeats every $r$ digits.\n",
"\n",
"The above is important for understanding the following: $10^q 10^r u - 10^q u$ *is an integer*. How do we know? Because by the logic above, the fractional part of both $10^q 10^r u$ and $10^q u$ only have the repetend, so they cancel out when we subtract, leaving us with an integer. It immediately follows that $d(10^q 10^r u - 10^q u) = 10^q 10^r - 10^q$ is also an integer that has $d$ as a factor. We can express this as\n",
"$$10^{q+r} \\equiv 10^q \\pmod{d}$$\n",
"Through the properties of [modular arithmetic](https://en.wikipedia.org/wiki/Modular_arithmetic), we can cancel common factors of 2 and 5 in $10^{q+r}$, $10^q$, and $d$ until we end up with\n",
"$$10^r \\equiv 1 \\pmod{n}$$\n",
"By definition, $r$ is the multiplicative order of 10 modulo $n$.\n",
"\n",
"As a concrete example of the above, consider $d = 2^4 \\times 5 \\times 63 = 5040$. The decimal representation of $u = \\frac{1}{d}$ is $0.0001(984126)$, where 984126 is the repetend. Therefore $q=4$ and $r=6$. Sure enough, $10^q 10^r u - 10^q u = 1984125$ is an integer; therefore, $10^{10} \\equiv 10^4 \\pmod{5040}$, and $10^6 \\equiv 1 \\pmod{63}$.\n",
"\n",
"## Relevant sequences\n",
"* Periods of reciprocals: [A007732](https://oeis.org/A007732)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

134
notebooks/problem0027.ipynb Normal file
View File

File diff suppressed because one or more lines are too long

50
notebooks/problem0028.ipynb Normal file
View File

@@ -0,0 +1,50 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "dee40daf",
"metadata": {},
"source": [
"# [Number Spiral Diagonals](https://projecteuler.net/problem=28)\n",
"\n",
"Whether you choose to program or just use pen and paper, there's a lot of different ways to tackle this problem. Here's one approach.\n",
"\n",
"If we have an $n \\times n$ spiral (note that $n$ must be odd!), it's pretty easy to get a formula for just the sum of the corners. The top right corner will always be $n^2$, and since it's an $n \\times n$ spiral (as mentioned one sentence ago), the top left corner will just be $n^2 - (n - 1)$. We can continue subtracting $n-1$ to get the values of the other corners. Adding these up, we get\n",
"$$f(n) = n^2 + (n^2 - (n-1)) + (n^2 - 2(n-1)) + (n^2 - 3(n-1)) = 4n^2 - 6n + 6$$\n",
"\n",
"Of course, we want the sums of the *diagonals*, not just the outermost corners. We can think about this as getting the sum of the corners of the $n \\times n$ spiral, plus the sum of the corners of the $(n-2) \\times (n-2)$ spiral one layer deeper, and so on until we reach the 1 at the center.\n",
"$$g(n) = f(n) + f(n-2) + f(n-4) + \\cdots + f(5) + f(3) + 1 = 1 + \\sum_{k=1}^{(n-1)/2} f(2k + 1) = 1 + \\sum_{k=1}^{(n-1)/2} (16k^2 + 4k + 4)$$\n",
"\n",
"Now we can apply [summation identities](https://en.wikipedia.org/wiki/Summation) (including the return of triangular numbers and square pyramidal numbers from [problem 6](https://projecteuler.net/problem=6)) to get a closed formula:\n",
"$$g(n) = \\frac{2}{3}n^3 + \\frac{1}{2}n^2 + \\frac{4}{3}n - \\frac{3}{2}$$\n",
"Plugging in 1001, we get our answer: $g(1001) = 669171001$.\n",
"\n",
"Side note: this practice of writing the natural numbers in a spiral, combined with marking the prime numbers, has been coined the [Ulam spiral](https://en.wikipedia.org/wiki/Ulam_spiral). Somewhat interestingly, lots of primes appear in vertical, horiztonal, and diagonal lines when laid out this way.\n",
"\n",
"## Relevant sequences\n",
"* Numbers on diagonals: [A200975](https://oeis.org/A200975)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

56
notebooks/problem0029.ipynb Normal file
View File

@@ -0,0 +1,56 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "3bbfb66e",
"metadata": {},
"source": [
"# [Distinct Powers](https://projecteuler.net/problem=29)\n",
"\n",
"Python offers a [built-in set type](https://docs.python.org/3/library/stdtypes.html) that will only keep track of distinct items for us. We also don't need to worry about overflows, since integers are unlimited precision."
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "9fa32dcf",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"9183"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"len({a ** b for a in range(2, 101) for b in range(2, 101)})"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

75
notebooks/problem0030.ipynb Normal file
View File

@@ -0,0 +1,75 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "cbbd0d3f",
"metadata": {},
"source": [
"# [Digit Fifth Powers](https://projecteuler.net/problem=30)\n",
"\n",
"These sorts of numbers are called [perfect digital invariants](https://en.wikipedia.org/wiki/Perfect_digital_invariant) (that page talks about the concept for arbitrary bases, but we obviously only need to worry about base 10).\n",
"\n",
"It's not difficult to implement the PDI function from the page linked above, and to check numbers one by one to see if they are a PDI, but how do we know when to stop searching? Well, let $f(n)$ be the fifth power PDI function, and consider the largest possible six digit number, 999999. The sum of the fifth powers of the digits is $f(999999) = 6 \\times 9^5 = 354294$. If we replaced any digits in 999999, they would have to be smaller than 9, so the power digit sum would be smaller than 354294. Therefore, for any $n \\leq 999999$, $f(n) \\leq 354294$, so for $n$ to *equal* $f(n)$, $n$ must be less than or equal to 354294, as well.\n",
"\n",
"(There's a formal proof on the page linked above that shows the upper bound can be reduced even further, but I think the logic is harder to follow than what's presented here, and using this bound solves the problem quickly enough.)"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "e8ab0179",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"443839\n"
]
}
],
"source": [
"def F(n, p):\n",
" total = 0\n",
" while n != 0:\n",
" total += (n % 10) ** p\n",
" n //= 10\n",
" \n",
" return total\n",
"\n",
"\n",
"print(sum(n for n in range(2, 354295) if n == F(n, 5)))"
]
},
{
"cell_type": "markdown",
"id": "f0754e9f",
"metadata": {},
"source": [
"## Relevant sequences\n",
"* Perfect digital invariants: [A252648](https://oeis.org/A252648)"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 9.5",
"language": "sage",
"name": "sagemath"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 3
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.11.2"
}
},
"nbformat": 4,
"nbformat_minor": 5
}