仕事の部屋: 2025

2025年11月29日土曜日

2025年10月25日土曜日

251025

Maxima

Zeilberger's Algorithm

MaximaではZeilberger's Algorithmが実装されているので、使ってみた。

load(zeilberger)$

/* Zeilberger 出力を自動でシフトして係数リストにする関数 */
shifted_zeilberger_coeffs(F_expr, k, n) := block(
    [Z, m],
    print("----------------------------------------"),
    print("Applying Zeilberger to F =", F_expr),
    Z : Zeilberger(F_expr, k, n)[1][2],
    m : length(Z) - 1,
    map(lambda([c], factor(subst(n = n-m, c))), Z)
)$


F(a) := a^k * binomial(n,k) * binomial(2*n,k)$

/* A005809 */
Z : shifted_zeilberger_coeffs(F(1), k, n)$
print("A005809 coefficients:", Z)$

/* A026000 */
Z : shifted_zeilberger_coeffs(F(2), k, n)$
print("A026000 coefficients:", Z)$

/* A387928 */
Z : shifted_zeilberger_coeffs(F(3), k, n)$
print("A387928 coefficients:", Z)$

出力結果

---------------------------------------- 
Applying Zeilberger to F = binomial(n, k) binomial(2 n, k) 
A005809 coefficients: [3 (3 n - 2) (3 n - 1), - 2 n (2 n - 1)] 
---------------------------------------- 
                            k
Applying Zeilberger to F = 2  binomial(n, k) binomial(2 n, k) 
A026000 coefficients: [(n - 1) (2 n - 3) (10 n - 3), 
                           3        2
                      220 n  - 506 n  + 334 n - 63, - n (2 n - 1) (10 n - 13)] 
---------------------------------------- 
                            k
Applying Zeilberger to F = 3  binomial(n, k) binomial(2 n, k) 
A387928 coefficients: [16 (n - 1) (2 n - 3) (11 n - 3), 
                       3         2
                  649 n  - 1475 n  + 964 n - 180, - 2 n (2 n - 1) (11 n - 14)]

2025年9月27日土曜日

250927

PARI

A386379

a(0) = 1; a(n) = Sum_{k=0..floor((n-1)/5)} a(5*k) * a(n-1-5*k)
を満たすとき、以下が成り立つことを確認してみた。
For k=0..4, a(5*n+k) = (k+1) * binomial(6*n+k+1,n)/(6*n+k+1).

N=30;

a(n) = if(n==0, 1, sum(k=0, (n-1)\5, a(5*k) * a(n-1-5*k)));
apr(n, p, r) = r*binomial(n*p+r, n)/(n*p+r);
b(n) = apr(n\5, 6, n%5+1);

\\ aとbの値が一致するかどうかを確認
for(n=0, N, if(a(n) != b(n), print(n)));
for(n=0, N, print1(a(n),", "));

出力結果

1, 1, 1, 1, 1, 1, 2, 3, 4, 5, 6, 13, 21, 30, 40, 51, 114, 190, 280, 385, 506, 1150, 1950, 2925, 4095, 5481, 12586, 21576, 32736, 46376, 62832,

2025年8月25日月曜日

250825

PARI

besseli

Transcendental functions にあるように

(x/2)^ν/Γ(ν+1) を割った値を返すので、使用の際に注意が必要。

K=10;

my_besseli(n, x) = besseli(n, x)*(x/2)^n/n!;

for(k=0, K, print([k, besseli(k, 2*x)]));

print;

for(k=0, K, print([k, my_besseli(k, 2*x)]));

出力結果

[0, 1 + x^2 + 1/4*x^4 + 1/36*x^6 + 1/576*x^8 + 1/14400*x^10 + 1/518400*x^12 + 1/25401600*x^14 + 1/1625702400*x^16 + O(x^17)]

[1, 1 + 1/2*x^2 + 1/12*x^4 + 1/144*x^6 + 1/2880*x^8 + 1/86400*x^10 + 1/3628800*x^12 + 1/203212800*x^14 + 1/14631321600*x^16 + O(x^17)]

[2, 1 + 1/3*x^2 + 1/24*x^4 + 1/360*x^6 + 1/8640*x^8 + 1/302400*x^10 + 1/14515200*x^12 + 1/914457600*x^14 + 1/73156608000*x^16 + O(x^17)]

[3, 1 + 1/4*x^2 + 1/40*x^4 + 1/720*x^6 + 1/20160*x^8 + 1/806400*x^10 + 1/43545600*x^12 + 1/3048192000*x^14 + 1/268240896000*x^16 + O(x^17)]

[4, 1 + 1/5*x^2 + 1/60*x^4 + 1/1260*x^6 + 1/40320*x^8 + 1/1814400*x^10 + 1/108864000*x^12 + 1/8382528000*x^14 + 1/804722688000*x^16 + O(x^17)]

[5, 1 + 1/6*x^2 + 1/84*x^4 + 1/2016*x^6 + 1/72576*x^8 + 1/3628800*x^10 + 1/239500800*x^12 + 1/20118067200*x^14 + 1/2092278988800*x^16 + O(x^17)]

[6, 1 + 1/7*x^2 + 1/112*x^4 + 1/3024*x^6 + 1/120960*x^8 + 1/6652800*x^10 + 1/479001600*x^12 + 1/43589145600*x^14 + 1/4881984307200*x^16 + O(x^17)]

[7, 1 + 1/8*x^2 + 1/144*x^4 + 1/4320*x^6 + 1/190080*x^8 + 1/11404800*x^10 + 1/889574400*x^12 + 1/87178291200*x^14 + 1/10461394944000*x^16 + O(x^17)]

[8, 1 + 1/9*x^2 + 1/180*x^4 + 1/5940*x^6 + 1/285120*x^8 + 1/18532800*x^10 + 1/1556755200*x^12 + 1/163459296000*x^14 + 1/20922789888000*x^16 + O(x^17)]

[9, 1 + 1/10*x^2 + 1/220*x^4 + 1/7920*x^6 + 1/411840*x^8 + 1/28828800*x^10 + 1/2594592000*x^12 + 1/290594304000*x^14 + 1/39520825344000*x^16 + O(x^17)]

[10, 1 + 1/11*x^2 + 1/264*x^4 + 1/10296*x^6 + 1/576576*x^8 + 1/43243200*x^10 + 1/4151347200*x^12 + 1/494010316800*x^14 + 1/71137485619200*x^16 + O(x^17)]

[0, 1 + x^2 + 1/4*x^4 + 1/36*x^6 + 1/576*x^8 + 1/14400*x^10 + 1/518400*x^12 + 1/25401600*x^14 + 1/1625702400*x^16 + O(x^17)]

[1, x + 1/2*x^3 + 1/12*x^5 + 1/144*x^7 + 1/2880*x^9 + 1/86400*x^11 + 1/3628800*x^13 + 1/203212800*x^15 + 1/14631321600*x^17 + O(x^18)]

[2, 1/2*x^2 + 1/6*x^4 + 1/48*x^6 + 1/720*x^8 + 1/17280*x^10 + 1/604800*x^12 + 1/29030400*x^14 + 1/1828915200*x^16 + 1/146313216000*x^18 + O(x^19)]

[3, 1/6*x^3 + 1/24*x^5 + 1/240*x^7 + 1/4320*x^9 + 1/120960*x^11 + 1/4838400*x^13 + 1/261273600*x^15 + 1/18289152000*x^17 + 1/1609445376000*x^19 + O(x^20)]

[4, 1/24*x^4 + 1/120*x^6 + 1/1440*x^8 + 1/30240*x^10 + 1/967680*x^12 + 1/43545600*x^14 + 1/2612736000*x^16 + 1/201180672000*x^18 + 1/19313344512000*x^20 + O(x^21)]

[5, 1/120*x^5 + 1/720*x^7 + 1/10080*x^9 + 1/241920*x^11 + 1/8709120*x^13 + 1/435456000*x^15 + 1/28740096000*x^17 + 1/2414168064000*x^19 + 1/251073478656000*x^21 + O(x^22)]

[6, 1/720*x^6 + 1/5040*x^8 + 1/80640*x^10 + 1/2177280*x^12 + 1/87091200*x^14 + 1/4790016000*x^16 + 1/344881152000*x^18 + 1/31384184832000*x^20 + 1/3515028701184000*x^22 + O(x^23)]

[7, 1/5040*x^7 + 1/40320*x^9 + 1/725760*x^11 + 1/21772800*x^13 + 1/958003200*x^15 + 1/57480192000*x^17 + 1/4483454976000*x^19 + 1/439378587648000*x^21 + 1/52725430517760000*x^23 + O(x^24)]

[8, 1/40320*x^8 + 1/362880*x^10 + 1/7257600*x^12 + 1/239500800*x^14 + 1/11496038400*x^16 + 1/747242496000*x^18 + 1/62768369664000*x^20 + 1/6590678814720000*x^22 + 1/843606888284160000*x^24 + O(x^25)]

[9, 1/362880*x^9 + 1/3628800*x^11 + 1/79833600*x^13 + 1/2874009600*x^15 + 1/149448499200*x^17 + 1/10461394944000*x^19 + 1/941525544960000*x^21 + 1/105450861035520000*x^23 + 1/14341317100830720000*x^25 + O(x^26)]

[10, 1/3628800*x^10 + 1/39916800*x^12 + 1/958003200*x^14 + 1/37362124800*x^16 + 1/2092278988800*x^18 + 1/156920924160000*x^20 + 1/15064408719360000*x^22 + 1/1792664637603840000*x^24 + 1/258143707814952960000*x^26 + O(x^27)]

2025年7月27日日曜日

250727

PARI

A386558

f_k(x) = x*Sum_{j>=0} binomial((k+1)*j+k-1, j)/(j+1) * x^j

がf_k(x) * (1 - f_k(x))^k = x を満たすことを確認してみた。

N=10;
K=10;

a(n, k) = binomial((k+1)*n+k-1, n)/(n+1);
b(n, k) = my(x='x+O('x^(n+1))); x*sum(j=0, n, a(j, k)*x^j);
for(k=0, K, f=b(N, k); if(f*(1-f)^k==x, print([k, Vec(f)])));

出力結果

[0, [1, 0, 0, 0, 0, 0, 0, 0, 0, 0]]

[1, [1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862]]

[2, [1, 2, 7, 30, 143, 728, 3876, 21318, 120175, 690690]]

[3, [1, 3, 15, 91, 612, 4389, 32890, 254475, 2017356, 16301164]]

[4, [1, 4, 26, 204, 1771, 16380, 158224, 1577532, 16112057, 167710664]]

[5, [1, 5, 40, 385, 4095, 46376, 548340, 6690585, 83615350, 1064887395]]

[6, [1, 6, 57, 650, 8184, 109668, 1533939, 22137570, 327203085, 4928006512]]

[7, [1, 7, 77, 1015, 14763, 228459, 3689595, 61474519, 1048927880, 18236463245]]

[8, [1, 8, 100, 1496, 24682, 433160, 7932196, 149846840, 2898753715, 57135036024]]

[9, [1, 9, 126, 2109, 38916, 763686, 15636192, 330237765, 7141879503, 157366449604]]

[10, [1, 10, 155, 2870, 58565, 1270752, 28765650, 671650110, 16057800980, 391139588190]]

2025年6月23日月曜日

250623

PARI

Fuss–Catalan number(2)

f_k(x) = Sum_{j>=0} binomial(k*j+1, j)/(k*j+1) * x^j

がf_k(x) = 1/f_k(-x*f_k(x)^(2*k-1)) を満たすことを確認してみた。

N=10;
K=10;

a(n, k) = binomial(k*n+1, n)/(k*n+1);
b(n, k) = my(x='x+O('x^(n+1))); sum(j=0, n, a(j, k)*x^j);
for(k=0, K, f=b(N, k); if(f*subst(f, x, -x*f^(2*k-1))==1, print([k, Vec(f)])));

出力結果

[0, [1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0]]

[1, [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]]

[2, [1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796]]

