1 Introduction↩︎

Throughout the paper, \(\mathbb{N}=\{1,2,\ldots\}\). Exponent parameters may be complex. Whenever a base is complex, powers are interpreted after fixing a branch of the logarithm, and no summand is allowed to contain a zero base at a positive integer. Infinite sums are used only under explicit convergence hypotheses.

Multiple zeta functions are among the central nested sums in experimental and computational mathematics. Their classical strict-chain form is \[\zeta(s_1,\ldots,s_d) = \sum_{n_1>\cdots>n_d\ge 1} \frac{1}{n_1^{s_1}\cdots n_d^{s_d}} . \label{eq:intro-mzv}\tag{1}\] For complex exponents, the ordinary convergence region for 1 is \[\Re(s_1+\cdots+s_k)>k, \qquad 1\le k\le d. \label{eq:intro-mzv-convergence}\tag{2}\] The case \(d=1\) gives the Riemann zeta function. Higher depth introduces the strict nesting that is responsible for the stuffle algebra of finite harmonic sums and for the rich system of relations among multiple zeta values [1]–[4].

The colored version inserts one decoration at each summation level: \[\mathop{\mathrm{Li}}_{s_1,\ldots,s_d}(z_1,\ldots,z_d) = \sum_{n_1>\cdots>n_d\ge 1} \frac{z_1^{n_1}\cdots z_d^{n_d}}{n_1^{s_1}\cdots n_d^{s_d}} . \label{eq:intro-mpl}\tag{3}\] We regard 3 both as a multiple polylogarithm and as a colored multiple zeta function. Root-of-unity colors give alternating and cyclotomic variants, while general complex colors give the usual polylogarithmic functions. For finite strict-chain colored sums we write \[\mathcal{H}_{\mathbf{s},\boldsymbol{z}}(N) = \sum_{N\ge n_1>\cdots>n_d\ge 1} \prod_{j=1}^{d} z_j^{n_j} n_j^{-s_j}, \qquad \mathop{\mathrm{Li}}_{\mathbf{s}}(\boldsymbol{z})=\lim_{N\to\infty}\mathcal{H}_{\mathbf{s},\boldsymbol{z}}(N) \label{eq:intro-finite-colored-harmonic}\tag{4}\] whenever the limit exists. The uncolored specialization is \(\zeta(\mathbf{s})=\mathop{\mathrm{Li}}_{\mathbf{s}}(1,\ldots,1)\).

The affine-base and polynomial-base zeta families used below, especially the multi-letter versions in which a single summation level carries several denominator letters, were introduced in [5]. The affine letters have the shape \[(\alpha+\beta n)^{-s}, \label{eq:intro-affine-letter}\tag{5}\] and the polynomial-base letters have the shape \[P(n)^{-s}, \qquad P(n)=c_0+c_1n+\cdots+c_m n^m. \label{eq:intro-polynomial-letter}\tag{6}\] The coefficients may be real or complex, subject to the zero-free conditions needed to define the summand and to keep the asymptotic tail on the chosen branch. These classes are natural numerical containers for shifted Hurwitz-type sums, residue-class decompositions, polynomial zeta functions, Mathieu-type series, Epstein-Hurwitz-type series, and related sums arising after partial fraction decomposition or changes of variables [6]–[11]. The present paper develops certified arbitrary-precision evaluation methods for these strict-chain families and introduces their weak-star analogues, where adjacent summation indices may be equal.

Multiple zeta functions and colored multiple zeta functions occur throughout the symbolic and numerical theory of special functions. They arise in Euler sums, iterated integrals, Mellin transforms, hypergeometric expansions, and perturbative quantum field theory. Vermaseren, Blümlein, Kurth, Broadhurst, Ablinger, Schneider, Phadikar, and others developed symbolic summation tools and structural frameworks for reducing, transforming, and organizing many of these objects [5], [12]–[16]. In particular, [5] develops a general sum formalism that includes colored, affine, and polynomial-base letters. Those methods address the symbolic and structural questions: when can a nested sum be organized, transformed, or reduced to a smaller basis of known functions or constants? The present paper is concerned with the complementary numerical question: given a convergent nested zeta-type sum, how can it be evaluated to arbitrary precision with reliable control of the infinite tail?

