お花見
q級数の展開計算をしたい、と思ったのです。以前に書いた記事
ではMaximaのプログラムで書いてみたわけですが、1000次まで正確に計算するのに8分もかかるという遅さでした。これでは色々と試すには遅すぎですよね。
そこで今回全部Lisp言語で書き直してみました。
さすがにそのプログラムが見たい、という人は少ないと思うので、ここには掲載しません。100行弱のCommon Lispのプログラムです。CLOS (Common Lisp Object System)を使って書いてあります。
まずは読み込んで、
(%i1) load("/Users/yasube/Dropbox/qexpand.lisp")$
(%i2) powerdisp:true;
$$ \tag{%o2} \mathbf{true} $$
tsujimotterさんの記事で出てきた無限積をやってみます。
(%i3) q*product((1-q^n)*(1-q^(23*n)),n,1,inf);
$$ \tag{%o3} q\, \prod_{n=1}^{\infty }{\left(1-q^{n}\right)\,\left(1-q^{23\,n}\right)} $$
今回作ったqsexpand()が高速展開関数です。第１引数には展開したいqの無限積を、第２引数には求めたい最高次数を指定します。
(%i4) qsexpand(%,1000);
$$ \tag{%o4} q-q^2-q^3+q^6+q^8-q^{13}-q^{16}+q^{23}-q^{24}+q^{25}+q^{26}+q^{27}-q^{29}-q^{31}+q^{39}-q^{41}-q^{46}-q^{47}+q^{48}+q^{49}-q^{50}-q^{54}+q^{58}+2\,q^{59}+q^{62}+q^{64}-q^{69}-q^{71}-q^{73}-q^{75}-q^{78}-q^{81}+q^{82}+q^{87}+q^{93}+q^{94}-q^{98}+2\,q^{101}-q^{104}-2\,q^{118}+q^{121}+q^{123}-q^{127}-q^{128}-q^{131}+q^{138}-q^{139}+q^{141}+q^{142}+q^{146}-q^{147}+q^{150}-q^{151}+q^{162}-q^{163}+2\,q^{167}+2\,q^{173}-q^{174}-2\,q^{177}-q^{179}+q^{184}-q^{186}-q^{192}-q^{193}-q^{197}+q^{200}-2\,q^{202}+q^{208}+2\,q^{211}+q^{213}+q^{216}+q^{219}+2\,q^{223}-q^{232}-q^{233}-q^{239}-q^{242}-q^{246}-q^{248}+q^{254}-q^{257}+q^{262}-q^{269}+2\,q^{271}-q^{277}+q^{278}-q^{282}+q^{289}+q^{294}-q^{299}+q^{302}-2\,q^{303}+2\,q^{307}-q^{311}+q^{312}+2\,q^{317}-q^{325}+q^{326}-q^{328}-q^{331}-2\,q^{334}-2\,q^{346}+2\,q^{347}-q^{349}-q^{351}-q^{353}+2\,q^{354}+q^{358}+q^{361}-q^{363}-q^{368}-q^{376}+q^{377}+q^{381}+q^{384}+q^{386}+q^{392}+q^{393}+q^{394}-q^{397}-q^{400}+q^{403}-q^{409}+q^{417}-2\,q^{422}-q^{426}-q^{432}-q^{438}-q^{439}-q^{443}-2\,q^{446}+2\,q^{449}+q^{453}-q^{461}+2\,q^{463}+q^{464}+q^{466}+2\,q^{472}+q^{478}-q^{487}+q^{489}-q^{491}+q^{496}-q^{499}-2\,q^{501}-q^{509}+q^{512}+q^{514}-2\,q^{519}+q^{529}+q^{533}+q^{537}+q^{538}-q^{541}-2\,q^{542}-q^{547}-q^{552}+q^{554}-q^{568}+q^{575}-q^{577}-q^{578}+q^{579}-q^{584}-q^{587}+q^{591}+2\,q^{593}+q^{598}+2\,q^{599}-q^{600}-q^{601}+2\,q^{606}+2\,q^{607}+q^{611}-2\,q^{614}+q^{621}+q^{622}-q^{624}+q^{625}-2\,q^{633}-2\,q^{634}-q^{637}-q^{647}-q^{648}+q^{650}-q^{653}+q^{656}+q^{662}-q^{667}-2\,q^{669}-q^{673}+q^{675}-q^{683}+2\,q^{691}-2\,q^{694}+q^{696}+q^{698}+q^{699}+q^{702}+q^{706}-q^{713}+q^{717}+2\,q^{719}-q^{722}-q^{725}+q^{726}+q^{729}-q^{739}+q^{744}+q^{752}-q^{754}-q^{761}-q^{762}-2\,q^{767}+q^{771}-q^{775}-q^{783}-q^{784}-q^{786}+q^{794}-q^{806}+q^{807}+2\,q^{808}+2\,q^{809}-q^{811}-2\,q^{813}+q^{818}+2\,q^{821}-q^{823}+2\,q^{829}+q^{831}-q^{832}-q^{834}-q^{837}+2\,q^{853}-q^{857}-q^{859}-q^{863}-q^{867}+2\,q^{877}+q^{878}+2\,q^{883}+q^{886}-q^{887}+q^{897}-2\,q^{898}+q^{899}-q^{906}-2\,q^{921}+q^{922}+q^{923}-2\,q^{926}-q^{929}+q^{933}-q^{943}-2\,q^{944}-q^{947}+q^{949}-2\,q^{951}-q^{967}+q^{968}+q^{974}+q^{975}-q^{978}+q^{982}+q^{984}+2\,q^{991}+q^{993}+2\,q^{997}+q^{998} $$
(%i5) time(%);
$$ \tag{%o5} \left[ 0.18 \right] $$
たった0.18秒で展開できました。maximaで書いたプログラムでは850秒位かかっていましたから実に4700倍くらいの高速化が出来ました！
次にラマヌジャンの保型形式を100次まで展開してみます。
(%i6) q*product((1-q^n)^24,n,1,inf);
$$ \tag{%o6} q\, \prod_{n=1}^{\infty }{\left(1-q^{n}\right)^{24}} $$
(%i7) qsexpand(%,100);
$$ \tag{%o7} q-24\,q^2+252\,q^3-1472\,q^4+4830\,q^5-6048\,q^6-16744\,q^7+84480\,q^8-113643\,q^9-115920\,q^{10}+534612\,q^{11}-370944\,q^{12}-577738\,q^{13}+401856\,q^{14}+1217160\,q^{15}+987136\,q^{16}-6905934\,q^{17}+2727432\,q^{18}+10661420\,q^{19}-7109760\,q^{20}-4219488\,q^{21}-12830688\,q^{22}+18643272\,q^{23}+21288960\,q^{24}-25499225\,q^{25}+13865712\,q^{26}-73279080\,q^{27}+24647168\,q^{28}+128406630\,q^{29}-29211840\,q^{30}-52843168\,q^{31}-196706304\,q^{32}+134722224\,q^{33}+165742416\,q^{34}-80873520\,q^{35}+167282496\,q^{36}-182213314\,q^{37}-255874080\,q^{38}-145589976\,q^{39}+408038400\,q^{40}+308120442\,q^{41}+101267712\,q^{42}-17125708\,q^{43}-786948864\,q^{44}-548895690\,q^{45}-447438528\,q^{46}+2687348496\,q^{47}+248758272\,q^{48}-1696965207\,q^{49}+611981400\,q^{50}-1740295368\,q^{51}+850430336\,q^{52}-1596055698\,q^{53}+1758697920\,q^{54}+2582175960\,q^{55}-1414533120\,q^{56}+2686677840\,q^{57}-3081759120\,q^{58}-5189203740\,q^{59}-1791659520\,q^{60}+6956478662\,q^{61}+1268236032\,q^{62}+1902838392\,q^{63}+2699296768\,q^{64}-2790474540\,q^{65}-3233333376\,q^{66}-15481826884\,q^{67}+10165534848\,q^{68}+4698104544\,q^{69}+1940964480\,q^{70}+9791485272\,q^{71}-9600560640\,q^{72}+1463791322\,q^{73}+4373119536\,q^{74}-6425804700\,q^{75}-15693610240\,q^{76}-8951543328\,q^{77}+3494159424\,q^{78}+38116845680\,q^{79}+4767866880\,q^{80}+1665188361\,q^{81}-7394890608\,q^{82}-29335099668\,q^{83}+6211086336\,q^{84}-33355661220\,q^{85}+411016992\,q^{86}+32358470760\,q^{87}+45164021760\,q^{88}-24992917110\,q^{89}+13173496560\,q^{90}+9673645072\,q^{91}-27442896384\,q^{92}-13316478336\,q^{93}-64496363904\,q^{94}+51494658600\,q^{95}-49569988608\,q^{96}+75013568546\,q^{97}+40727164968\,q^{98}-60754911516\,q^{99} $$
(%i8) time(%);
$$ \tag{%o8} \left[ 0.01 \right] $$
このくらいの速度になれば色々と実験してみたくなります。