[3, [1, 1, 3, 12, 55, 273, 1428, 7752, 43263, 246675, 1430715]]

[4, [1, 1, 4, 22, 140, 969, 7084, 53820, 420732, 3362260, 27343888]]

[5, [1, 1, 5, 35, 285, 2530, 23751, 231880, 2330445, 23950355, 250543370]]

[6, [1, 1, 6, 51, 506, 5481, 62832, 749398, 9203634, 115607310, 1478314266]]

[7, [1, 1, 7, 70, 819, 10472, 141778, 1997688, 28989675, 430321633, 6503352856]]

[8, [1, 1, 8, 92, 1240, 18278, 285384, 4638348, 77652024, 1329890705, 23190029720]]

[9, [1, 1, 9, 117, 1785, 29799, 527085, 9706503, 184138713, 3573805950, 70625252863]]

[10, [1, 1, 10, 145, 2470, 46060, 910252, 18730855, 397089550, 8612835715, 190223180840]]

2025年5月2日金曜日

250502

Ruby

An identity motivated by an amazing identity of Ramanujan

James Mc Laughlin氏が得た結果を確認してみた。

# a(n) = 99*a(n-1) - 99*a(n-2) + a(n-3)
def A(n, a, b, c)
  ary = [a, b, c]
  (3..n).each{|i| ary << 99*ary[-1] - 99*ary[-2] + ary[-3]}
  ary
end

n = 10

# A261004
a_ary = A(n, -3, -461, -45343)
# A269548
b_ary = A(n, -1, -233, -22961)
# A269549
c_ary = A(n,  1, -199, -19799)
# A269550
d_ary = A(n,  7,  465,  45347)
# A269551
e_ary = A(n,  5,  237,  22965)
# A269552
f_ary = A(n,  3,  203,  19803)
# A269553
p_ary = A(n, -3, -435, -42763)
# A269554
q_ary = A(n, -1, -343, -33861)
# A269555
r_ary = A(n,  7,  439,  42767)
# A269556
s_ary = A(n,  5,  347,  33865)
# A010701
t_ary = [3] * (n + 1)

a = [
  a_ary,
  b_ary,
  c_ary,
  d_ary,
  e_ary,
  f_ary,
  p_ary,
  q_ary,
  r_ary,
  s_ary,
  t_ary
]

(0..n).each{|i|
  (1..5).each{|j|
    ary = (0..10).map{|k| a[k][i]}

    terms = ary.each_with_index.map{|x, idx|
      display = x < 0 ? "(#{x})^#{j}" : "#{x}^#{j}"
      operator = idx == 0 ? "" : (idx < 6 ? " + " : " - ")
      operator + display
    }
    expression = terms.join

    # 実際の数値計算
    sum = ary.each_with_index.map{|x, idx| idx < 6 ? x ** j : -x ** j}.sum

    puts "#{expression} = #{sum}"
  }
  p ""
}

出力結果

(-3)^1 + (-1)^1 + 1^1 + 7^1 + 5^1 + 3^1 - (-3)^1 - (-1)^1 - 7^1 - 5^1 - 3^1 = 1

(-3)^2 + (-1)^2 + 1^2 + 7^2 + 5^2 + 3^2 - (-3)^2 - (-1)^2 - 7^2 - 5^2 - 3^2 = 1

(-3)^3 + (-1)^3 + 1^3 + 7^3 + 5^3 + 3^3 - (-3)^3 - (-1)^3 - 7^3 - 5^3 - 3^3 = 1

(-3)^4 + (-1)^4 + 1^4 + 7^4 + 5^4 + 3^4 - (-3)^4 - (-1)^4 - 7^4 - 5^4 - 3^4 = 1

(-3)^5 + (-1)^5 + 1^5 + 7^5 + 5^5 + 3^5 - (-3)^5 - (-1)^5 - 7^5 - 5^5 - 3^5 = 1

(-461)^1 + (-233)^1 + (-199)^1 + 465^1 + 237^1 + 203^1 - (-435)^1 - (-343)^1 - 439^1 - 347^1 - 3^1 = 1

