In [ ]:
A Probabilist's Approach to the Legendre Conjecture¶
Legendre's Conjecture: For every $n \in \mathbb{N}$, we have that there exists a prime number $p$ such that $n^2 < p < (n+1)^2 $.
Empirically, what is the probability that this conjecture is true? (Inspired by Feynman's work on Fermat's last theorem)
Prime Number Theorem (Corollary): For any integer N, the probability that N is prime tends to $\frac{1}{\log{N}}$.
$$ P(\text{no prime between}\; n^2 \; \text{and} \; (n+1)^2) = \prod_{i=n^2}^{(n+1)^2} 1 - \frac{1}{\log{i}} $$$$ \implies P(\text{Legendre's conjecture is false}) = \sum_{n=N}^{\infty} \prod_{i=n^2}^{(n+1)^2} 1-\frac{1}{\log{i}} $$Estimating this product¶
We have $$\log{\left(\prod_{i=n^2}^{(n+1)^2} 1-\frac{1}{\log{i}}\right)} = \sum_{i=n^2}^{(n+1)^2} \left( \log{(1-1/\log{i})} \right) \approx \int_{n^2}^{(n+1)^2} \log{(1-1/\log(x))}dx$$
Amazingly, this integral can be expressed in closed form - $$\int \log(1-1/\log(x)) = -e \text{Ei}(\log(x)-1) + \text{li}(x) + x \log(1-1/\log(x))$$
We therefore have that $$P := \prod_{i=n^2}^{(n+1)^2} (1-\frac{1}{\log{i}}) \approx P' := \exp(e \; [\text{Ei}(2\log n-1)-\text{Ei}(2\log(n+1) - 1)]) + \text{li}((n+1)^2) - \text{li}(n^2))\left(1-\frac{1}{2\log n}\right)^{-n^2}\left(1-\frac{1}{2\log(n+1)}\right)^{(n+1)^2} $$
Using the definitions of the exponential integral and the logarithmic integral, we get
$$P \approx \exp\left(e \left(\int_{2\log n-1}^{2\log(n+1)-1} \frac{e^t}{t}dt\right)+\int_{n^2}^{(n+1)^2} \frac{dt}{\log{t}}\right) \left(1-\frac{1}{2\log n}\right)^{-n^2}\left(1-\frac{1}{2\log(n+1)}\right)^{(n+1)^2} $$Now, since we have that $$\left(1-\frac{1}{2\log n}\right)^{-n^2}\left(1-\frac{1}{2\log(n+1)}\right)^{(n+1)^2} = \left(\frac{2-\frac{1}{\log{n}+1}}{2-\frac{1}{\log{n}}}\right)^{n^2} \left(1-\frac{1}{2\log(n+1)}\right)^{2n+1} $$
we have that $$P \approx P'' := \exp\left(e \left(\int_{2\log n-1}^{2\log(n+1)-1} \frac{e^t}{t}dt\right)+\int_{n^2}^{(n+1)^2} \frac{dt}{\log{t}}\right) \left(\frac{2-\frac{1}{\log{n}+1}}{2-\frac{1}{\log{n}}}\right)^{n^2} \left(1-\frac{1}{2\log(n+1)}\right)^{2n+1}$$
In [9]:
# We sum up to n=10000 terms from N=10
from math import log
from functools import reduce
prod = lambda arr: reduce(lambda a,b: a*b, arr)
for N in np.linspace(10,500,100):
N = int(N)
print(N, sum([prod([1-1/log(i) for i in range(n**2, (n+1)**2+1)]) for n in range(10,N+1)]))
10 0.005196273509080254 14 0.01707694784532844 19 0.022336237875580096 24 0.024011500993202848 29 0.02456434731521764 34 0.024752536009686448 39 0.02481832217024484 44 0.024841852658789763 49 0.024850438741033355 54 0.024853627248810017 59 0.024854829959129385 64 0.024855290021491377 69 0.024855468247351763 74 0.024855538091494375 79 0.02485556575286258 84 0.024855576815001702 89 0.024855581278869075 94 0.02485558309528063 99 0.02485558384019228 104 0.02485558414792023 108 0.02485558425840319 113 0.024855584322151594 118 0.024855584348970977 123 0.024855584360320163 128 0.024855584365149487 133 0.02485558436721532 138 0.02485558436810346 143 0.02485558436848712 148 0.02485558436865361 153 0.024855584368726175 158 0.024855584368757934 163 0.02485558436877189 168 0.02485558436877805 173 0.024855584368780777 178 0.024855584368781988 183 0.02485558436878253 188 0.024855584368782775 193 0.02485558436878288 198 0.024855584368782928 203 0.02485558436878295 207 0.02485558436878296 212 0.02485558436878296 217 0.02485558436878296 222 0.02485558436878296 227 0.02485558436878296 232 0.02485558436878296 237 0.02485558436878296 242 0.02485558436878296 247 0.02485558436878296 252 0.02485558436878296 257 0.02485558436878296 262 0.02485558436878296 267 0.02485558436878296 272 0.02485558436878296 277 0.02485558436878296 282 0.02485558436878296 287 0.02485558436878296 292 0.02485558436878296 297 0.02485558436878296 302 0.02485558436878296 306 0.02485558436878296 311 0.02485558436878296 316 0.02485558436878296 321 0.02485558436878296 326 0.02485558436878296 331 0.02485558436878296 336 0.02485558436878296 341 0.02485558436878296 346 0.02485558436878296 351 0.02485558436878296 356 0.02485558436878296 361 0.02485558436878296 366 0.02485558436878296 371 0.02485558436878296 376 0.02485558436878296 381 0.02485558436878296 386 0.02485558436878296 391 0.02485558436878296 396 0.02485558436878296 401 0.02485558436878296 405 0.02485558436878296 410 0.02485558436878296 415 0.02485558436878296 420 0.02485558436878296 425 0.02485558436878296 430 0.02485558436878296 435 0.02485558436878296 440 0.02485558436878296 445 0.02485558436878296 450 0.02485558436878296 455 0.02485558436878296 460 0.02485558436878296 465 0.02485558436878296 470 0.02485558436878296 475 0.02485558436878296 480 0.02485558436878296 485 0.02485558436878296 490 0.02485558436878296 495 0.02485558436878296 500 0.02485558436878296
In [18]:
# Let's try more digits with mpmath
from mpmath import *
mp.dps = 100; mp.pretty = False
from functools import reduce
prod = lambda arr: reduce(lambda a,b: a*b, arr)
out = 0
for n in arange(10,500):
out += prod([1-1/log(i) for i in arange(n**2, (n+1)**2+1)])
if n%10==0: print(n, out)
10.0 0.005196273509080254847009590716693511629156854414891029893227564618764245586915259508555271288886037229 20.0 0.02283633772007478906453499354964824319445739465678895633554391445354217680655266559930421962661168564 30.0 0.02461943671806240319455749447494419794474767922492881886186933920324224216587415777024306304327835058 40.0 0.02482510994832935906756473138775584192306845291447143013687629679057455281538092674198926981837073502 50.0 0.02485134813750385803193277518136496302651787755103058126933492722085466771840821556852906146235361709 60.0 0.02485495998347029652553789504642911134401849449715201595733620101148680745419005392777890765650893874 70.0 0.02485548783865370192630023288144958099548644845242217717159009035492083736733978705073085767955632548 80.0 0.02485556883603585736021763751138588059718511252840316126847251030479101424809247121786786732338851662 90.0 0.02485558178234395784529517357529406904465937423595310706280644169017614133444712565222593287323389048 100.0 0.02485558392507732060983508735185638268072280168441674799606005172754484893783619675740412400609012339 110.0 0.02485558429063797506046611852446301954467393321373143548819071339991474387893565992072728072451110966 120.0 0.02485558435469368838663720216459312742356526897962543737780450010252054667592094385326652204048223392 130.0 0.02485558436618853019920593800848866840518218015863491681773861340273790158193764272496648640553880588 140.0 0.02485558436829598947161021346819686921800584493306563776275317857362986750351132381726514140647474299 150.0 0.02485558436868995871972984341231987553881304653711349613983971726592266239791944338065411858235209484 160.0 0.02485558436876492773784426303107499350258331359086899275904019737646833961032062364417405049786809053 170.0 0.02485558436877942857004897417274029120010286962029723941544690634234547558898050261004567537734002481 180.0 0.02485558436878227601937318098542907928496013841841400885854800370715489421098241654964969777212932198 190.0 0.02485558436878284304257228806495529752591317762282809840704277736366689458548260201413360365384093825 200.0 0.02485558436878295743952119265153786370288925833295563432293173434502189483612392393481037805717383335 210.0 0.02485558436878298080268405896356098275840125073963200675424353605785910759337724916943831061090390724 220.0 0.0248555843687829856291294845476600171532661442942590894801116453456963501652306536368427614547472179 230.0 0.02485558436878298663701416408598478615066490043302972027686653460800093593504946339455108677520118941 240.0 0.02485558436878298684964207780502173583336809297076639656474054001422573554203212421431423464918917783 250.0 0.02485558436878298689493383501911274407585704999109631869836019136660999609037160936559371055829122727 260.0 0.02485558436878298690467010947321250270656277328595461546926879628282608107793245700733296125791023963 270.0 0.02485558436878298690678138789151264458743918300078967323628354647683678506817220948486556344329442668 280.0 0.02485558436878298690724302242197172986585734839134400150514524040175040640830850581783992229091473178 290.0 0.0248555843687829869073447612837619312746851237133079164756892554244494048425902483429898824244542621 300.0 0.0248555843687829869073673535625463569635124079275209331567353223643125198616489447164097705291100464 310.0 0.02485558436878298690737240689813083240530335015750070230478758909120905621990606498502560589651015055 320.0 0.02485558436878298690737354508061153516177173035163910564012123415764621129973729074331887936598439299 330.0 0.02485558436878298690737380315286201184390026982869958881926378260496185684274599593821390378958169488 340.0 0.02485558436878298690737386204464877849898740393135598217664806733906593587540118265170217817840613305 350.0 0.02485558436878298690737387556687056051635606378526402769677080088741160446110703942885688868199490435 360.0 0.02485558436878298690737387869024461814724778249139378922615465837161270661372289971011663044807450711 370.0 0.02485558436878298690737387941583493910287026920840295378631506501275230520321536726095525721804904191 380.0 0.02485558436878298690737387958533312832766922912543380041543083703523181884516825090641723241191903441 390.0 0.02485558436878298690737387962514058773597734830810704373334677673189988864331722823497122125309455387 400.0 0.02485558436878298690737387963453812609851444370709638390683772544610644803590314107795198069817343439 410.0 0.02485558436878298690737387963676780724402174568666464256605142708361861973478283961352815396440548881 420.0 0.02485558436878298690737387963729940502712013107018893895366603714392487369076062496653818849626741491 430.0 0.02485558436878298690737387963742674703977950267747377260654331967357931694999026085672015358193084394 440.0 0.02485558436878298690737387963745739123314293517205895790724007393248463218686854957060038037678757089 450.0 0.02485558436878298690737387963746479845804914951258577066726670186900615682970747684452970816294979339 460.0 0.02485558436878298690737387963746659666444936063176781803891896182442094633120610346787779944101500427 470.0 0.0248555843687829869073738796374670350424251190028244865005731797210343195222705527434844586800247167 480.0 0.02485558436878298690737387963746714235092112465185290682019418379284527799790506679036182968404602501 490.0 0.02485558436878298690737387963746716872328301350544613075203123556312014710804961879854997149073995345
In [21]:
# And trying P'
out = 0
for n in arange(10,500):
if n%10==0: print(n, out)
out += exp(e*ei(2*log(n)-1)-e*ei(2*log(n+1)-1)-li(n**2)+li((n+1)**2))*(1-1/(2*log(n)))**(-n**2)*(1-1/(2*log(n+1)))**((n+1)**2)
10.0 0 20.0 0.02779157631077098261120826272861472447750611577299816730649958813651946021331165684877380054328373331 30.0 0.03044196208687949895358984879034041305237755952661883705302101417852895871874869686395512749275890839 40.0 0.03073822394868029145237338630214509367265160069871876196998066641209331548258017754134731819797512408 50.0 0.03077525365839126402154597908408678745691110356043537196245094966764204878398809134435028097994055592 60.0 0.03078027625585673627631731421973103599793889902757214847173270599293795169562957485379858768748251073 70.0 0.03078100198009662091653604658743779618898450152407443293033375025997828815184148303435338090203025156 80.0 0.03078111232542702337266153705444182060330908489170796064334854297777315247434373021396768190105316503 90.0 0.03078112982926019662066535288923382888994770319379927759674908836183066953106732734792273060013842107 100.0 0.03078113270770176547266556825687225333805912337712388632198813352500789954072259222780481784952632691 110.0 0.03078113319605011268428729927484415518381868460042475865356842451996106432842116614701885453321244693 120.0 0.03078113328120429127762017075259152365091918817916072651196919230540897417190589824506542062624704773 130.0 0.03078113329641911353259852249855604249705360014775836108258795558393021871067121811524824445180262326 140.0 0.03078113329919776042464028847867355221984750609590166846503785005269699625577493350532127814516446898 150.0 0.03078113329971537717906087626317710529095034915167283553108592849131091308701242069032568297601850734 160.0 0.03078113329981356004283130467510828948775305785557342629183870699006172665089219593244699699355568913 170.0 0.03078113329983249527347352601984604472018145354464551399026449840269906388309177824252515084550220121 180.0 0.03078113329983620342348735703013976631579745244680442563124144702146769133695658866418955283603159269 190.0 0.03078113329983693999085300994484358415472601148038852836057838702466310020019988726324174480247436535 200.0 0.03078113329983708824689956635986327587150699966720286491980029621307199607594954940015271705513331248 210.0 0.03078113329983711845914250845079635392276569764825014253980034145766127559374189386243428530913885774 220.0 0.03078113329983712468775929355415334924940710726478900113892247985105326680033006082257678936101556863 230.0 0.03078113329983712598595778237680040677600458399479484346366797239411840919461337212758181270183768794 240.0 0.03078113329983712625933625518467116231934748349535401102313152139794556030456126201484449823110651421 250.0 0.03078113329983712631746889634229915531637472888846881356271729059541644567289924926384217818724067823 260.0 0.03078113329983712632994530444185493443278055986122499225394239055758467230092407206259328054278135676 270.0 0.03078113329983712633264660972062199109492329833536507582424425747220935484039810688338389463873541498 280.0 0.03078113329983712633323639003957548768930574570432314157065203420558358257437532377895791461432228956 290.0 0.03078113329983712633336618936043763513938072092812890047563642962764083919507867530085060069328681841 300.0 0.03078113329983712633339497437611645165220200842383298372190901969810618749534472712150265131871777161 310.0 0.0307811332998371263334014046676039942752059697195107827249876075459156577049191977565499086298196776 320.0 0.03078113329983712633340285121880043611646457196037922832189930649215047743946578451864949767310035188 330.0 0.03078113329983712633340317882658930587290977507741667740142144210947199364382947920332628323899528428 340.0 0.03078113329983712633340325350218849407539897190733087560346682036169722576409355940734976587777124034 350.0 0.03078113329983712633340327062999773495186616855811888539158091031826828760507565896309236646825196437 360.0 0.03078113329983712633340327458207108763861798020545521441643287262946562539562011666920398799685825755 370.0 0.03078113329983712633340327549925418503546300754376195277588659512142518117749980731687125884446137695 380.0 0.03078113329983712633340327571330086223414370510339052290389080237096132448390877958874906814636236318 390.0 0.0307811332998371263334032757635236850901626431219938644790603501756787618765534288313906108828487321 400.0 0.03078113329983712633340327577536929715407152634553916120984062482987750071571905332613616677769570413 410.0 0.03078113329983712633340327577817735224738355522819413504597039666640251624749596303678057882199041712 420.0 0.03078113329983712633340327577884627734724381258470389700835445434525191333120095773324362027578405989 430.0 0.0307811332998371263334032757790063838049410599640486394683911138904409851600019467532034009494206969 440.0 0.03078113329983712633340327577904488186715540380336401476959781823881187402104769164927591486810667157 450.0 0.0307811332998371263334032757790541803031655494381843335328043681735029253314294195556147201727034552 460.0 0.03078113329983712633340327577905643593000168252586849066769265485558869798903957031292566069956790136 470.0 0.03078113329983712633340327577905698541891306151740400178696901495187761396910444310992237726027409048 480.0 0.03078113329983712633340327577905711983012850041458342725379411239290998736435549900510098193158324317 490.0 0.03078113329983712633340327577905715284046137901772681056407668687405251351418317905397606424547046497
In [ ]:
Questions:
How to get more digits?
- Euler Maclaurin formula to get the error of the integral approximation
- convergence acceleration on the original sum
- original product does not converge if infinite (see here)
In [ ]: