generalized_multiple_zeta_values
#
Generalized multiple zeta values
TODO:
turns warning for convergence check (ie off/warning/error)
introduce linear combinations of divergent multiple zeta values that we possibly simplify at the end.
- exception surface_dynamics.misc.generalized_multiple_zeta_values.DivergentZetaError#
- surface_dynamics.misc.generalized_multiple_zeta_values.Z2(a, b, c, check_convergence=True)#
- surface_dynamics.misc.generalized_multiple_zeta_values.Z3(a, b, c, d, e, f, g, check_convergence=True)#
The function
Z3(a,b,c,d,e,f,g)
.\[\sum_{x,y,z \geq 1} x^{-a} y^{-b} z^{-c} (x+y)^{-d} (x+z)^{-e} (y+z)^{-f} (x+y+z)^{-g}\]The reduction algorithm was designed by Bill Allombert.
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import Z3 sage: M = Multizetas(QQ) sage: Z3(1,1,1,1,1,1,1) 21/2*ζ(1,1,5) + 9/2*ζ(1,2,4) - 3/2*ζ(1,3,3) - 3/2*ζ(1,4,2) + 9/2*ζ(1,6) sage: Z3(3,0,0,0,0,3,0) 6*ζ(1,4) - 12*ζ(1,5) + 3*ζ(2,3) - 6*ζ(2,4) + ζ(3,2) - 2*ζ(3,3) sage: assert Z3(2,3,4,0,0,0,0) == M((2,)) * M((3,)) * M((4,)) sage: assert Z3(1,0,0,2,0,0,3) == M((1,2,3)) sage: assert Z3(0,0,0,2,0,1,1) == M((4,)) / 2 sage: assert Z3(1,0,1,1,0,0,1) == 3 * M((1,1,2)) sage: assert Z3(0,0,0,0,0,0,4) == 1/2 * M((2,)) - 3/2 * M((3,)) + M((4,)) sage: assert Z3(0,0,0,2,0,1,1) == 2 * M((1,1,2)) - M((2,2)) - 3 *M((1,3))
- surface_dynamics.misc.generalized_multiple_zeta_values.Z3_sort_abc(a, b, c, d, e, f, g)#
- surface_dynamics.misc.generalized_multiple_zeta_values.apply_relation(n, den_tuple, i, relation)#
- surface_dynamics.misc.generalized_multiple_zeta_values.clean_term(n, den_tuple)#
- surface_dynamics.misc.generalized_multiple_zeta_values.convergent_multizeta(t)#
Multizeta value at a convergent index
t
.TESTS:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import convergent_multizeta sage: assert all(convergent_multizeta(t) == Multizeta(*t) for t in [(2,),(3,),(1,2),(3,2),(1,1,2)]) sage: convergent_multizeta((0,3)) ζ(2) - ζ(3) sage: convergent_multizeta((0,2,2)) ζ(1,2) - ζ(2,2) sage: convergent_multizeta((1,0,3)) ζ(1,2) - ζ(1,3) - ζ(2) + ζ(3) sage: convergent_multizeta((0,1,3)) -ζ(1,3) + ζ(2) - ζ(3) sage: convergent_multizeta((0, 4)) ζ(3) - ζ(4) sage: convergent_multizeta((-1, 5)) 1/2*ζ(3) - 1/2*ζ(4) sage: convergent_multizeta((-2, 5)) 1/3*ζ(2) - 1/2*ζ(3) + 1/6*ζ(4) sage: convergent_multizeta((-1, 3, 4)) 1/2*ζ(1,4) - 1/2*ζ(2,4) sage: convergent_multizeta((-1, -1, 8)) 1/8*ζ(4) - 5/12*ζ(5) + 3/8*ζ(6) - 1/12*ζ(7) sage: convergent_multizeta((-2,-2,10)) 1/18*ζ(4) - 4/15*ζ(5) + 31/72*ζ(6) - 1/4*ζ(7) + 1/72*ζ(8) + 1/60*ζ(9) sage: convergent_multizeta((-1, -2, 10)) 1/10*ζ(5) - 3/8*ζ(6) + 5/12*ζ(7) - 1/8*ζ(8) - 1/60*ζ(9) sage: convergent_multizeta((4,-4,10)) -1/3*ζ(1,10) + 1/30*ζ(3,10) + 1/5*ζ(4,5) - 1/2*ζ(4,6) + 1/3*ζ(4,7) - 1/30*ζ(4,9) - 1/10*ζ(8) - 2/5*ζ(9) + 1/2*ζ(10) sage: convergent_multizeta((2,-1,8)) -1/2*ζ(1,8) + 1/2*ζ(2,6) - 1/2*ζ(2,7) - 1/2*ζ(7) + 1/2*ζ(8) sage: convergent_multizeta((0,3,2)) ζ(2,2) - ζ(3,2)
Divergent cases:
sage: convergent_multizeta((0,2)) Traceback (most recent call last): ... DivergentZetaError: divergent multizeta value (0, 2) sage: convergent_multizeta((0,0,3)) Traceback (most recent call last): ... DivergentZetaError: divergent multizeta value (0, 0, 3) sage: convergent_multizeta((1,0,2)) Traceback (most recent call last): ... DivergentZetaError: divergent multizeta value (1, 0, 2) sage: convergent_multizeta((0,1,2)) Traceback (most recent call last): ... DivergentZetaError: divergent multizeta value (0, 1, 2)
- surface_dynamics.misc.generalized_multiple_zeta_values.handle_term(n, den_tuple)#
Entry point for reduction of generalized multiple zeta values into standard multiple zeta values.
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import handle_term, is_convergent sage: M = Multizetas(QQ) sage: V1 = FreeModule(ZZ, 1) sage: v = V1((1,)); v.set_immutable() sage: dt = ((v,3),) sage: assert is_convergent(1, dt) and handle_term(1, dt) == M((3,)) sage: V2 = FreeModule(ZZ, 2) sage: va = V2((1,0)); va.set_immutable() sage: vb = V2((0,1)); vb.set_immutable() sage: vc = V2((1,1)); vc.set_immutable() sage: dt = ((va,2), (vc,3)) sage: assert is_convergent(2, dt) and handle_term(2, dt) == M((2,3)) sage: dt1 = ((va,2),(vb,3)) sage: dt2 = ((va,3),(vb,2)) sage: assert is_convergent(2,dt1) and is_convergent(2,dt2) sage: assert handle_term(2, ((va,2), (vb,3))) == handle_term(2, ((va,3), (vb,2))) == M((2,)) * M((3,)) sage: V3 = FreeModule(ZZ, 3) sage: va = V3((1,0,0)); va.set_immutable() sage: vb = V3((0,1,0)); vb.set_immutable() sage: vc = V3((0,0,1)); vc.set_immutable() sage: vd = V3((1,1,0)); vd.set_immutable() sage: ve = V3((1,0,1)); ve.set_immutable() sage: vf = V3((0,1,1)); vf.set_immutable() sage: vg = V3((1,1,1)); vg.set_immutable() sage: assert handle_term(3, ((va,2), (vd,3), (vg,4))) == M((2,3,4)) sage: assert handle_term(3, ((va,2), (vb,3), (vc,4))) == handle_term(3, ((va,3), (vb,2), (vc,4))) # optional mzv sage: assert handle_term(3, ((va,2), (vb,3), (vc,4))) == M((2,)) * M((3,)) * M((4,)) sage: assert handle_term(3, ((va,1), (vc,2), (vd,3))) == handle_term(3, ((va,1), (vb,2), (ve,3))) sage: assert handle_term(3, ((va,1), (vc,2), (vd,3))) == handle_term(3, ((va,2), (vb,1), (vf,3))) sage: assert handle_term(3, ((va,1), (vc,2), (vd,3))) == M((2,)) * M((1,3))
- surface_dynamics.misc.generalized_multiple_zeta_values.has_term_sum_of_smaller_terms(n, den_tuple)#
Look for a vector v_i and {v_j} with sum(v_j) = v_i
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import linear_forms, has_term_sum_of_smaller_terms sage: va,vb,vc = linear_forms(2) sage: has_term_sum_of_smaller_terms(2, ((va,1),(vb,1),(vc,1))) [2, [1, 1, 0]] sage: assert has_term_sum_of_smaller_terms(2, ((va,1),(vc,1))) is None sage: assert has_term_sum_of_smaller_terms(2, ((vb,1),(vc,1))) is None sage: assert has_term_sum_of_smaller_terms(2, ((va,1),(vb,1))) is None sage: va, vb, vd, vc, ve, vf, vg = linear_forms(3) sage: has_term_sum_of_smaller_terms(3, ((va,1),(vb,1),(vc,1),(vg,1))) [3, [1, 1, 1, 0]] sage: has_term_sum_of_smaller_terms(3, ((va,1),(vb,1),(ve,1),(vf,1),(vg,1))) [4, [0, 1, 1, 0, 0]] sage: has_term_sum_of_smaller_terms(3, ((vd,1),(ve,1),(vf,1),(vg,1))) [3, [1/2, 1/2, 1/2, 0]] sage: has_term_sum_of_smaller_terms(3, ((va,1),(vd,1),(ve,1),(vg,1))) [3, [-1, 1, 1, 0]]
- surface_dynamics.misc.generalized_multiple_zeta_values.is_Z2_convergent(a, b, c)#
TESTS:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import is_Z2_convergent
Convergent examples:
sage: assert is_Z2_convergent(1,1,1) sage: assert is_Z2_convergent(2,2,0) sage: assert is_Z2_convergent(0,0,3)
Divergent examples:
sage: assert not is_Z2_convergent(0,0,2) sage: assert not is_Z2_convergent(1,2,0)
- surface_dynamics.misc.generalized_multiple_zeta_values.is_Z3_convergent(a, b, c, d, e, f, g)#
TESTS:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import is_Z3_convergent
Convergent examples:
sage: assert is_Z3_convergent(2,0,2,2,0,0,0) sage: assert is_Z3_convergent(0,0,0,0,0,0,4) sage: assert is_Z3_convergent(1,0,0,1,0,0,2) sage: assert is_Z3_convergent(0,1,0,1,0,0,2) sage: assert is_Z3_convergent(0,1,0,0,0,1,2)
Divergent examples:
sage: assert not is_Z3_convergent(0,0,0,1,1,1,0) sage: assert not is_Z3_convergent(0,0,0,0,0,0,3)
- surface_dynamics.misc.generalized_multiple_zeta_values.is_convergent(n, den)#
TESTS:
sage: import itertools sage: from surface_dynamics.misc.generalized_multiple_zeta_values import is_convergent sage: V = FreeModule(ZZ, 3) sage: va = V((1,0,0)); va.set_immutable() sage: vb = V((0,1,0)); vb.set_immutable() sage: vc = V((0,0,1)); vc.set_immutable() sage: vd = V((1,1,0)); vd.set_immutable() sage: ve = V((1,0,1)); ve.set_immutable() sage: vf = V((0,1,1)); vf.set_immutable() sage: vg = V((1,1,1)); vg.set_immutable() sage: gens = [va,vb,vc,vd,ve,vf,vg] sage: N = 0 sage: for p in itertools.product([0,1,2], repeat=7): ....: if sum(map(bool,p)) == 3 and is_convergent(3, list(zip(gens,p))): ....: print(p) ....: N += 1 (0, 0, 0, 0, 1, 1, 2) (0, 0, 0, 0, 1, 2, 1) (0, 0, 0, 0, 1, 2, 2) (0, 0, 0, 0, 2, 1, 1) (0, 0, 0, 0, 2, 1, 2) (0, 0, 0, 0, 2, 2, 1) (0, 0, 0, 0, 2, 2, 2) (0, 0, 0, 1, 0, 1, 2) ... (2, 0, 2, 0, 0, 2, 0) (2, 0, 2, 2, 0, 0, 0) (2, 1, 0, 0, 0, 0, 2) (2, 1, 0, 0, 0, 2, 0) (2, 2, 0, 0, 0, 0, 2) (2, 2, 0, 0, 0, 2, 0) (2, 2, 0, 0, 2, 0, 0) (2, 2, 2, 0, 0, 0, 0) sage: print(N) 125
- surface_dynamics.misc.generalized_multiple_zeta_values.is_multizeta(n, den_tuple)#
Return the corresponding zeta values if it is one.
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import linear_forms, is_multizeta sage: va, vb, vd, vc, ve, vf, vg = linear_forms(3) sage: is_multizeta(3, ((vb, 2), (vf, 5), (vg, 2))) ζ(2,5,2) sage: is_multizeta(3, [(vg, 5)]) 1/2*ζ(3) - 3/2*ζ(4) + ζ(5) sage: assert is_multizeta(3, ((va,3), (vb,3), (vc,3))) is None sage: assert is_multizeta(3, ((vb,2), (ve,5), (vg,1))) is None
- surface_dynamics.misc.generalized_multiple_zeta_values.is_product(n, den_tuple)#
INPUT:
n
- number of variablesden_tuple
- tuple of pairs(vector, power)
TESTS:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import is_product sage: is_product(3, [((1,0,0),2), ((0,1,0),5), ((1,1,0),1), ((0,0,1),3)]) [(2, (((0, 1), 5), ((1, 0), 2), ((1, 1), 1))), (1, (((1), 3),))] sage: is_product(3, [((1,0,0),2), ((0,1,0),3), ((1,0,1),1), ((0,0,1),5)]) [(2, (((0, 1), 5), ((1, 0), 2), ((1, 1), 1))), (1, (((1), 3),))] sage: is_product(3, [((1,1,1),3)]) is None True
- surface_dynamics.misc.generalized_multiple_zeta_values.is_reducible(n, den_tuple)#
If (x1+x2+…+xd) is not present, use a linear relation to create it. Then try to kill using other forms. Then try to write as P(x1,x2,…,x_{d-1}) (x1+x2+…+xd) and recurse.
Should solve all Tonrheim
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import linear_forms, is_reducible sage: va, vb, vd, vc, ve, vf, vg = linear_forms(3) sage: is_reducible(3, ((va,3),(vb,3),(vc,3))) [(1, (((0, 0, 1), 3), ((0, 1, 1), 3), ((1, 1, 1), 3))), (1, (((0, 1, 0), 3), ((0, 1, 1), 3), ((1, 1, 1), 3))), (3, (((0, 0, 1), 2), ((0, 1, 1), 4), ((1, 1, 1), 3))), (3, (((0, 1, 0), 2), ((0, 1, 1), 4), ((1, 1, 1), 3))), ... (90, (((1, 0, 0), 1), ((1, 0, 1), 1), ((1, 1, 1), 7))), (90, (((0, 1, 0), 1), ((1, 1, 0), 1), ((1, 1, 1), 7))), (90, (((1, 0, 0), 1), ((1, 1, 0), 1), ((1, 1, 1), 7)))]
- surface_dynamics.misc.generalized_multiple_zeta_values.is_stufflisable(n, den_tuple)#
- surface_dynamics.misc.generalized_multiple_zeta_values.kill_relation(n, den_tuple, i, relation)#
Make calls to apply_relation until we get rid of term
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import kill_relation sage: V = FreeModule(ZZ, 2) sage: va = V((1,0)); va.set_immutable() sage: vb = V((0,1)); vb.set_immutable() sage: vc = V((1,1)); vc.set_immutable() sage: den_tuple = ((va,2),(vb,3),(vc,4)) sage: kill_relation(2, den_tuple, 2, [1,1,0]) [(1, (((0, 1), 3), ((1, 1), 6))), (2, (((0, 1), 2), ((1, 1), 7))), (1, (((1, 0), 2), ((1, 1), 7))), (3, (((0, 1), 1), ((1, 1), 8))), (3, (((1, 0), 1), ((1, 1), 8)))]
- surface_dynamics.misc.generalized_multiple_zeta_values.linear_form(R, v)#
- surface_dynamics.misc.generalized_multiple_zeta_values.linear_forms(d)#
EXAMPLES:
sage: from surface_dynamics.misc.generalized_multiple_zeta_values import linear_forms sage: list(linear_forms(1)) [(1)] sage: list(linear_forms(2)) [(1, 0), (0, 1), (1, 1)] sage: L1 = list(linear_forms(3)) sage: L2 = L1[:] sage: L2.sort(key = lambda x: x[::-1]) sage: assert L1 == L2
- surface_dynamics.misc.generalized_multiple_zeta_values.negative_rays(n)#
- surface_dynamics.misc.generalized_multiple_zeta_values.stuffle(n, den_tuple, i, j)#
- surface_dynamics.misc.generalized_multiple_zeta_values.to_Z2(den_tuple)#
- surface_dynamics.misc.generalized_multiple_zeta_values.to_Z3(den_tuple, sort=True)#
- surface_dynamics.misc.generalized_multiple_zeta_values.write_array(level, explanations)#