2017年2月19日日曜日

170219

Ruby


Eisenstein series(2)

Q とR で表してみた。

require 'prime'

def power0(a, n)
  return 1 if n == 0
  k = power0(a, n >> 1)
  k *= k
  return k if n & 1 == 0
  return k * a
end

# x > 0
def sigma(x, i)
  sum = 1
  pq = i.prime_division
  pq.each{|a, n| sum *= (power0(a, (n + 1) * x) - 1) / (power0(a, x) - 1)}
  sum
end

def bernoulli(n)
  ary = []
  a = []
  (0..n).each{|i|
    a << 1r / (i + 1)
    i.downto(1){|j| a[j - 1] = j * (a[j - 1] - a[j])}
    ary << a[0] # Bn = a[0]
  }
  ary
end

def E_2k(k, n)
  a = -4 * k / bernoulli(2 * k)[-1]
  b = a.denominator
  c = a.numerator
  [b] + (1..n).map{|i| c * sigma(2 * k - 1, i)}
end

# m次以下を取り出す
def mul(f_ary, b_ary, m)
  s1, s2 = f_ary.size, b_ary.size
  ary = Array.new(s1 + s2 - 1, 0)
  (0..s1 - 1).each{|i|
    (0..s2 - 1).each{|j|
      ary[i + j] += f_ary[i] * b_ary[j]
    }
  }
  ary[0..m]
end

# m次以下を取り出す
def power(ary, n, m)
  return [1] if n == 0
  k = power(ary, n >> 1, m)
  k = mul(k, k, m)
  return k if n & 1 == 0
  return mul(k, ary, m)
end

def f(i, j, n)
  mul(power(@q_ary, i, n), power(@r_ary, j, n), n)
end

def show(e, ary0, ary1)
  print "#{e} = "
  (0..ary0.size - 1).each{|i|
    j = ary1[i]
    if i == 0
      if j < 0
        print "-#{ary0[i]} * #{-ary1[i]}"
      else
        print "#{ary0[i]} * #{ary1[i]}"
      end
    else
      if j < 0
        print " - #{ary0[i]} * #{-ary1[i]}"
      else
        print " + #{ary0[i]} * #{ary1[i]}"
      end
    end
  }
  puts
end

n = 10
@q_ary = E_2k(2, n)
@r_ary = E_2k(3, n)

e_ary = E_2k(6, n)
c_ary = [441, 250]
f_ary0, f_ary1 = f(3, 0, n), f(0, 2, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i]])}
p ""
e_ary = E_2k(7, n)
c_ary = [1]
f_ary0 = f(2, 1, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i]])}
p ""
e_ary = E_2k(8, n)
c_ary = [1617, 2000]
f_ary0, f_ary1 = f(4, 0, n), f(1, 2, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i]])}
p ""
e_ary = E_2k(9, n)
c_ary = [38367, 5500]
f_ary0, f_ary1 = f(3, 1, n), f(0, 3, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i]])}
p ""
e_ary = E_2k(10, n)
c_ary = [53361, 121250]
f_ary0, f_ary1 = f(5, 0, n), f(2, 2, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i]])}
p ""
e_ary = E_2k(11, n)
c_ary = [57183, 20500]
f_ary0, f_ary1 = f(4, 1, n), f(1, 3, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i]])}
p ""
e_ary = E_2k(12, n)
c_ary = [49679091, 176400000, 10285000]
f_ary0, f_ary1, f_ary2 = f(6, 0, n), f(3, 2, n), f(0, 4, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i], f_ary2[i]])}
p ""
e_ary = E_2k(13, n)
c_ary = [392931, 265000]
f_ary0, f_ary1 = f(5, 1, n), f(2, 3, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i]])}
p ""
e_ary = E_2k(14, n)
c_ary = [489693897, 2507636250, 395450000]
f_ary0, f_ary1, f_ary2 = f(7, 0, n), f(4, 2, n), f(1, 4, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i], f_ary2[i]])}
p ""
e_ary = E_2k(15, n)
c_ary = [815806500201, 881340705000, 26021050000]
f_ary0, f_ary1, f_ary2 = f(6, 1, n), f(3, 3, n), f(0, 5, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i], f_ary2[i]])}
p ""
e_ary = E_2k(16, n)
c_ary = [764412173217, 5323905468000, 1621003400000]
f_ary0, f_ary1, f_ary2 = f(8, 0, n), f(5, 2, n), f(2, 4, n)
(0..n).each{|i| show(e_ary[i], c_ary, [f_ary0[i], f_ary1[i], f_ary2[i]])}
p ""

