ore_algebra.examples.iint

Iterated integrals

After J. Ablinger, J. Blümlein, C. G. Raab, and C. Schneider, Iterated Binomial Sums and their Associated Iterated Integrals, Journal of Mathematical Physics 55(11), 2014.

Thanks to Jakob Ablinger and Clemens Raab for the data.

sage: from ore_algebra.examples.iint import *

sage: dop = diffop([f[1], w[3]]) # (A.15)
sage: dop.local_basis_monomials(1)
[log(x - 1), 1, sqrt(x - 1)]
sage: iint_value(dop, [0, 0, 4*i])
[-4.934802200544679...]
sage: n(pi^2/2)
4.934802200544...

sage: dop = diffop([f[1/4], w[1], f[1]]) # (A.8)
sage: dop.local_basis_monomials(1)
[1, x - 1, (x - 1)^(3/2)*log(x - 1), (x - 1)^(3/2)]
sage: iint_value(dop, [0, 0, 16/9*i, -16/27*i*(-3*i*pi+8)])
[-3.445141853366646...]
sage: n(-pi^3/9)
-3.445141853366...

sage: dop = diffop([w[29], w[8]]) # (A.23)
sage: dop.local_basis_monomials(1)
[log(x - 1), 1, x - 1]
sage: iint_value(dop, [0, 0, -4/3])
[1.562604825792972...]
sage: n(4/9*(psi(1,1/3)-2/3*pi^2))
1.562604825792...

sage: dop = diffop([h[0], w[8], w[8], f[1], f[1]]) # (A.69)
sage: dop.local_basis_monomials(1)
[log(x - 1),
1,
x - 1,
1/2*(x - 1)^2*log(x - 1)^2,
(x - 1)^2*log(x - 1),
(x - 1)^2]
sage: myini = [0, 0, 0, 1/3, -2/3-i*pi/3, 11/12+2*i*pi/3-pi^2/6]
sage: iint_value(dop, myini)
[0.97080469562493...]
sage: iint_value(dop, myini, algorithm="binsplit")
[0.97080469...]
sage: iint_value(dop, myini, 1e-500) # long time (3.4 s)
[0.97080469562493...0383420...]

Known exact values:

sage: psi1_a = RBF(10.09559712542709408179200409989251636051890, rad=1e-35)
sage: Li_a = RBF(0.67795750683172255251567721037981473544402459, rad=1e-35)
sage: Li_b = RBF(-0.4484142069236462024430644059157743208342699, rad=1e-35)
sage: Li_c = RBF(-0.4725978446588968746186231931265476736649961, rad=1e-35)
sage: Li_d = RBF(0.34882786115484008421425193538699285536624537, rad=1e-35)
sage: Li_e = RBF(0.51747906167389938633075816189886294562237747, rad=1e-35)

sage: ref = dict()
sage: ref[1] = 4*catalan
sage: ref[2] = pi/2*zeta(2)
sage: ref[3] = -2*pi*ln(2)
sage: ref[4] = pi*(2*ln(2)^2 + zeta(2))
sage: ref[5] = 2*pi*ln(2)
sage: ref[6] = pi*(4*ln(2)^2 - zeta(2))
sage: ref[7] = 0
sage: ref[8] = -2*pi/3*zeta(2)
sage: ref[9] = -pi*(2/9*zeta(3) - 4/3*ln(2)*zeta(2)
....:               + 2*pi/(9*sqrt(3))*(4*zeta(2) - psi1_a))
sage: ref[10] = -2*pi/3*zeta(3)
sage: ref[11] = -2*pi*ln((sqrt(5)-1)/2)
sage: ref[12] = 2*pi*(2*Li_a - 6/5*zeta(3) - 6/5*ln((sqrt(5)-1)/2)*zeta(2)
....:                 + 2/3*ln((sqrt(5)-1)/2)^3)
sage: ref[13] = -pi*zeta(2)
sage: ref[14] = -2*pi*(4/5*zeta(3) + 9/5*ln((sqrt(5)-1)/2)*zeta(2)
....:                  - 2/3*ln((sqrt(5)-1)/2)^3)
sage: ref[15] = 3*zeta(2)
sage: ref[16] = zeta(2)/3
sage: ref[17] = 2*arccot(sqrt(7))^2
sage: ref[18] = 1/4*hypergeometric([1]*5, [3/2, 2, 2, 2], 1/8)
sage: ref[19] = (1/3*ln(3/2)^3 + ln(3/2)*zeta(2) + ln(3/2)*Li_b
....:            - Li_c - 2*Li_d)
sage: ref[20] = -1/3*ln(2)^2 + 4/9*zeta(2) - 2/3*Li_b
sage: ref[21] = 1/2*ln(2)^2 + zeta(2) + Li_b
sage: ref[22] = 1/24*ln(2)^4 - ln(2)^2*zeta(2) - 4/5*zeta(2)^2 + Li_e
sage: ref[23] = 4/9*(psi1_a - 4*zeta(2))
sage: ref[24] = sqrt(2)*(2/3*zeta(2) - 2*Li_b - ln(2)^2)

