daviddoji/pytudes
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Add files via upload

Browse Source
This commit is contained in:
Peter Norvig
2017-07-06 01:06:25 -04:00
committed by GitHub
parent 601228e6e8
commit 50c7fdb519
1 changed files with 224 additions and 43 deletions
Download Patch File Download Diff File Expand all files Collapse all files

265
Beal.ipynb
View File

@@ -239,24 +239,23 @@
},
"outputs": [
{
"ename": "NameError",
"evalue": "name 'esimpsons' is not defined",
"output_type": "error",
"traceback": [
"\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
"\u001b[0;31mNameError\u001b[0m Traceback (most recent call last)",
"\u001b[0;32m<ipython-input-5-c5d9ee9dd94c>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mesimpsons\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mrange\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m1000\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m2000\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0;36m11\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m12\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m13\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
"\u001b[0;31mNameError\u001b[0m: name 'esimpsons' is not defined"
"data": {
"text/plain": [
"(1782, 1841, 1922, 12)"
]
},
"execution_count": 5,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"esimpsons(range(1000, 2000), [11, 12, 13])"
"simpsons(range(1000, 2000), [11, 12, 13])"
]
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 6,
"metadata": {
"button": false,
"new_sheet": false,
@@ -264,7 +263,18 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"(3987, 4365, 4472, 12)"
]
},
"execution_count": 6,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"simpsons(range(3000, 5000), [12])"
]
@@ -284,7 +294,7 @@
"In October 2015 I looked back at my original [program from 2000](http://norvig.com/beal2000.html).\n",
"I ported it from Python 1.5 to 3.5 (by putting parens around the argument to `print` and adding `long = int`). The program runs 250 times faster today than it did in 2000, a tribute to both computer hardware engineers and the developers of the Python interpreter.\n",
"\n",
"I found that I was a bit confused about the definition of [the problem in 2000](https://web.archive.org/web/19991127081319/http://bealconjecture.com/). I thought then that, *by definition*, $A$ and $B$ could not have a common factor, but actually, the definition of the conjecture only rules out examples where all three of $A, B, C$ share a common factor. Mark Tiefenbruck (and later Edward P. Berlin and Shen Lixing) wrote to point out that my thought was actually correct, not by definition, but by derivation: if $A$ and $B$ have a commmon prime factor $p$, then the sum of $A^x + B^y$ must also have that factor $p$, and since $A^x + B^y = C^z$, then $C^z$ and hence $C$ must have the factor $p$. So I was wrong twice, and in this case two wrongs made a right.\n",
"I found that I was a bit confused about the definition of [the problem in 2000](https://web.archive.org/web/19991127081319/http://bealconjecture.com/). I thought then that, *by definition*, $A$ and $B$ could not have a common factor, but actually, the definition of the conjecture only rules out examples where all three of $A, B, C$ share a common factor. Mark Tiefenbruck (and later Edward P. Berlin and Shen Lixing) wrote to point out that my thought was actually correct, not by definition, but by derivation: if $A$ and $B$ have a commmon prime factor $p$, then the sum of $A^x + B^y$ must also have that factor $p$, and since $A^x + B^y = C^z$, then $C^z$ and hence $C$ must have the factor $p$. So I was wrong twice, and in this case two wrongs did make a right.\n",
"\n",
"Mark Tiefenbruck also suggested an optimization: only consider exponents that are odd primes, or 4. The idea is that a number like 512 can be expressed as either $2^9$ or $8^3$, and my program doesn't need to consider both. In general, any time we have a composite exponent, such as $b^{qp}$, where $p$ is prime, we should ignore $b^{(qp)}$, and instead consider only $(b^q)^p$. There's one complication to this scheme: 2 is a prime, but 2 is not a valid exponent for a Beal counterexample. So we will allow 4 as an exponent, as well as all odd primes up to `max_x`.\n",
"\n",
@@ -293,7 +303,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 7,
"metadata": {
"button": false,
"new_sheet": false,
@@ -358,7 +368,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 8,
"metadata": {
"button": false,
"new_sheet": false,
@@ -366,7 +376,16 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"CPU times: user 256 ms, sys: 1.44 ms, total: 257 ms\n",
"Wall time: 256 ms\n"
]
}
],
"source": [
"%time beal(100, 100)"
]
@@ -381,12 +400,12 @@
}
},
"source": [
"The execution time goes up roughly with the square of `max_A`, so with 5 times more `A` values, this computation takes about 25 times longer:"
"The execution time goes up roughly with the square of `max_A`, so the following should take about 100 times longer:"
]
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 9,
"metadata": {
"button": false,
"new_sheet": false,
@@ -394,9 +413,18 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"CPU times: user 29.1 s, sys: 25.2 ms, total: 29.2 s\n",
"Wall time: 29.2 s\n"
]
}
],
"source": [
"%time beal(500, 100)"
"%time beal(1000, 100)"
]
},
{
@@ -421,7 +449,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 10,
"metadata": {
"button": false,
"new_sheet": false,
@@ -429,7 +457,23 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"{1: [1],\n",
" 2: [8, 16, 32, 128],\n",
" 3: [27, 81, 243, 2187],\n",
" 4: [64, 256, 1024, 16384],\n",
" 5: [125, 625, 3125, 78125],\n",
" 6: [216, 1296, 7776, 279936]}"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Apowers = make_Apowers(6, 10)\n",
"Apowers"
@@ -450,7 +494,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 11,
"metadata": {
"button": false,
"new_sheet": false,
@@ -458,7 +502,38 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"{1: 1,\n",
" 8: 2,\n",
" 16: 2,\n",
" 27: 3,\n",
" 32: 2,\n",
" 64: 4,\n",
" 81: 3,\n",
" 125: 5,\n",
" 128: 2,\n",
" 216: 6,\n",
" 243: 3,\n",
" 256: 4,\n",
" 625: 5,\n",
" 1024: 4,\n",
" 1296: 6,\n",
" 2187: 3,\n",
" 3125: 5,\n",
" 7776: 6,\n",
" 16384: 4,\n",
" 78125: 5,\n",
" 279936: 6}"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Czroots = make_Czroots(Apowers)\n",
"Czroots"
@@ -466,9 +541,20 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 12,
"metadata": {},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"3"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Czroots[243]"
]
@@ -493,7 +579,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 13,
"metadata": {
"button": false,
"new_sheet": false,
@@ -501,7 +587,25 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"3 ** 3 + 6 ** 3 == 3 ** 5 == 243\n",
"7 ** 7 + 49 ** 3 == 98 ** 3 == 941192\n",
"8 ** 4 + 16 ** 3 == 2 ** 13 == 8192\n",
"8 ** 5 + 32 ** 3 == 16 ** 4 == 65536\n",
"9 ** 3 + 18 ** 3 == 9 ** 4 == 6561\n",
"16 ** 5 + 32 ** 4 == 8 ** 7 == 2097152\n",
"17 ** 4 + 34 ** 4 == 17 ** 5 == 1419857\n",
"19 ** 4 + 38 ** 3 == 57 ** 3 == 185193\n",
"27 ** 3 + 54 ** 3 == 3 ** 11 == 177147\n",
"28 ** 3 + 84 ** 3 == 28 ** 4 == 614656\n",
"34 ** 5 + 51 ** 4 == 85 ** 4 == 52200625\n"
]
}
],
"source": [
"def gcd(a, b): return 1\n",
"\n",
@@ -523,7 +627,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 14,
"metadata": {
"button": false,
"new_sheet": false,
@@ -531,7 +635,18 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"{True}"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"{3 ** 3 + 6 ** 3 == 3 ** 5 == 243,\n",
" 7 ** 7 + 49 ** 3 == 98 ** 3 == 941192,\n",
@@ -561,7 +676,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 15,
"metadata": {
"button": false,
"new_sheet": false,
@@ -591,7 +706,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 16,
"metadata": {
"button": false,
"new_sheet": false,
@@ -599,7 +714,18 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"'tests pass'"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"def tests():\n",
" assert make_Apowers(6, 10) == {\n",
@@ -700,7 +826,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 17,
"metadata": {
"button": false,
"collapsed": true,
@@ -715,7 +841,7 @@
"from itertools import combinations, product\n",
"from collections import defaultdict\n",
" \n",
"def beal_modp(max_A, max_x, m=2**31-1):\n",
"def beal_modm(max_A, max_x, m=2**31-1):\n",
" \"\"\"See if any A ** x + B ** y equals some C ** z (mod p), with gcd(A, B) == 1.\n",
" If so, verify that the equation works without the (mod m).\n",
" Consider any 1 <= A,B <= max_A and x,y <= max_x, with x,y prime or 4.\"\"\"\n",
@@ -765,9 +891,25 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 18,
"metadata": {},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"{1: [(1, 3)],\n",
" 2: [(8, 3), (16, 4), (32, 5), (128, 7)],\n",
" 3: [(27, 3), (81, 4), (243, 5), (187, 7)],\n",
" 4: [(64, 3), (256, 4), (24, 5), (384, 7)],\n",
" 5: [(125, 3), (625, 4), (125, 5), (125, 7)],\n",
" 6: [(216, 3), (296, 4), (776, 5), (936, 7)]}"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"Apowers = make_Apowers_modm(6, 10, 1000)\n",
"Apowers"
@@ -789,7 +931,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 19,
"metadata": {
"button": false,
"new_sheet": false,
@@ -797,7 +939,37 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"data": {
"text/plain": [
"defaultdict(list,\n",
" {1: [(1, 3)],\n",
" 8: [(2, 3)],\n",
" 16: [(2, 4)],\n",
" 24: [(4, 5)],\n",
" 27: [(3, 3)],\n",
" 32: [(2, 5)],\n",
" 64: [(4, 3)],\n",
" 81: [(3, 4)],\n",
" 125: [(5, 3), (5, 5), (5, 7)],\n",
" 128: [(2, 7)],\n",
" 187: [(3, 7)],\n",
" 216: [(6, 3)],\n",
" 243: [(3, 5)],\n",
" 256: [(4, 4)],\n",
" 296: [(6, 4)],\n",
" 384: [(4, 7)],\n",
" 625: [(5, 4)],\n",
" 776: [(6, 5)],\n",
" 936: [(6, 7)]})"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"make_Czroots_modm(Apowers)"
]
@@ -817,7 +989,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 20,
"metadata": {
"button": false,
"new_sheet": false,
@@ -825,9 +997,18 @@
"read_only": false
}
},
"outputs": [],
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"CPU times: user 35.5 s, sys: 44.1 ms, total: 35.5 s\n",
"Wall time: 35.5 s\n"
]
}
],
"source": [
"%time beal_modp(500, 100)"
"%time beal_modm(1000, 100)"
]
},
{
@@ -840,7 +1021,7 @@
}
},
"source": [
"We don't see a big speedup here, but the idea is that as we start dealing with much larger integers, this version will be faster. I could improve this version by caching certain computations, managing the memory layout better, moving some computations out of loops, considering using multiple different numbers as the modulus (as in a Bloom filter), finding a way to parallelize the program, and re-coding in a faster compiled language (such as C++ or Go or Julia). Then I could invest thousands (or millions) of CPU hours searching for counterexamples. \n",
"We don't see a speedup here, but the idea is that as we start dealing with much larger integers, this version should be faster. I could improve this version by caching certain computations, managing the memory layout better, moving some computations out of loops, considering using multiple different numbers as the modulus (as in a Bloom filter), finding a way to parallelize the program, and re-coding in a faster compiled language (such as C++ or Go or Julia). Then I could invest thousands (or millions) of CPU hours searching for counterexamples. \n",
"\n",
"But [Witold Jarnicki](https://plus.sandbox.google.com/+WitoldJarnicki/posts) and [David Konerding](http://www.konerding.com/~dek/) already did that: they wrote a C++ program that, in parallel across thousands of machines, searched for $A, B$ up to 200,000 and $x, y$ up to 5,000, but found no counterexamples. So I don't think it is worthwhile to continue on that path.\n",
"\n",