Direct summation is usually inefficient. In depth \(d\), the strict-chain or weak-star region contains many lattice points, and convergence may be only polynomial. High-precision numerical evaluation has therefore been essential in the subject. Borwein, Bradley, and Broadhurst used extensive high-precision computations in their compendium of Euler-Zagier sums [3]. Crandall developed fast evaluation methods for multiple zeta sums [17]. Vollinga and Weinzierl gave algorithms for numerical evaluation of multiple polylogarithms and implemented them using arbitrary-precision arithmetic in GiNaC [18]. For rigorous special-function computation, Johansson used Euler-Maclaurin summation for high-precision evaluation of the Hurwitz zeta function and its derivatives, and the Arb library provided an efficient midpoint-radius framework for propagating numerical error bounds [19], [20].

The strategy developed here keeps the nested-sum structure explicit. A finite prefix is accumulated by dynamic recurrences, and the omitted infinite part is then controlled recursively. In certified regimes, acceptance is based on an analytic inequality of the form \(|Z-V|\le E\), not on agreement between two numerical truncations alone. One branch expands the tail by Euler-Maclaurin, reducing retained terms to Hurwitz tails and to Lerch-type tails in strict-disk colored cases. A second branch uses direct absolute majorants of the omitted region. The two branches are complementary: Euler-Maclaurin is often much faster for regular power-like tails, while direct absolute bounds are robust for multi-letter and branch-sensitive inputs where the reduction to suffix constants is less favorable.

The term certified is used conservatively. A numerical value is called certified only when it is accompanied by an analytic radius that bounds the omitted infinite tail. Boundary-color calculations that rely only on unit-modulus oscillatory cancellation are therefore separated from the certification claim unless an independent analytic remainder bound is available.

The paper is organized as follows. Section 2 defines the strict and weak-star families, including colored, affine, polynomial-base, and multi-letter levels. Section 3 gives the finite-prefix recurrences, the two tail mechanisms, and the scope of the certification claim. Section 4 describes the implementation, its public interface, method selection, and returned certification data. Section 5 concludes with limitations and future directions.

2 The Family of Generalized Sums↩︎

Let \(w_j:\mathbb{N}\to\mathbb{C}\) be a one-variable level weight. The strict affine and polynomial-base variants in this section are the affine and polynomial-base multi-letter families discussed in the introduction; the weak-star sums are the star analogues introduced here. The strict and weak-star chain sums associated with \(\mathbf{w}=(w_1,\ldots,w_d)\) are \[Z^{>}(\mathbf{w}) = \sum_{n_1>\cdots>n_d\ge 1} \prod_{j=1}^{d} w_j(n_j), \label{eq:strict-general-weight}\tag{7}\] and \[Z^{\star}(\mathbf{w}) = \sum_{n_1\ge\cdots\ge n_d\ge 1} \prod_{j=1}^{d} w_j(n_j). \label{eq:star-general-weight}\tag{8}\] The strict sum is the usual multiple-zeta chain. The weak-star sum is not computed by converting it to a linear combination of strict sums; it is a separate chain geometry with its own prefix recurrence and tail bounds.

The ordinary and colored power-letter cases are obtained from \[w_j(n)=x_j^n n^{-s_j}. \label{eq:power-level-weight}\tag{9}\] Here \(x_j=1\) for ordinary multiple zeta functions, and arbitrary complex colors \(x_j\) give colored multiple zeta functions or multiple polylogarithms. We use “strict-disk” for colors with modulus less than one and “boundary-color” for unit-modulus colors. A boundary color different from one supplies oscillatory cancellation unless it is neutralized by a cumulative color product, or by a collapsed block in the weak-star case.

The affine multi-letter case allows a product of affine denominator letters at a single level: \[w_j(n) = x_j^n \prod_{\ell=1}^{m_j} (\alpha_{j\ell}+\beta_{j\ell}n)^{-s_{j\ell}}. \label{eq:affine-multiletter-weight}\tag{10}\] There is one color \(x_j\) per level, not one color per denominator letter. The single-letter affine case is \(m_j=1\) for every \(j\). Positive real affine data are included as the special case \(\alpha_{j\ell}\ge0\) and \(\beta_{j\ell}>0\). For complex affine data one assumes \[\beta_{j\ell}\ne 0, \qquad \alpha_{j\ell}+\beta_{j\ell}n\ne0 \quad(n\in\mathbb{N}), \label{eq:affine-zero-free}\tag{11}\] so that every summand is defined.

