redo implementation of fundamental pell solver
This commit is contained in:
parent
ccbb6027c8
commit
2c41f67244
|
@ -19,10 +19,10 @@
|
|||
{
|
||||
"data": {
|
||||
"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",
|
||||
" -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",
|
||||
" 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)]"
|
||||
"[(-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",
|
||||
" (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)]"
|
||||
]
|
||||
},
|
||||
"execution_count": 1,
|
||||
|
@ -95,9 +95,9 @@
|
|||
"## Solving Pell equations\n",
|
||||
"Lagrange proved that if $(x_0, y_0)$ is a solution to\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 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",
|
||||
"then $\\frac{x_0}{y_0}$ is a [convergent of the continued fraction](https://en.wikipedia.org/wiki/Simple_continued_fraction) of $\\sqrt{d}$. Specifically, if $p$ is the period of the continued fraction, then the first solution will be the $(p-1)$th convergent if $p$ is even, and the $(2p-1)$th convergent if $p$ is odd.\n",
|
||||
"\n",
|
||||
"Here, we iterate over each convergent to see if its numerator and denominator are a solution to the Pell equation."
|
||||
"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 (using SageMath's constructions adds overhead that makes this implementation a little slow, but it makes the code easier to read - and it's still considerably faster than using `solve_diophantine`)."
|
||||
]
|
||||
},
|
||||
{
|
||||
|
@ -108,10 +108,16 @@
|
|||
"outputs": [],
|
||||
"source": [
|
||||
"def pell_fundamental_solution(d):\n",
|
||||
" for f in continued_fraction(sqrt(d)).convergents():\n",
|
||||
" x, y = f.as_integer_ratio()\n",
|
||||
" if x^2 - d*y^2 == 1:\n",
|
||||
" return (x, y)"
|
||||
" K.<s> = QuadraticField(d)\n",
|
||||
" f = continued_fraction(s)\n",
|
||||
" p = f.period_length()\n",
|
||||
" if p % 2 == 0:\n",
|
||||
" n = p - 1\n",
|
||||
" else:\n",
|
||||
" n = 2 * p - 1\n",
|
||||
" \n",
|
||||
" x, y = f.convergent(n).as_integer_ratio()\n",
|
||||
" return x, y"
|
||||
]
|
||||
},
|
||||
{
|
||||
|
|
Loading…
Reference in New Issue