(-461)^2 + (-233)^2 + (-199)^2 + 465^2 + 237^2 + 203^2 - (-435)^2 - (-343)^2 - 439^2 - 347^2 - 3^2 = 1

(-461)^3 + (-233)^3 + (-199)^3 + 465^3 + 237^3 + 203^3 - (-435)^3 - (-343)^3 - 439^3 - 347^3 - 3^3 = 1

(-461)^4 + (-233)^4 + (-199)^4 + 465^4 + 237^4 + 203^4 - (-435)^4 - (-343)^4 - 439^4 - 347^4 - 3^4 = 1

(-461)^5 + (-233)^5 + (-199)^5 + 465^5 + 237^5 + 203^5 - (-435)^5 - (-343)^5 - 439^5 - 347^5 - 3^5 = 1

(-45343)^1 + (-22961)^1 + (-19799)^1 + 45347^1 + 22965^1 + 19803^1 - (-42763)^1 - (-33861)^1 - 42767^1 - 33865^1 - 3^1 = 1

(-45343)^2 + (-22961)^2 + (-19799)^2 + 45347^2 + 22965^2 + 19803^2 - (-42763)^2 - (-33861)^2 - 42767^2 - 33865^2 - 3^2 = 1

(-45343)^3 + (-22961)^3 + (-19799)^3 + 45347^3 + 22965^3 + 19803^3 - (-42763)^3 - (-33861)^3 - 42767^3 - 33865^3 - 3^3 = 1

(-45343)^4 + (-22961)^4 + (-19799)^4 + 45347^4 + 22965^4 + 19803^4 - (-42763)^4 - (-33861)^4 - 42767^4 - 33865^4 - 3^4 = 1

(-45343)^5 + (-22961)^5 + (-19799)^5 + 45347^5 + 22965^5 + 19803^5 - (-42763)^5 - (-33861)^5 - 42767^5 - 33865^5 - 3^5 = 1

(-4443321)^1 + (-2250073)^1 + (-1940399)^1 + 4443325^1 + 2250077^1 + 1940403^1 - (-4190475)^1 - (-3318283)^1 - 4190479^1 - 3318287^1 - 3^1 = 1

(-4443321)^2 + (-2250073)^2 + (-1940399)^2 + 4443325^2 + 2250077^2 + 1940403^2 - (-4190475)^2 - (-3318283)^2 - 4190479^2 - 3318287^2 - 3^2 = 1

(-4443321)^3 + (-2250073)^3 + (-1940399)^3 + 4443325^3 + 2250077^3 + 1940403^3 - (-4190475)^3 - (-3318283)^3 - 4190479^3 - 3318287^3 - 3^3 = 1

(-4443321)^4 + (-2250073)^4 + (-1940399)^4 + 4443325^4 + 2250077^4 + 1940403^4 - (-4190475)^4 - (-3318283)^4 - 4190479^4 - 3318287^4 - 3^4 = 1

(-4443321)^5 + (-2250073)^5 + (-1940399)^5 + 4443325^5 + 2250077^5 + 1940403^5 - (-4190475)^5 - (-3318283)^5 - 4190479^5 - 3318287^5 - 3^5 = 1

(-435400283)^1 + (-220484321)^1 + (-190139599)^1 + 435400287^1 + 220484325^1 + 190139603^1 - (-410623923)^1 - (-325158121)^1 - 410623927^1 - 325158125^1 - 3^1 = 1

(-435400283)^2 + (-220484321)^2 + (-190139599)^2 + 435400287^2 + 220484325^2 + 190139603^2 - (-410623923)^2 - (-325158121)^2 - 410623927^2 - 325158125^2 - 3^2 = 1

(-435400283)^3 + (-220484321)^3 + (-190139599)^3 + 435400287^3 + 220484325^3 + 190139603^3 - (-410623923)^3 - (-325158121)^3 - 410623927^3 - 325158125^3 - 3^3 = 1

(-435400283)^4 + (-220484321)^4 + (-190139599)^4 + 435400287^4 + 220484325^4 + 190139603^4 - (-410623923)^4 - (-325158121)^4 - 410623927^4 - 325158125^4 - 3^4 = 1