出力結果
691 = 441 * 1 + 250 * 1
65520 = 441 * 720 - 250 * 1008
134250480 = 441 * 179280 + 250 * 220752
11606736960 = 441 * 16954560 + 250 * 16519104
274945048560 = 441 * 396974160 + 250 * 399517776
3199218815520 = 441 * 4632858720 + 250 * 4624512480
23782204031040 = 441 * 34413301440 + 250 * 34423752384
129554448266880 = 441 * 187477879680 + 250 * 187506813312
563087459516400 = 441 * 814940600400 + 250 * 814794618960
2056098632318640 = 441 * 2975469665040 + 250 * 2975666040144
6555199353000480 = 441 * 9486467837280 + 250 * 9486668147040
""
1 = 1 * 1
-24 = -1 * 24
-196632 = -1 * 196632
-38263776 = -1 * 38263776
-1610809368 = -1 * 1610809368
-29296875024 = -1 * 29296875024
-313495116768 = -1 * 313495116768
-2325336249792 = -1 * 2325336249792
-13195750342680 = -1 * 13195750342680
-61004818143672 = -1 * 61004818143672
-240029297071632 = -1 * 240029297071632
""
3617 = 1617 * 1 + 2000 * 1
16320 = 1617 * 960 - 2000 * 768
534790080 = 1617 * 354240 - 2000 * 19008
234174178560 = 1617 * 61543680 + 2000 * 67329024
17524001357760 = 1617 * 4858169280 + 2000 * 4834170816
498046875016320 = 1617 * 137745912960 + 2000 * 137655866880
7673653657232640 = 1617 * 2120861041920 + 2000 * 2122110676224
77480203842286080 = 1617 * 21423820362240 + 2000 * 21418943158272
574226476491096000 = 1617 * 158753769048000 + 2000 * 158760815970240
3360143509958850240 = 1617 * 928983317334720 + 2000 * 928988742914304
16320498047409790080 = 1617 * 4512174992346240 + 2000 * 4512155542392960
""
43867 = 38367 * 1 + 5500 * 1
-28728 = 38367 * 216 - 5500 * 1512
-3765465144 = -38367 * 200232 + 5500 * 712152
-3709938631392 = -38367 * 85500576 - 5500 * 78097824
-493547047383096 = -38367 * 11218984488 - 5500 * 11474230824
-21917724609403728 = -38367 * 499862636784 - 5500 * 498089967984
-486272786232443616 = -38367 * 11084671590048 - 5500 * 11088580243104
-6683009405824511424 = -38367 * 152346382155072 - 5500 * 152351956669248
-64690198594597187640 = -38367 * 1474691273530920 - 5500 * 1474676091461160
-479102079577959825624 = -38367 * 10921720940625672 - 5500 * 10921529499813576
-2872821917728374840144 = -38367 * 65489246355989232 - 5500 * 65490182325115632
""
174611 = 53361 * 1 + 121250 * 1
13200 = 53361 * 1200 - 121250 * 528
6920614800 = 53361 * 586800 - 121250 * 201168
15341851377600 = 53361 * 148641600 + 121250 * 61114944
3628395292275600 = 53361 * 20400279600 + 121250 * 20946935856
251770019531263200 = 53361 * 1439038231200 + 121250 * 1443146395680
8043563916910526400 = 53361 * 46093334702400 + 121250 * 46053422547264
150465416446925500800 = 53361 * 861697555612800 + 121250 * 861726789128832
1902324110996589786000 = 53361 * 10894180752126000 + 121250 * 10894843149545520
17831242688625346952400 = 53361 * 102121497049868400 + 121250 * 102119072037503664
132000251770026451864800 = 53361 * 755966260027216800 + 121250 * 755968133350219680
""
77683 = 57183 * 1 + 20500 * 1
-552 = 57183 * 456 - 20500 * 1272
-1157628456 = -57183 * 146232 + 20500 * 351432
-5774114968608 = -57183 * 133082976 + 20500 * 89559456
-2427722831757864 = -57183 * 32170154808 - 20500 * 28689603384
-263214111328125552 = -57183 * 3378441902544 - 20500 * 3415837464144
-12109202528761173024 = -57183 * 155862776255328 - 20500 * 155926897275744
-308317316973972772416 = -57183 * 3969266446940352 - 20500 * 3967939206760128
-5091303792066668003880 = -57183 * 65538944782146360 - 20500 * 65540990858009400
-60399282006368937251976 = -57183 * 777506848190979672 - 20500 * 777517458842153496
-552000263214112485753456 = -57183 * 7105808014591457232 - 20500 * 7105797244669716432
""
236364091 = 49679091 * 1 + 176400000 * 1 + 10285000 * 1
131040 = 49679091 * 1440 - 176400000 * 288 - 10285000 * 2016
1099243323360 = 49679091 * 876960 - 176400000 * 325728 + 10285000 * 1457568
12336522153621120 = 49679091 * 292072320 + 176400000 * 11700864 - 10285000 * 411997824
9221121336284413920 = 49679091 * 57349833120 + 176400000 * 35176468896 + 10285000 * 16227967392
1562118530273437631040 = 49679091 * 6660135541440 + 176400000 * 6601058210880 + 10285000 * 6497071680960
103486260766565509822080 = 49679091 * 436536302762880 + 176400000 * 438061091013504 + 10285000 * 440015323483008
3586400651444203277717760 = 49679091 * 15172132360815360 + 176400000 * 15173572442740992 + 10285000 * 15172068869975808
77352372210526124884754400 = 49679091 * 327295477379498400 + 176400000 * 327251435243536800 + 10285000 * 327221898778968480
1161399411211600265764157280 = 49679091 * 4913576699608450080 + 176400000 * 4913611331706352224 + 10285000 * 4913597307075535008
13104001562118531372680823360 = 49679091 * 55439481453769056960 + 176400000 * 55439979246339307200 + 10285000 * 55440561879404210880
""
657931 = 392931 * 1 + 265000 * 1
-24 = 392931 * 696 - 265000 * 1032
-805306392 = -392931 * 34632 + 265000 * 48312
-20334926626656 = -392931 * 167186976 + 265000 * 171162336
-27021598569529368 = -392931 * 64422848328 - 265000 * 6444771144
-7152557373046875024 = -392931 * 11387712944304 - 265000 * 10105554483504
-682326933054044766048 = -392931 * 1037073232984608 - 265000 * 1037089473751584
-32185646871935157619392 = -392931 * 48892286706157632 - 265000 * 48959817978105408
-906694391732570450559000 = -392931 * 1378097272692189000 - 265000 * 1378102838778701640
-17229551704624797057112632 = -392931 * 26188038166214133672 - 265000 * 26186640301645703016
-240000007152557373852181392 = -392931 * 364779879415169299632 - 265000 * 364779940958775418032
""
3392780147 = 489693897 * 1 + 2507636250 * 1 + 395450000 * 1
6960 = 489693897 * 1680 - 2507636250 * 48 - 395450000 * 1776
934155393840 = 489693897 * 1224720 - 2507636250 * 392688 + 395450000 * 975888
53074158495516480 = 489693897 * 505659840 - 2507636250 * 67089216 - 395450000 * 66529344
125380214560150002480 = 489693897 * 129351117840 + 2507636250 * 37279185936 - 395450000 * 79516693488
51856040954589843756960 = 489693897 * 21060890131680 + 2507636250 * 15066490704480 + 395450000 * 9511628122080
7123493021854278627673920 = 489693897 * 2160822606183360 + 2507636250 * 2098369148842944 + 395450000 * 2031621786790848
457358042050198589771226240 = 489693897 * 134717272385473920 + 2507636250 * 134803101024250752 + 395450000 * 134911299030780288
16828247534415852672059972400 = 489693897 * 4957295423282269200 + 2507636250 * 4960096515113176080 + 395450000 * 4962883791154433040
404722169541211889603611092720 = 489693897 * 119288258695393463760 + 2507636250 * 119289357755096403984 + 395450000 * 119289719378991436368
6960000051856040955523999143840 = 489693897 * 2051465861242156554720 + 2507636250 * 2051412780505054295520 + 395450000 * 2051366007318600561120
""
1723168255201 = 815806500201 * 1 + 881340705000 * 1 + 26021050000 * 1
-171864 = 815806500201 * 936 - 881340705000 * 792 - 26021050000 * 2520
-92268782591832 = 815806500201 * 134568 - 881340705000 * 197208 + 26021050000 * 2457000
-11795091175438423776 = -815806500201 * 173988576 + 881340705000 * 180534816 - 26021050000 * 1113204960
-49536425459206569762648 = -815806500201 * 104617833048 + 881340705000 * 34731625896 + 26021050000 * 199879986600
-32012164592742919922046864 = -815806500201 * 27210540914064 - 881340705000 * 11282282306064 + 26021050000 * 4992350445936
-6332441368275869747902027488 = -815806500201 * 3910401774129888 - 881340705000 * 3475192229286624 - 26021050000 * 3054519828108000
-553385882817076320573218661312 = -815806500201 * 322823174243838912 - 881340705000 * 319729598062193088 - 26021050000 * 316433406335739840
-26594665913504249904864455466840 = -815806500201 * 15429983442476298840 - 881340705000 * 15436589476561121880 - 26021050000 * 15444821445342229080
-809501558423540417400141934830072 = -815806500201 * 469709326015243815672 - 881340705000 * 469831003553540798136 - 26021050000 * 469944493113793897080
-17186400032012164592835188704466832 = -815806500201 * 9973673112569954220432 - 881340705000 * 9973761497118317484432 - 26021050000 * 9973874479528786860432
""
7709321041217 = 764412173217 * 1 + 5323905468000 * 1 + 1621003400000 * 1
32640 = 764412173217 * 1920 + 5323905468000 * 192 - 1621003400000 * 1536
70093866303360 = 764412173217 * 1630080 - 5323905468000 * 402048 + 1621003400000 * 551808
20160859654708062720 = 764412173217 * 803228160 - 5323905468000 * 161431296 + 1621003400000 * 163854336
150525431711563807489920 = 764412173217 * 253366181760 + 5323905468000 * 20329262976 - 1621003400000 * 93387735168
151991844177246093750032640 = 764412173217 * 53205643249920 + 5323905468000 * 23865942948480 - 1621003400000 * 9709554816000
43295116458269350559666465280 = 764412173217 * 7498254194403840 + 5323905468000 * 5794392238723584 + 1621003400000 * 4142226444876288
5149788469617367127914995164160 = 764412173217 * 699684356363412480 + 5323905468000 * 671204645516954112 + 1621003400000 * 642510156233453568
323250903208723929093223124860800 = 764412173217 * 42100628403784982400 + 5323905468000 * 41947216018774335360 + 1621003400000 * 41792421673548259200
12452826654927532105218787405188480 = 764412173217 * 1614922125605880493440 + 5323905468000 * 1615253348424607402944 + 1621003400000 * 1615606968766288470528
326400000151991844177316187616303360 = 764412173217 * 42332208491309728078080 + 5323905468000 * 42337765240473386384640 + 1621003400000 * 42343208407470359036160
""

0 件のコメント:

コメントを投稿

注: コメントを投稿できるのは、このブログのメンバーだけです。