refer to relevant problems instead of re-implementing
This commit is contained in:
parent
275bfae648
commit
fc7c74c323
|
@ -36,7 +36,7 @@
|
||||||
"id": "396468b7",
|
"id": "396468b7",
|
||||||
"metadata": {},
|
"metadata": {},
|
||||||
"source": [
|
"source": [
|
||||||
"If we wanted to compute the convergents ourselves, we could first make a generator for the partial denominators of the continued fraction of $e$."
|
"To compute the convergents ourselves, we'll first make a generator for the partial denominators of the continued fraction of $e$."
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
|
@ -48,14 +48,9 @@
|
||||||
"source": [
|
"source": [
|
||||||
"from itertools import count, chain\n",
|
"from itertools import count, chain\n",
|
||||||
"\n",
|
"\n",
|
||||||
"def partial_denominators_e(n):\n",
|
"def partial_denominators_e():\n",
|
||||||
" yield 2\n",
|
" yield 2\n",
|
||||||
" denominators = chain.from_iterable((1, 2 * k, 1) for k in count(1))\n",
|
" yield from chain.from_iterable((1, 2 * k, 1) for k in count(1))"
|
||||||
" for (i, b) in enumerate(denominators):\n",
|
|
||||||
" if i >= n:\n",
|
|
||||||
" break\n",
|
|
||||||
" \n",
|
|
||||||
" yield b"
|
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
|
@ -63,7 +58,7 @@
|
||||||
"id": "33288cd0",
|
"id": "33288cd0",
|
||||||
"metadata": {},
|
"metadata": {},
|
||||||
"source": [
|
"source": [
|
||||||
"Then write a function for computing a continued fraction from a sequence of partial denominators (outside of SageMath, you might want to use a [fraction type](https://docs.python.org/3/library/fractions.html))."
|
"Then we'll apply a simple algorithm for computing [convergents using the partial denominators](https://en.wikipedia.org/wiki/Simple_continued_fraction) (outside of SageMath, you might want to use a [fraction type](https://docs.python.org/3/library/fractions.html))."
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
|
@ -73,13 +68,21 @@
|
||||||
"metadata": {},
|
"metadata": {},
|
||||||
"outputs": [],
|
"outputs": [],
|
||||||
"source": [
|
"source": [
|
||||||
"def cf(denominators):\n",
|
"def convergents(partial_denoms):\n",
|
||||||
" a = next(denominators)\n",
|
" h, hprev = 1, 0\n",
|
||||||
" \n",
|
" k, kprev = 0, 1\n",
|
||||||
" try:\n",
|
" for b in partial_denoms:\n",
|
||||||
" return a + 1 / cf(denominators)\n",
|
" h, hprev = b * h + hprev, h\n",
|
||||||
" except StopIteration:\n",
|
" k, kprev = b * k + kprev, k\n",
|
||||||
" return a"
|
" yield h/k"
|
||||||
|
]
|
||||||
|
},
|
||||||
|
{
|
||||||
|
"cell_type": "markdown",
|
||||||
|
"id": "99b790f3",
|
||||||
|
"metadata": {},
|
||||||
|
"source": [
|
||||||
|
"Now just iterate until we reach the 100th convergent."
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
|
@ -100,7 +103,11 @@
|
||||||
}
|
}
|
||||||
],
|
],
|
||||||
"source": [
|
"source": [
|
||||||
"sum(cf(partial_denominators_e(99)).numerator().digits())"
|
"for (i, c) in enumerate(convergents(partial_denominators_e())):\n",
|
||||||
|
" if i == 99:\n",
|
||||||
|
" break\n",
|
||||||
|
"\n",
|
||||||
|
"sum(c.numerator().digits())"
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
|
|
|
@ -19,10 +19,10 @@
|
||||||
{
|
{
|
||||||
"data": {
|
"data": {
|
||||||
"text/plain": [
|
"text/plain": [
|
||||||
"[(-sqrt(2)*(2*sqrt(2) + 3)^t + sqrt(2)*(-2*sqrt(2) + 3)^t - 3/2*(2*sqrt(2) + 3)^t - 3/2*(-2*sqrt(2) + 3)^t,\n",
|
"[(sqrt(2)*(2*sqrt(2) + 3)^t - sqrt(2)*(-2*sqrt(2) + 3)^t + 3/2*(2*sqrt(2) + 3)^t + 3/2*(-2*sqrt(2) + 3)^t,\n",
|
||||||
" 3/4*sqrt(2)*(2*sqrt(2) + 3)^t - 3/4*sqrt(2)*(-2*sqrt(2) + 3)^t + (2*sqrt(2) + 3)^t + (-2*sqrt(2) + 3)^t),\n",
|
" -3/4*sqrt(2)*(2*sqrt(2) + 3)^t + 3/4*sqrt(2)*(-2*sqrt(2) + 3)^t - (2*sqrt(2) + 3)^t - (-2*sqrt(2) + 3)^t),\n",
|
||||||
" (sqrt(2)*(2*sqrt(2) + 3)^t - sqrt(2)*(-2*sqrt(2) + 3)^t + 3/2*(2*sqrt(2) + 3)^t + 3/2*(-2*sqrt(2) + 3)^t,\n",
|
" (-sqrt(2)*(2*sqrt(2) + 3)^t + sqrt(2)*(-2*sqrt(2) + 3)^t - 3/2*(2*sqrt(2) + 3)^t - 3/2*(-2*sqrt(2) + 3)^t,\n",
|
||||||
" -3/4*sqrt(2)*(2*sqrt(2) + 3)^t + 3/4*sqrt(2)*(-2*sqrt(2) + 3)^t - (2*sqrt(2) + 3)^t - (-2*sqrt(2) + 3)^t)]"
|
" 3/4*sqrt(2)*(2*sqrt(2) + 3)^t - 3/4*sqrt(2)*(-2*sqrt(2) + 3)^t + (2*sqrt(2) + 3)^t + (-2*sqrt(2) + 3)^t)]"
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
"execution_count": 1,
|
"execution_count": 1,
|
||||||
|
@ -97,51 +97,20 @@
|
||||||
"## Solving Pell equations\n",
|
"## Solving Pell equations\n",
|
||||||
"Lagrange proved that if $(x_0, y_0)$ is a solution to\n",
|
"Lagrange proved that if $(x_0, y_0)$ is a solution to\n",
|
||||||
"$$x^2 - dy^2 = 1$$\n",
|
"$$x^2 - dy^2 = 1$$\n",
|
||||||
"then $\\frac{x_0}{y_0}$ is a [convergent of the continued fraction](https://en.wikipedia.org/wiki/Simple_continued_fraction) of $\\sqrt{d}$. This is great for us, since we can write generators for computing convergents of square roots (FYI, SageMath can do this with built-in methods: `continued_fraction(sqrt(d)).convergents()`)."
|
"then $\\frac{x_0}{y_0}$ is a [convergent of the continued fraction](https://en.wikipedia.org/wiki/Simple_continued_fraction) of $\\sqrt{d}$. This is great for us, since there are algorithms to compute these convergents. We'll use SageMath here; see [problem 64](https://projecteuler.net/problem=64) for how to compute the partial denominators of the continued fraction of a square root, and see [problem 65](https://projecteuler.net/problem=65) for an algorithm that uses partial denominators to compute convergents of continued fractions.\n",
|
||||||
|
"\n",
|
||||||
|
"Here, we iterate over each convergent to see if its numerator and denominator are a solution to the Pell equation."
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
"cell_type": "code",
|
"cell_type": "code",
|
||||||
"execution_count": 4,
|
"execution_count": 4,
|
||||||
"id": "c51f9b7e",
|
|
||||||
"metadata": {},
|
|
||||||
"outputs": [],
|
|
||||||
"source": [
|
|
||||||
"def continued_fraction_sqrt(d):\n",
|
|
||||||
" x = sqrt(d)\n",
|
|
||||||
" while True:\n",
|
|
||||||
" b = floor(x)\n",
|
|
||||||
" yield b\n",
|
|
||||||
" x = (x - b)^-1\n",
|
|
||||||
" \n",
|
|
||||||
" \n",
|
|
||||||
"def convergents(partial_denoms):\n",
|
|
||||||
" h, hprev = 1, 0\n",
|
|
||||||
" k, kprev = 0, 1\n",
|
|
||||||
" for b in partial_denoms:\n",
|
|
||||||
" h, hprev = b * h + hprev, h\n",
|
|
||||||
" k, kprev = b * k + kprev, k\n",
|
|
||||||
" yield h/k"
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "markdown",
|
|
||||||
"id": "d5cbf5eb",
|
|
||||||
"metadata": {},
|
|
||||||
"source": [
|
|
||||||
"Then we can just iterate over each convergent to see if its numerator and denominator are a solution to the Pell equation."
|
|
||||||
]
|
|
||||||
},
|
|
||||||
{
|
|
||||||
"cell_type": "code",
|
|
||||||
"execution_count": 5,
|
|
||||||
"id": "5d95125c",
|
"id": "5d95125c",
|
||||||
"metadata": {},
|
"metadata": {},
|
||||||
"outputs": [],
|
"outputs": [],
|
||||||
"source": [
|
"source": [
|
||||||
"def pell_fundamental_solution(d):\n",
|
"def pell_fundamental_solution(d):\n",
|
||||||
" partial_denoms = continued_fraction_sqrt(d)\n",
|
" for f in continued_fraction(sqrt(d)).convergents():\n",
|
||||||
" for f in convergents(partial_denoms):\n",
|
|
||||||
" x, y = f.as_integer_ratio()\n",
|
" x, y = f.as_integer_ratio()\n",
|
||||||
" if x^2 - d*y^2 == 1:\n",
|
" if x^2 - d*y^2 == 1:\n",
|
||||||
" return (x, y)"
|
" return (x, y)"
|
||||||
|
@ -157,7 +126,7 @@
|
||||||
},
|
},
|
||||||
{
|
{
|
||||||
"cell_type": "code",
|
"cell_type": "code",
|
||||||
"execution_count": 6,
|
"execution_count": 5,
|
||||||
"id": "03d7ba9d",
|
"id": "03d7ba9d",
|
||||||
"metadata": {},
|
"metadata": {},
|
||||||
"outputs": [
|
"outputs": [
|
||||||
|
@ -167,7 +136,7 @@
|
||||||
"661"
|
"661"
|
||||||
]
|
]
|
||||||
},
|
},
|
||||||
"execution_count": 6,
|
"execution_count": 5,
|
||||||
"metadata": {},
|
"metadata": {},
|
||||||
"output_type": "execute_result"
|
"output_type": "execute_result"
|
||||||
}
|
}
|
||||||
|
|
Loading…
Reference in New Issue