The polynomial-base multi-letter case similarly has \[w_j(n) = x_j^n \prod_{\ell=1}^{m_j}P_{j\ell}(n)^{-s_{j\ell}}, \qquad P_{j\ell}(n)=\sum_{r=0}^{M_{j\ell}}c_{j\ell r}n^r. \label{eq:polynomial-multiletter-weight}\tag{12}\] The coefficient convention is increasing in powers of \(n\). Positive polynomial bases, \(P_{j\ell}(n)>0\) for \(n\in\mathbb{N}\), give a real-denominator subfamily. Complex polynomial bases are allowed when the leading coefficient is nonzero and no positive integer is a root: \[c_{j\ell M_{j\ell}}\ne0, \qquad P_{j\ell}(n)\ne0 \quad(n\in\mathbb{N}). \label{eq:polynomial-zero-free}\tag{13}\] Degree-one polynomial bases are affine bases.

The convergence tests are most compactly stated in terms of effective level exponents. For power and affine letters set \[\sigma_j= \begin{cases} s_j, & \text{power letter},\\[2mm] \sum_{\ell=1}^{m_j}s_{j\ell}, & \text{affine multi-letter level}, \end{cases} \label{eq:effective-power-affine-exponents}\tag{14}\] and for polynomial-base levels set \[\sigma_j=\sum_{\ell=1}^{m_j}M_{j\ell}s_{j\ell}. \label{eq:effective-polynomial-exponents}\tag{15}\] These exponents describe the power-law decay of the level weights. For uncolored strict chains, the ordinary convergence condition is \[\Re(\sigma_1+\cdots+\sigma_k)>k, \qquad 1\le k\le d. \label{eq:effective-strict-convergence}\tag{16}\] For colored strict chains, let \[X_k=x_1x_2\cdots x_k, \qquad S_k=\sigma_1+\cdots+\sigma_k, \qquad r_k=\#\{1\le j\le k: X_j=1\}. \label{eq:colored-cumulative-data}\tag{17}\] A standard sufficient-and-necessary convergence test for the colored nested series in this setting is \[|X_k|<1 \quad\text{or}\quad \bigl(|X_k|=1\;\text{and}\;\Re(S_k)>r_k\bigr), \qquad 1\le k\le d, \label{eq:colored-strict-convergence}\tag{18}\] together with exclusion of any prefix satisfying \(|X_k|>1\). The uncolored criterion 16 is recovered from 18 when all \(x_j=1\).

For weak-star chains, one must also consider the strata where adjacent indices are equal. Let \(\pi=(p_1,\ldots,p_r)\) be a composition of \(d\), and collapse successive blocks of lengths \(p_1,\ldots,p_r\). The collapsed level exponent is the sum of the \(\sigma_j\) in that block, and the collapsed color is the product of the colors in that block. The weak-star chain is convergent when every such contiguous block collapse satisfies the corresponding strict-chain criterion. In the uncolored case this means that every collapsed exponent vector satisfies 16 ; in the colored case it means that every collapsed colored vector satisfies 18 .

3 Finite Prefixes and Tail Certification↩︎

This section describes the numerical core at the mathematical level. The common pattern is to compute a finite prefix directly at the chosen working precision and then control the omitted infinite tail. The strict and weak-star prefixes use different dynamic recurrences, while the tails are handled by either Euler-Maclaurin expansion or direct absolute majorants.

3.1 Finite strict and weak-star prefixes↩︎

For a strict chain define \[H_i^{>}(N) = \sum_{N\ge m_i>\cdots>m_d\ge1} \prod_{j=i}^{d}w_j(m_j), \qquad H_{d+1}^{>}(N)=1. \label{eq:strict-prefix-suffixes}\tag{19}\] Then \[H_i^{>}(N)=H_i^{>}(N-1)+w_i(N)H_{i+1}^{>}(N-1), \qquad 1\le i\le d, \label{eq:strict-prefix-recurrence}\tag{20}\] with \(H_i^{>}(0)=0\) for \(i\le d\). Thus the strict prefix with top index at most \(N\) is \(H_1^{>}(N)\).

For a weak-star chain define \[H_i^{\star}(N) = \sum_{N\ge m_i\ge\cdots\ge m_d\ge1} \prod_{j=i}^{d}w_j(m_j), \qquad H_{d+1}^{\star}(N)=1. \label{eq:star-prefix-suffixes}\tag{21}\] The update is \[H_i^{\star}(N)=H_i^{\star}(N-1)+w_i(N)H_{i+1}^{\star}(N), \qquad 1\le i\le d. \label{eq:star-prefix-recurrence}\tag{22}\] The occurrence of \(H_{i+1}^{\star}(N)\), rather than \(H_{i+1}^{\star}(N-1)\), is precisely the weak-chain correction allowing \(m_i=m_{i+1}\).

3.2 Euler-Maclaurin tails↩︎

The basic one-variable tail is the Hurwitz tail \[\sum_{n>N}n^{-a}=\zeta(a,N+1).\] For \(a\ne1\) the Euler-Maclaurin expansion gives, for an integer \(K\ge1\), \[\zeta(a,N+1) = \frac{(N+1)^{1-a}}{a-1} +\frac{(N+1)^{-a}}{2} +\sum_{j=1}^{K} \frac{B_{2j}}{(2j)!}(a)_{2j-1}(N+1)^{-a-2j+1} +\mathcal{R}_K(a,N), \label{eq:hurwitz-em-tail}\tag{23}\] where \(B_{2j}\) are Bernoulli numbers. A useful mixed form is obtained by multiplying an inner Hurwitz tail by an outer power. For \[S_{p,a}(N)=\sum_{n>N}n^{-p}\zeta(a,n+1),\] one obtains the approximation \[\begin{align} S_{p,a}(N) &\approx \frac{1}{a-1}\zeta(p+a-1,N+1) +\frac{1}{2}\zeta(p+a,N+1) \notag\\ &\quad+ \sum_{j=1}^{K} \frac{B_{2j}}{(2j)!}(a)_{2j-1} \zeta(p+a+2j-1,N+1). \label{eq:power-hurwitz-tail} \end{align}\tag{24}\] When \[\Re(a)+2K-1>0, \qquad \Re(p)+\Re(a)+2K-1>1, \label{eq:power-hurwitz-bound-conditions}\tag{25}\] the discarded part is bounded by \[|\mathcal{R}_K(p,a,N)| \le \frac{2\zeta(2K)}{(2\pi)^{2K}} \frac{|(a)_{2K}|}{\Re(a)+2K-1} \zeta(\Re(p)+\Re(a)+2K-1,N+1). \label{eq:power-hurwitz-bound}\tag{26}\] This estimate is the basic analytic certificate for ordinary power tails.

The nested recursion may be described in terms of suffix expansions. A suffix finite sum has the form \[H_{i+1}(n-1) = C_{i+1} +\sum_r c_r n^{-\alpha_r} +R_{i+1}(n), \label{eq:suffix-expansion}\tag{27}\] where \(C_{i+1}\) is the limiting suffix constant and \(R_{i+1}(n)\) is bounded by a finite list of positive power or power-logarithmic majorants. Substitution into the next level gives \[\sum_{n>N}w_i(n)H_{i+1}(n-1) = C_{i+1}\sum_{n>N}w_i(n) +\sum_r c_r\sum_{n>N}w_i(n)n^{-\alpha_r} +\sum_{n>N}w_i(n)R_{i+1}(n). \label{eq:tail-recursion-general}\tag{28}\] For power letters, the retained sums in 28 are Hurwitz or Lerch tails. For affine and polynomial-base letters, the level weight is first expanded as a finite asymptotic power series and the remainder is bounded by a positive power majorant.

For an affine letter, the large-\(n\) expansion is \[(\alpha+\beta n)^{-s}n^{-\gamma} = \beta^{-s} \sum_{r=0}^{R}\binom{-s}{r} \left(\frac{\alpha}{\beta}\right)^r n^{-s-\gamma-r} +\mathcal{E}_R(n), \label{eq:affine-binomial-expansion}\tag{29}\] valid beyond a cutoff at which the chosen logarithm branch is stable. In a multi-letter level, the finite expansions of the letters are multiplied before the one-variable tail is evaluated. The omitted product terms and the individual remainders are collected into positive power majorants.

For a polynomial base of degree \(m\), \[P(n)=c_m n^m\left(1+q_1n^{-1}+\cdots+q_mn^{-m}\right),\] and hence, beyond a branch-safe cutoff, \[P(n)^{-s}n^{-\gamma} = \sum_{r=0}^{R}a_r n^{-ms-\gamma-r}+\mathcal{E}_R(n). \label{eq:polynomial-asymptotic}\tag{30}\] The coefficients \(a_r\) come from the Taylor expansion of \((1+q_1x+\cdots+q_mx^m)^{-s}\) at \(x=0\). A Cauchy estimate on a disk in the \(x=1/n\) plane gives a computable remainder bound of the form \[|\mathcal{E}_R(n)|\le C_R n^{-\Re(ms+\gamma)-R-1}, \qquad n\ge A, \label{eq:polynomial-cauchy-majorant}\tag{31}\] provided the disk avoids zeros of the normalized polynomial factor and remains within the selected logarithm branch. Multi-letter polynomial levels are handled by multiplying the retained level expansions and summing the resulting power tails.

When an intermediate exponent in the Euler-Maclaurin recursion approaches a pole, the expression can pass through quantities such as \((a-1)^{-1}\) even though the original nested sum is convergent. A finite-part construction resolves this situation. The exponents are perturbed in a fixed generic direction, \[(s_1,\ldots,s_d) \mapsto (s_1,\ldots,s_d)+\varepsilon(2^{-1},2^{-2},\ldots,2^{-d}), \label{eq:finite-part-direction}\tag{32}\] the tail expression is expanded as a Laurent series in \(\varepsilon\), and the constant term is retained: \[\operatorname{FP}_{\varepsilon=0}F(\varepsilon)=[\varepsilon^0]F(\varepsilon). \label{eq:finite-part}\tag{33}\] The same principle applies to composite levels by perturbing at the level-block scale. In the certified Euler-Maclaurin branches, the retained finite part is accompanied by propagated remainder bounds. Parameter-stability comparisons are useful diagnostics, but they are not used as substitutes for analytic remainder bounds in the certification statements.

3.3 Direct absolute tail majorants↩︎

The direct-tail method avoids suffix constants and symbolic reductions. It uses positive majorants for the summand and then bounds the omitted chain region directly. Suppose each level admits a bound \[|w_j(n)|\le C_j q_j^n n^{-\rho_j}, \qquad 0\le q_j\le1, \label{eq:level-majorant}\tag{34}\] possibly after replacing a complex power by the elementary estimate \[|z^{-s}|\le e^{\pi|\Im s|}|z|^{-\Re s} \label{eq:complex-power-majorant}\tag{35}\] for the principal logarithm. Affine and polynomial levels satisfy 34 beyond a computable cutoff, with a finite correction for the initial range.

The suffix of a weak or strict chain with power-type majorants can be bounded recursively by a finite list of power-logarithmic terms \[C n^{\lambda}(\log n)^\ell, \qquad C\ge0, \quad \ell\in\mathbb{N}_0. \label{eq:power-log-term}\tag{36}\] The first omitted level is then bounded by sums of the form \[\sum_{n>N} n^{-\rho+\lambda}(\log n)^\ell. \label{eq:power-log-tail}\tag{37}\] When \(\rho-\lambda>1\), these tails are controlled by a finite initial correction plus the integral \[\int_N^\infty x^{-\rho+\lambda}(\log x)^\ell\,dx. \label{eq:power-log-integral}\tag{38}\] The integral has an explicit finite expression in powers of \(\log N\) divided by powers of \(\rho-\lambda-1\).

If a color satisfies \(|q|<1\), geometric decay replaces 38 by a Lerch-type tail. A basic bound is \[\sum_{n>N}|q|^n n^{-\rho} =|q|^{N+1}\Phi(|q|,\rho,N+1), \label{eq:lerch-tail-bound}\tag{39}\] and logarithmic factors are absorbed by inequalities of the form \((\log n)^\ell\le (\ell/(e\beta))^\ell n^{\beta}\) with a small positive \(\beta\). This gives certified absolute bounds in the strict disk. If \(|q|=1\), the direct-tail certificate deliberately ignores oscillatory cancellation and therefore succeeds only when the remaining power-logarithmic absolute majorant is summable. Thus an alternating example such as \(\sum_{n\ge1}(-1)^n n^{-2}\) can be certified through the absolute majorant \(\sum n^{-2}\), whereas a conditionally convergent boundary example such as \(\sum_{n\ge1}(-1)^n n^{-1/2}\) is outside the direct absolute-tail certified class. For mixed color data the geometric and power-logarithmic majorants are combined level by level; certification never depends on unproved cancellation between terms.

The direct-tail certificate is usually less asymptotically sharp than a successful Euler-Maclaurin expansion, but it is structurally robust. It works without assigning values to divergent suffix constants, applies naturally to weak-star chains, and treats multi-letter affine or polynomial levels as a single composite weight.

3.4 Scope of certification↩︎

A branch is called certified in this paper only when it produces a value \(V_N\) and a positive radius \(E_N\) satisfying \[|Z-V_N|\le E_N. \label{eq:certified-output}\tag{40}\] The radius must come from a proved bound for the omitted analytic tail: Euler-Maclaurin remainders, affine binomial remainders, polynomial Cauchy remainders, direct absolute tail majorants, or certified finite-part combinations of these bounds. Floating-point roundoff is controlled by guard precision rather than by a full directed interval-arithmetic proof, so the certificate asserted here is an analytic truncation certificate for the mathematical tail.

For colored sums, the distinction is not whether a color lies on the unit circle, but whether the selected branch proves a tail bound. Strict-disk data are certified by absolute Lerch-type majorants. Boundary-color data with non-one unit-modulus factors are also certified when oscillation can be ignored and the absolute majorant is summable, or when another explicit Euler-Maclaurin/Lerch remainder bound is available. This includes any case in which the relevant cumulative color products, and in the weak-star case the collapsed block products, still lead to a summable certified majorant.

The cases excluded from the certification claim are the conditional boundary cases in which the remaining tail is known only through cancellation from non-one unit-modulus factors. For such inputs, oscillatory Hurwitz/Lerch-type formulas may be evaluated as numerical extensions: two nearby truncation or expansion choices are compared, and the observed difference is reported as a stability indicator. This diagnostic quantity is not an \(E_N\) satisfying 40 .

4 Wolfram Language Implementation↩︎

The implementation is supplied as a Wolfram Language package file. From an examples notebook in the same directory, it can be loaded by evaluating

The package exposes twelve public numerical evaluators: six strict-chain functions and six weak-star functions.

The suffix S denotes the weak-star chain \(n_1\ge\cdots\ge n_d\ge1\). Star sums are accumulated directly; they are not converted into a sum of strict-chain values.

The input conventions are as follows.

A singleton level should use the ordinary single-letter form; a genuine multi-letter block is used when a level contains two or more denominator letters.

All public evaluators accept an optional requested decimal precision; when it is omitted, the default is 10 digits. With "ReturnCertificationData" -> True, the result is an association rather than just the numerical value. The common keys include "Value", "ErrorRadius", "N", "R", "K", "Depth", "WorkingPrecision", "Branch", and "Certification". Depending on the family, additional metadata may include "SelectedMethod", "MethodSelection", "LevelLetterCounts", "Degrees", "CoefficientDomain", and "CoefficientConvention". By default only the value stored under "Value" is returned.

The public method option has four choices.

For Automatic calls, successful data record the chosen branch through "SelectedMethod" and "MethodSelection". The "Certification" field is part of the contract: only entries described as certified analytic truncation-error data use "ErrorRadius" as an analytic tail radius, whereas entries described as “internal stability check only” use "ErrorRadius" as a stability difference and should not be used as rigorous enclosures.

The principal control options are summarized below.

For certified branches, the adaptive search uses the accepted controls \(N\), \(R\), and \(K\) according to the selected method until the analytic error radius is below the requested tolerance or the specified maxima are exhausted. Failure is a controlled outcome: it means that the selected branch and parameter bounds did not prove the requested precision, not necessarily that the mathematical series is divergent. For boundary-color diagnostic branches, the reported radius has the separate non-certified meaning described in Section 3.4.

The ordinary multiple zeta evaluator first checks the convergence region 2 . The Euler-Maclaurin method then chooses between a non-singular recursive Hurwitz-Euler-Maclaurin branch and a finite-part branch when intermediate resonances are detected. The direct-tail method uses an absolute power-log tail radius. The weak-star ordinary evaluator uses the corresponding weak-chain recurrence and weak-tail bounds.

The colored power-letter evaluators detect the all-colors-one case and redirect it to the appropriate zeta evaluator. Strict-disk color inputs use absolute Lerch majorants and are certified. Boundary inputs with non-one unit-modulus colors are split according to the proof available for the tail: when an absolute or analytic majorant proves the remainder, the result is certified; when convergence depends only on oscillatory cancellation, the Euler-Maclaurin/Lerch-type branch is returned only as an explicitly non-certified diagnostic. In the latter case the reported radius is a comparison estimate, not an analytic enclosure.

The affine evaluators support positive real data and complex coefficient data. For complex affine coefficients, the validator requires nonzero slopes and no positive-integer zero of \(a+bn\). The Euler-Maclaurin path uses branch-safe binomial tails when every needed suffix is convergent; otherwise the certified direct-tail fallback can be used. Multi-letter affine input is handled as a composite product at each level, with the number of letters per level recorded in the returned metadata.

The polynomial-base evaluators use increasing-power coefficient lists. Complex polynomial coefficients are supported when leading coefficients are nonzero and there are no positive-integer roots. Degree-one polynomial inputs are routed through the affine machinery when this is sharper. Higher-degree polynomial Euler-Maclaurin tails use branch-safe asymptotic expansions and Cauchy remainder bounds; the direct-tail path uses complex polynomial power-log majorants. For low requested digit counts in complex polynomial-base cases, an inexpensive low-order Euler-Maclaurin prepass is tried before the heavier parameter grid. Multi-letter polynomial levels are treated as composite level products rather than expanded into separate one-letter zeta calls.

The supplementary examples notebook exercises the full public interface: strict and weak-star sums, ordinary and colored power-letter cases, affine and polynomial-base families, complex coefficients, multi-letter levels, explicit method choices, certification-data output, and convergence-gate failures.

4.1 Sample evaluations↩︎

The following basic calls illustrate the strict-chain evaluators at 10 decimal digits.

The corresponding weak-star evaluators use the same input conventions and add the suffix S to the function name.

Options can be fixed explicitly. For example, the following call requests the Euler-Maclaurin branch and fixes both the product-expansion and shifted-tail orders.

Certification metadata can be requested when one wants the analytic radius and the selected branch rather than only the numerical value.

5 Conclusion↩︎

We presented a certified arbitrary-precision framework for a family of generalized multiple zeta functions. The family includes ordinary and colored multiple zeta values, affine and polynomial-base multi-letter strict-chain families, the weak-star analogues introduced here, and complex coefficient data under zero-free and branch-safe hypotheses.

The numerical method combines dynamic finite-prefix accumulation with two complementary analytic-tail strategies. Euler-Maclaurin tails are fast in regular situations and reduce retained terms to Hurwitz or Lerch-type one-variable objects. Direct absolute tails are often more robust for branch-sensitive, multi-letter, weak-star, or divergent-suffix inputs. The certification claim is confined to branches with analytic tail radii; boundary-color diagnostic outputs are deliberately separated from certified output.

The approach is intended as a numerical complement to symbolic reduction methods. It supports testing identities, validating conjectural reductions, performing integer-relation searches, and exploring constants for which no closed reduction is known. Future work includes certified cancellation-aware bounds for conditionally convergent unit-modulus colored tails, stronger directed roundoff control, and tighter integration with stuffle reduction and Shintani-type decompositions.

5.0.0.1 AI-assisted preparation disclosure

AI-assisted tools were used in a limited way to improve the clarity and presentation of selected parts of the manuscript. The mathematical ideas, definitions, results, proofs, computations, and overall substance are the author’s own. All AI-assisted suggestions were carefully reviewed, and the author remains fully responsible for the accuracy, originality, and final form of the article.

Declaration of competing interest↩︎

The author declares no competing interests.

References↩︎