sage: for k in sorted(word.keys()): # long time (16 s)
....:     dop = diffop(word[k])
....:     val = iint_value(dop, ini[k])
....:     ok = "ok" if k in ref and val in RBF(ref[k]).add_error(1e-10) else ""
....:     print("(A.{})\t{}\t{}".format(k, val, ok))
(A.1)   [3.66386237670887...]  ok
(A.2)   [2.58385639002498...]  ok
(A.3)   [-4.3551721806072...]  ok
(A.4)   [8.18648809789096...]  ok
(A.5)   [4.35517218060720...]  ok
(A.6)   [0.86983785563201...]  ok
(A.7)   [+/- ...e-1...]        ok
(A.8)   [-3.4451418533666...]  ok
(A.9)   [8.38881875897492...]  ok
(A.10)  [-2.5175820907753...]  ok
(A.11)  [3.02354306885557...]  ok
(A.12)  [4.9576404218637...]   ok
(A.13)  [-5.1677127800499...]  ok
(A.14)  [2.44339103708582...]  ok
(A.15)  [4.93480220054467...]  ok
(A.16)  [0.54831135561607...]  ok
(A.17)  [0.26117239648121...]  ok
(A.18)  [0.25268502737327...]  ok
(A.19)  [0.28230892870993...]  ok
(A.20)  [0.86987360746446...]  ok
(A.21)  [1.43674636688368...]  ok
(A.22)  [-2.4278628067547...]  ok
(A.23)  [1.56260482579297...]  ok
(A.24)  [2.13970244864910...]  ok
(A.27)  [-0.5342475125153...]
(A.28)  [0.14083166505781...]
(A.29)  [0.00671775899937...]
(A.30)  [9.29746030760470...]
(A.31)  [-1.2068741628578...]
(A.32)  [0.76813058387865...]
(A.34)  [0.02493993273621...]
(A.35)  [0.25758589176283...]
(A.36)  [0.01419229300628...]
(A.37)  [0.01752989857942...]
(A.38)  [0.00931804685640...]
(A.39)  [0.00269874086778...]
(A.41)  [0.00909837539280...]
(A.42)  [0.12315273017782...]
(A.43)  [0.13779105121530...]
(A.44)  [0.03130068493969...]
(A.45)  [-0.0016512243669...]
(A.46)  [0.93615701249494...]
(A.47)  [-0.3836443021183...]
(A.48)  [-3.9987993660442...]
(A.49)  [-0.0145924441289...]
(A.50)  [0.00042865207521...]
(A.51)  [0.00722568081669...]
(A.52)  [-0.0017644712040...]
(A.53)  [5.86168742725...e-5 ...]
(A.54)  [-0.0005385608328...]
(A.55)  [-0.0007156119070...]
(A.56)  [0.00319604460177...]
(A.57)  [0.00241736336932...]
(A.58)  [0.06418646378008...]
(A.59)  [0.07359150737673...]
(A.60)  [0.01646793507128...]
(A.61)  [-0.0002687223080...]
(A.62)  [-0.0003463692070...]
(A.63)  [-0.0001223539098...]
(A.64)  [-6.9281364974...e-5 ...]
(A.65)  [0.00029942509011...]
(A.66)  [0.00076351201665...]
(A.67)  [0.34902731697236...]
(A.68)  [-0.1404863030270...]
(A.69)  [-0.9708046956249...]
(A.70)  [-0.0018467475220...]

Functions

diffop(word)
get_ini_Hstar()
iint_value(dop, ini[, eps])