(-435400283)^5 + (-220484321)^5 + (-190139599)^5 + 435400287^5 + 220484325^5 + 190139603^5 - (-410623923)^5 - (-325158121)^5 - 410623927^5 - 325158125^5 - 3^5 = 1

(-42664784581)^1 + (-21605213513)^1 + (-18631740599)^1 + 42664784585^1 + 21605213517^1 + 18631740603^1 - (-40236954115)^1 - (-31862177823)^1 - 40236954119^1 - 31862177827^1 - 3^1 = 1

(-42664784581)^2 + (-21605213513)^2 + (-18631740599)^2 + 42664784585^2 + 21605213517^2 + 18631740603^2 - (-40236954115)^2 - (-31862177823)^2 - 40236954119^2 - 31862177827^2 - 3^2 = 1

(-42664784581)^3 + (-21605213513)^3 + (-18631740599)^3 + 42664784585^3 + 21605213517^3 + 18631740603^3 - (-40236954115)^3 - (-31862177823)^3 - 40236954119^3 - 31862177827^3 - 3^3 = 1

(-42664784581)^4 + (-21605213513)^4 + (-18631740599)^4 + 42664784585^4 + 21605213517^4 + 18631740603^4 - (-40236954115)^4 - (-31862177823)^4 - 40236954119^4 - 31862177827^4 - 3^4 = 1

(-42664784581)^5 + (-21605213513)^5 + (-18631740599)^5 + 42664784585^5 + 21605213517^5 + 18631740603^5 - (-40236954115)^5 - (-31862177823)^5 - 40236954119^5 - 31862177827^5 - 3^5 = 1

(-4180713488823)^1 + (-2117090440081)^1 + (-1825720439399)^1 + 4180713488827^1 + 2117090440085^1 + 1825720439403^1 - (-3942810879483)^1 - (-3122168268781)^1 - 3942810879487^1 - 3122168268785^1 - 3^1 = 1

(-4180713488823)^2 + (-2117090440081)^2 + (-1825720439399)^2 + 4180713488827^2 + 2117090440085^2 + 1825720439403^2 - (-3942810879483)^2 - (-3122168268781)^2 - 3942810879487^2 - 3122168268785^2 - 3^2 = 1

(-4180713488823)^3 + (-2117090440081)^3 + (-1825720439399)^3 + 4180713488827^3 + 2117090440085^3 + 1825720439403^3 - (-3942810879483)^3 - (-3122168268781)^3 - 3942810879487^3 - 3122168268785^3 - 3^3 = 1

(-4180713488823)^4 + (-2117090440081)^4 + (-1825720439399)^4 + 4180713488827^4 + 2117090440085^4 + 1825720439403^4 - (-3942810879483)^4 - (-3122168268781)^4 - 3942810879487^4 - 3122168268785^4 - 3^4 = 1

(-4180713488823)^5 + (-2117090440081)^5 + (-1825720439399)^5 + 4180713488827^5 + 2117090440085^5 + 1825720439403^5 - (-3942810879483)^5 - (-3122168268781)^5 - 3942810879487^5 - 3122168268785^5 - 3^5 = 1

(-409667257120241)^1 + (-207453257914553)^1 + (-178901971320799)^1 + 409667257120245^1 + 207453257914557^1 + 178901971320803^1 - (-386355229235355)^1 - (-305940628162963)^1 - 386355229235359^1 - 305940628162967^1 - 3^1 = 1

(-409667257120241)^2 + (-207453257914553)^2 + (-178901971320799)^2 + 409667257120245^2 + 207453257914557^2 + 178901971320803^2 - (-386355229235355)^2 - (-305940628162963)^2 - 386355229235359^2 - 305940628162967^2 - 3^2 = 1

(-409667257120241)^3 + (-207453257914553)^3 + (-178901971320799)^3 + 409667257120245^3 + 207453257914557^3 + 178901971320803^3 - (-386355229235355)^3 - (-305940628162963)^3 - 386355229235359^3 - 305940628162967^3 - 3^3 = 1

(-409667257120241)^4 + (-207453257914553)^4 + (-178901971320799)^4 + 409667257120245^4 + 207453257914557^4 + 178901971320803^4 - (-386355229235355)^4 - (-305940628162963)^4 - 386355229235359^4 - 305940628162967^4 - 3^4 = 1

(-409667257120241)^5 + (-207453257914553)^5 + (-178901971320799)^5 + 409667257120245^5 + 207453257914557^5 + 178901971320803^5 - (-386355229235355)^5 - (-305940628162963)^5 - 386355229235359^5 - 305940628162967^5 - 3^5 = 1

(-40143210484294963)^1 + (-20328302185186241)^1 + (-17530567468999199)^1 + 40143210484294967^1 + 20328302185186245^1 + 17530567468999203^1 - (-37858869654185443)^1 - (-29979059391701841)^1 - 37858869654185447^1 - 29979059391701845^1 - 3^1 = 1

(-40143210484294963)^2 + (-20328302185186241)^2 + (-17530567468999199)^2 + 40143210484294967^2 + 20328302185186245^2 + 17530567468999203^2 - (-37858869654185443)^2 - (-29979059391701841)^2 - 37858869654185447^2 - 29979059391701845^2 - 3^2 = 1

(-40143210484294963)^3 + (-20328302185186241)^3 + (-17530567468999199)^3 + 40143210484294967^3 + 20328302185186245^3 + 17530567468999203^3 - (-37858869654185443)^3 - (-29979059391701841)^3 - 37858869654185447^3 - 29979059391701845^3 - 3^3 = 1

(-40143210484294963)^4 + (-20328302185186241)^4 + (-17530567468999199)^4 + 40143210484294967^4 + 20328302185186245^4 + 17530567468999203^4 - (-37858869654185443)^4 - (-29979059391701841)^4 - 37858869654185447^4 - 29979059391701845^4 - 3^4 = 1

(-40143210484294963)^5 + (-20328302185186241)^5 + (-17530567468999199)^5 + 40143210484294967^5 + 20328302185186245^5 + 17530567468999203^5 - (-37858869654185443)^5 - (-29979059391701841)^5 - 37858869654185447^5 - 29979059391701845^5 - 3^5 = 1

(-3933624960203786301)^1 + (-1991966160890337193)^1 + (-1717816709990600999)^1 + 3933624960203786305^1 + 1991966160890337197^1 + 1717816709990601003^1 - (-3709782870880938195)^1 - (-2937641879758617703)^1 - 3709782870880938199^1 - 2937641879758617707^1 - 3^1 = 1

(-3933624960203786301)^2 + (-1991966160890337193)^2 + (-1717816709990600999)^2 + 3933624960203786305^2 + 1991966160890337197^2 + 1717816709990601003^2 - (-3709782870880938195)^2 - (-2937641879758617703)^2 - 3709782870880938199^2 - 2937641879758617707^2 - 3^2 = 1

(-3933624960203786301)^3 + (-1991966160890337193)^3 + (-1717816709990600999)^3 + 3933624960203786305^3 + 1991966160890337197^3 + 1717816709990601003^3 - (-3709782870880938195)^3 - (-2937641879758617703)^3 - 3709782870880938199^3 - 2937641879758617707^3 - 3^3 = 1

(-3933624960203786301)^4 + (-1991966160890337193)^4 + (-1717816709990600999)^4 + 3933624960203786305^4 + 1991966160890337197^4 + 1717816709990601003^4 - (-3709782870880938195)^4 - (-2937641879758617703)^4 - 3709782870880938199^4 - 2937641879758617707^4 - 3^4 = 1

(-3933624960203786301)^5 + (-1991966160890337193)^5 + (-1717816709990600999)^5 + 3933624960203786305^5 + 1991966160890337197^5 + 1717816709990601003^5 - (-3709782870880938195)^5 - (-2937641879758617703)^5 - 3709782870880938199^5 - 2937641879758617707^5 - 3^5 = 1

(-385455102889486762703)^1 + (-195192355465067858801)^1 + (-168328507011609898999)^1 + 385455102889486762707^1 + 195192355465067858805^1 + 168328507011609899003^1 - (-363520862476677757803)^1 - (-287858925156952833301)^1 - 363520862476677757807^1 - 287858925156952833305^1 - 3^1 = 1

(-385455102889486762703)^2 + (-195192355465067858801)^2 + (-168328507011609898999)^2 + 385455102889486762707^2 + 195192355465067858805^2 + 168328507011609899003^2 - (-363520862476677757803)^2 - (-287858925156952833301)^2 - 363520862476677757807^2 - 287858925156952833305^2 - 3^2 = 1

(-385455102889486762703)^3 + (-195192355465067858801)^3 + (-168328507011609898999)^3 + 385455102889486762707^3 + 195192355465067858805^3 + 168328507011609899003^3 - (-363520862476677757803)^3 - (-287858925156952833301)^3 - 363520862476677757807^3 - 287858925156952833305^3 - 3^3 = 1

(-385455102889486762703)^4 + (-195192355465067858801)^4 + (-168328507011609898999)^4 + 385455102889486762707^4 + 195192355465067858805^4 + 168328507011609899003^4 - (-363520862476677757803)^4 - (-287858925156952833301)^4 - 363520862476677757807^4 - 287858925156952833305^4 - 3^4 = 1

(-385455102889486762703)^5 + (-195192355465067858801)^5 + (-168328507011609898999)^5 + 385455102889486762707^5 + 195192355465067858805^5 + 168328507011609899003^5 - (-363520862476677757803)^5 - (-287858925156952833301)^5 - 363520862476677757807^5 - 287858925156952833305^5 - 3^5 = 1

2025年4月28日月曜日

250428

PARI

A074143等

A(1,k) = 1; A(n,k) = k + n * Sum_{j=1..n-1} A(j,k)

が

A(n,k) = (n+2) * A(n-1, k) - (n-1) * A(n-2, k) for n > 3

を満たすことを確認してみた。

M=20;
N=10;

b(n, k) = if(n==1, 1, k+n*sum(j=1, n-1, b(j, k)));
a(n, k) = if(n<4, b(n, k), (n+2)*a(n-1, k)-(n-1)*a(n-2, k));

\\ aとbの値が一致するかどうかを確認
for(k=-5, 5, for(n=1, M, if(a(n, k)!=b(n, k), print([n, k]))));
for(k=-5, 5, print1(k,": "); for(n=1, N, print1(a(n, k),", ")); print);

出力結果

-5: 1, -3, -11, -57, -355, -2555, -20865, -190765, -1931495, -21461055,

-4: 1, -2, -7, -36, -224, -1612, -13164, -120356, -1218604, -13540044,

-3: 1, -1, -3, -15, -93, -669, -5463, -49947, -505713, -5619033,

-2: 1, 0, 1, 6, 38, 274, 2238, 20462, 207178, 2301978,

-1: 1, 1, 5, 27, 169, 1217, 9939, 90871, 920069, 10222989,

0: 1, 2, 9, 48, 300, 2160, 17640, 161280, 1632960, 18144000,

1: 1, 3, 13, 69, 431, 3103, 25341, 231689, 2345851, 26065011,

2: 1, 4, 17, 90, 562, 4046, 33042, 302098, 3058742, 33986022,

3: 1, 5, 21, 111, 693, 4989, 40743, 372507, 3771633, 41907033,

4: 1, 6, 25, 132, 824, 5932, 48444, 442916, 4484524, 49828044,

5: 1, 7, 29, 153, 955, 6875, 56145, 513325, 5197415, 57749055,

2025年3月23日日曜日

250323

PARI

Fuss–Catalan number(1)

f_k(x) = Sum_{j>=0} binomial(k*j+1, j)/(k*j+1) * x^j

がf_k(x) = 1 + x*f_k(x)^k を満たすことを確認してみた。

N=10;
K=10;

a(n, k) = binomial(k*n+1, n)/(k*n+1);
b(n, k) = my(x='x+O('x^(n+1))); sum(j=0, n, a(j, k)*x^j);
for(k=0, K, f=b(N, k); if(f==1+x*f^k, print([k, Vec(f)])));