Skip to content

Instantly share code, notes, and snippets.

@isuruf

isuruf/mkarakoc Secret

Last active May 4, 2016 13:50
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save isuruf/03fe806581a23c61f9c8 to your computer and use it in GitHub Desktop.
Save isuruf/03fe806581a23c61f9c8 to your computer and use it in GitHub Desktop.
Download ZIP
Raw
mkarakoc
from symengine import symbols, Integer, exp, RealDouble, sympify
import numpy as np
import sympy as sym
En = symbols("x")
nmax = 40
prec = 500
d = [[0,
-(-192 + 270*exp(Integer(1)/15))*exp(Integer(-1)/15)/5,
(-2886 + 4050*exp(Integer(1)/15))*exp(Integer(-1)/15)/25,
-(-43292 + 60750*exp(Integer(1)/15))*exp(Integer(-1)/15)/125,
(-1298761 + 1822500*exp(Integer(1)/15))*exp(Integer(-1)/15)/1250,
-(-48703538 + 68343750*exp(Integer(1)/15))*exp(Integer(-1)/15)/15625,
(-8766636841 + 12301875000*exp(Integer(1)/15))*exp(Integer(-1)/15)/937500,
-(-460248434153 + 645848437500*exp(Integer(1)/15))*exp(Integer(-1)/15)/16406250,
(-110459624196721 + 155003625000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/1312500000,
-(-1864006158319667 + 2615686171875000*exp(Integer(1)/15))*exp(Integer(-1)/15)/7382812500,
(-2236807389983600401 + 3138823406250000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/2953125000000,
-(-184536609673647033083 + 258952931015625000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/81210937500000,
(-66433179482512931909881 + 93223055165625000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/9745312500000000,
-(-3238617499772505430606699 + 4544623939324218750000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/158361328125000000,
(-2720438699808904561709627161 + 3817484109032343750000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/44341171875000000000,
-(-23542257979115520245564081201 + 33035920174318359375000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/127907226562500000000,
(-146903689789680846332319866694241 + 206144141887746562500000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/266047031250000000000000,
-(-1170638778011519244210673937719733 + 1642711130667980419921875000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/706687426757812500000000,
(-778024541693809713075401755530653317 + 1091771089920873140625000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/156558445312500000000000000,
-(-1441290463487782493472181752120535269743 + 2022505944078417493007812500000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/96674839980468750000000000000,
(-864774278092669496083309051272321161845801 + 1213503566447050495804687500000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/19334967996093750000000000000000,
-(-6190997672708883892414598889790481045032439 + 8687582350700475140419921875000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/46140264536132812500000000000000,
(-5287840365160764359874116228221046163168883193 + 7420217396010052884641015625000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/13136404726757812500000000000000000,
-(-15506591870833941485330845839258217873492749963473 + 21759787513799480084209778320312500000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/12840835620405761718750000000000000000,
(-11164746147000437869438209004265916868914779973700561 + 15667047009935625660631040390625000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/3081800548897382812500000000000000000000,
-(-27544603981086606585785055109208676485809490066695463 + 38652254136354339623267369384765625000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/2534375451395874023437500000000000000000,
(-1632844123998814038405338066873890342078786571153707046641 + 2291305625203085252867289657128906250000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/50079258919582470703125000000000000000000000,
-(-330650935109759842777080958541962794270954280658625676944803 + 463989389103624763705626155568603515625000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/3380349977071816772460937500000000000000000000,
(-21365137345553712917903692705788365168277045827172736048741117 + 29980852834388061654825074667509765625000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/72807537967700668945312500000000000000000000000,
-(-30204962922276561637686345562808301256651673538165455588907754159 + 42385430694616122164508949311191931152343750000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/34310552267278940240478515625000000000000000000000,
(-54368933260097810947835422013054942261973012368697820060033957486201 + 76293775250309019896116108760145476074218750000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/20586331360367364144287109375000000000000000000000000,
-(-972367460228672388105518124464251851992978875055557166458299624272441 + 1364484826592065163526691945133371014404296875000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/122726206186805440090942382812500000000000000000000000,
(-1103195082150348309414260563028533010261125123699395767036325391901823971 + 1548070057806270294619374134114951841796875000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/46412819794282784616210937500000000000000000000000000000,
-(-375431076394290409047540547855647640054489143658950621969549484931589470131 + 526827591547196359637655760015994548661499023437500000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/5264954245413953379901428222656250000000000000000000000000,
(-3063517583377409737827930870502084742844631412257037075271523797041770076268961 + 4298913147025122294643271001730515517077832031250000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/14320675547525953193331884765625000000000000000000000000000000,
-(-34964059375503045920862254500295532391161554161629227489511956379281071522634881 + 49063682656264982710602549476272187966649169921875000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/54480830887326995844197387695312500000000000000000000000000000,
(-45710907099341876877590442199333738136655421335519474254709315603228516664328970739 + 64144309325348535291124596262663218373240283203125000000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/23742172618266711873155493164062500000000000000000000000000000000,
-(-120505378840640022918547803247993567162757854495763214003977433259011177056337249110689 + 169100435458950076161227216897445909436454696594238281250000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/20863434188301873058535389617919921875000000000000000000000000000000,
(-274752263756659252254288991405425333131087908250340127929068547830545483688448927972370921 + 385548992846406173647598054526176673515116708234863281250000000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/15856209983109423524486896109619140625000000000000000000000000000000000,
-(-127766354767603865316978585033524499111038494695110472844598649030897542096774739955355317 + 179289476005681726219272545228786449925551092462158203125000000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/2457838590386595856339383737182617187500000000000000000000000000000000,
(-96438044578587397541255435983304291929011855795869384903103060288521464774645573718302193271601 + 135327696489088566950306917138688012403805964590437011718750000000000000000000000000000000000000*exp(Integer(1)/15))*exp(Integer(-1)/15)/618392189341267517454988948275146484375000000000000000000000000000000000000]
]
c = [[0, 18, -54, 162, -486, 1458, -4374, 13122, -39366, 118098, -354294, 1062882, -3188646, 9565938, -28697814, 86093442, -258280326, 774840978, -2324522934, 6973568802, -20920706406, 62762119218, -188286357654, 564859072962, -1694577218886, 5083731656658, -15251194969974, 45753584909922, -137260754729766, 411782264189298, -1235346792567894, 3706040377703682, -11118121133111046, 33354363399333138, -100063090197999414, 300189270593998242, -900567811781994726, 2701703435345984178, -8105110306037952534, 24315330918113857602, -72945992754341572806]
]
c[0] = [sympify(i).expand().n(prec, real=True) for i in c[0]]
d[0] = [sympify(i).expand().n(prec, real=True) for i in d[0]]
d[0][0] = (-(En*exp(Integer(1)/15).n(prec, real=True) - 9*exp(Integer(1)/15).n(prec, real=True) + 12)*exp(Integer(-1)/15).n(prec, real=True)).expand()
for n in range(1, nmax+2):
c += [[]]
d += [[]]
for i in range(nmax+1 - n):
c[n] += [((i+1)*c[n-1][i+1]).expand() + d[n-1][i] + sum([(c[0][k]*c[n-1][i-k]).expand() for k in range(i+1)])]
d[n] += [((i+1)*d[n-1][i+1]).expand() + sum([(d[0][k]*c[n-1][i-k]).expand() for k in range(i+1)])]
n = nmax
t = d[n][0]*c[n-1][0]- d[n-1][0]*c[n][0]
np.set_printoptions(precision=prec)
sol = sym.poly(t, En._sympy_())
solc = sol.coeffs()
np.roots(solc)
@mkarakoc
Copy link 

c[0] = [0, 18, -54, 162, -486, 1458, -4374, 13122, -39366, 118098, -354294, 1062882, -3188646, 9565938, -28697814, 86093442, -258280326, 774840978, -2324522934, 6973568802, -20920706406, 62762119218, -188286357654, 564859072962, -1694577218886, 5083731656658, -15251194969974, 45753584909922, -137260754729766, 411782264189298, -1235346792567894, 3706040377703682, -11118121133111046, 33354363399333138, -100063090197999414, 300189270593998242, -900567811781994726, 2701703435345984178, -8105110306037952534, 24315330918113857602, -72945992754341572806]

@mkarakoc
Copy link 

d[0] = [-(En*exp(1/15) - 9*exp(1/15) + 12)*exp(-1/15),
 -(-192 + 270*exp(1/15))*exp(-1/15)/5,
 (-2886 + 4050*exp(1/15))*exp(-1/15)/25,
 -(-43292 + 60750*exp(1/15))*exp(-1/15)/125,
 (-1298761 + 1822500*exp(1/15))*exp(-1/15)/1250,
 -(-48703538 + 68343750*exp(1/15))*exp(-1/15)/15625,
 (-8766636841 + 12301875000*exp(1/15))*exp(-1/15)/937500,
 -(-460248434153 + 645848437500*exp(1/15))*exp(-1/15)/16406250,
 (-110459624196721 + 155003625000000*exp(1/15))*exp(-1/15)/1312500000,
 -(-1864006158319667 + 2615686171875000*exp(1/15))*exp(-1/15)/7382812500,
 (-2236807389983600401 + 3138823406250000000*exp(1/15))*exp(-1/15)/2953125000000,
 -(-184536609673647033083 + 258952931015625000000*exp(1/15))*exp(-1/15)/81210937500000,
 (-66433179482512931909881 + 93223055165625000000000*exp(1/15))*exp(-1/15)/9745312500000000,
 -(-3238617499772505430606699 + 4544623939324218750000000*exp(1/15))*exp(-1/15)/158361328125000000,
 (-2720438699808904561709627161 + 3817484109032343750000000000*exp(1/15))*exp(-1/15)/44341171875000000000,
 -(-23542257979115520245564081201 + 33035920174318359375000000000*exp(1/15))*exp(-1/15)/127907226562500000000,
 (-146903689789680846332319866694241 + 206144141887746562500000000000000*exp(1/15))*exp(-1/15)/266047031250000000000000,
 -(-1170638778011519244210673937719733 + 1642711130667980419921875000000000*exp(1/15))*exp(-1/15)/706687426757812500000000,
 (-778024541693809713075401755530653317 + 1091771089920873140625000000000000000*exp(1/15))*exp(-1/15)/156558445312500000000000000,
 -(-1441290463487782493472181752120535269743 + 2022505944078417493007812500000000000000*exp(1/15))*exp(-1/15)/96674839980468750000000000000,
 (-864774278092669496083309051272321161845801 + 1213503566447050495804687500000000000000000*exp(1/15))*exp(-1/15)/19334967996093750000000000000000,
 -(-6190997672708883892414598889790481045032439 + 8687582350700475140419921875000000000000000*exp(1/15))*exp(-1/15)/46140264536132812500000000000000,
 (-5287840365160764359874116228221046163168883193 + 7420217396010052884641015625000000000000000000*exp(1/15))*exp(-1/15)/13136404726757812500000000000000000,
 -(-15506591870833941485330845839258217873492749963473 + 21759787513799480084209778320312500000000000000000*exp(1/15))*exp(-1/15)/12840835620405761718750000000000000000,
 (-11164746147000437869438209004265916868914779973700561 + 15667047009935625660631040390625000000000000000000000*exp(1/15))*exp(-1/15)/3081800548897382812500000000000000000000,
 -(-27544603981086606585785055109208676485809490066695463 + 38652254136354339623267369384765625000000000000000000*exp(1/15))*exp(-1/15)/2534375451395874023437500000000000000000,
 (-1632844123998814038405338066873890342078786571153707046641 + 2291305625203085252867289657128906250000000000000000000000*exp(1/15))*exp(-1/15)/50079258919582470703125000000000000000000000,
 -(-330650935109759842777080958541962794270954280658625676944803 + 463989389103624763705626155568603515625000000000000000000000*exp(1/15))*exp(-1/15)/3380349977071816772460937500000000000000000000,
 (-21365137345553712917903692705788365168277045827172736048741117 + 29980852834388061654825074667509765625000000000000000000000000*exp(1/15))*exp(-1/15)/72807537967700668945312500000000000000000000000,
 -(-30204962922276561637686345562808301256651673538165455588907754159 + 42385430694616122164508949311191931152343750000000000000000000000*exp(1/15))*exp(-1/15)/34310552267278940240478515625000000000000000000000,
 (-54368933260097810947835422013054942261973012368697820060033957486201 + 76293775250309019896116108760145476074218750000000000000000000000000*exp(1/15))*exp(-1/15)/20586331360367364144287109375000000000000000000000000,
 -(-972367460228672388105518124464251851992978875055557166458299624272441 + 1364484826592065163526691945133371014404296875000000000000000000000000*exp(1/15))*exp(-1/15)/122726206186805440090942382812500000000000000000000000,
 (-1103195082150348309414260563028533010261125123699395767036325391901823971 + 1548070057806270294619374134114951841796875000000000000000000000000000000*exp(1/15))*exp(-1/15)/46412819794282784616210937500000000000000000000000000000,
 -(-375431076394290409047540547855647640054489143658950621969549484931589470131 + 526827591547196359637655760015994548661499023437500000000000000000000000000*exp(1/15))*exp(-1/15)/5264954245413953379901428222656250000000000000000000000000,
 (-3063517583377409737827930870502084742844631412257037075271523797041770076268961 + 4298913147025122294643271001730515517077832031250000000000000000000000000000000*exp(1/15))*exp(-1/15)/14320675547525953193331884765625000000000000000000000000000000,
 -(-34964059375503045920862254500295532391161554161629227489511956379281071522634881 + 49063682656264982710602549476272187966649169921875000000000000000000000000000000*exp(1/15))*exp(-1/15)/54480830887326995844197387695312500000000000000000000000000000,
 (-45710907099341876877590442199333738136655421335519474254709315603228516664328970739 + 64144309325348535291124596262663218373240283203125000000000000000000000000000000000*exp(1/15))*exp(-1/15)/23742172618266711873155493164062500000000000000000000000000000000,
 -(-120505378840640022918547803247993567162757854495763214003977433259011177056337249110689 + 169100435458950076161227216897445909436454696594238281250000000000000000000000000000000*exp(1/15))*exp(-1/15)/20863434188301873058535389617919921875000000000000000000000000000000,
 (-274752263756659252254288991405425333131087908250340127929068547830545483688448927972370921 + 385548992846406173647598054526176673515116708234863281250000000000000000000000000000000000*exp(1/15))*exp(-1/15)/15856209983109423524486896109619140625000000000000000000000000000000000,
 -(-127766354767603865316978585033524499111038494695110472844598649030897542096774739955355317 + 179289476005681726219272545228786449925551092462158203125000000000000000000000000000000000*exp(1/15))*exp(-1/15)/2457838590386595856339383737182617187500000000000000000000000000000000,
 (-96438044578587397541255435983304291929011855795869384903103060288521464774645573718302193271601 + 135327696489088566950306917138688012403805964590437011718750000000000000000000000000000000000000*exp(1/15))*exp(-1/15)/618392189341267517454988948275146484375000000000000000000000000000000000000]

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment