Optimal Mittag–Lefﬂer Summation

: A novel method of an optimal summation is developed that allows for calculating from small-variable asymptotic expansions the characteristic amplitudes for variables tending to inﬁnity. The method is developed in two versions, as the self-similar Borel–Leroy or Mittag–Lefﬂer summations. It is based on optimized self-similar iterated roots approximants applied to the Borel–Leroy and Mittag–Lefﬂer- transformed series with the subsequent inverse transformations. As a result, simple and transparent expressions for the critical amplitudes are obtained in explicit form. The control parameters come into play from the Borel–Leroy and Mittag–Lefﬂer transformations. They are determined from the optimization procedure, either from the minimal derivative or minimal difference conditions, imposed on the analytically expressed critical amplitudes. After diff-log transformation, virtually the same procedure can be applied to critical indices at inﬁnity. The results are obtained for a number of various examples. The examples vary from a rapid growth of the coefﬁcients to a fast decay, as well as intermediate cases. The methods give good estimates for the large-variable critical amplitudes and exponents. The Mittag–Lefﬂer summation works uniformly well for a wider variety of examples.


Introduction
A method for determining from small-variable asymptotic expansions the characteristic exponents for variables tending to infinity, called a self-similar Borel-Leroy summation, was suggested recently [1]. It is based on optimized self-similar factor approximation applied to the Borel-Leroy transformed series. The celebrated method of Borel summation [2] was also explored in [1], but it should be noted that it does not allow for controlling the accuracy of approximations. The results can be improved only by adding more and more terms to the expansions, but even then there are some important physical cases when the technique could fail because of the poles of gamma-functions entering the formulas. In addition, the factor approximants are difficult to handle analytically in high-orders. Therefore, application of some remedial techniques is recommended in the current paper.
First, to solve some debilitating technical problems, instead of the factor approximants, we apply the so-called iterated roots. Such innovation allows for extending calculation systematically into high orders of perturbation theory.
Second, in order to broaden the spectrum of problems to be tackled, in addition to Borel and Borel-Leroy summations, we consider also Mittag-Leffler summation.
Third, optimization ideas are applied now in a consistent manner in rather high orders, as a part of the summation procedures, in order to determine control parameters introduced by the transformations.
As is commonly understood, the most often applied Padé-Borel summation has the following meaning. Suppose we know that the coefficients a n of the power series for the sought function change as n!. Let us get rid of the factorial dependence that complicates direct calculation with the original series, by performing Borel-transform. The transformation basically smoothed the coefficients by excluding n! growth. To the transformed series, one can apply Padé approximations with a justifiable hope to obtain better results in a finite domain because the transformed coefficients could be decreasing in magnitude, or grow much slower than originally. The whole table of Padé approximants could be used for resummation applied to the transformed series [3]. After the summation is accomplished, we have to get back to the original series by performing inverse Borel transformation, by means of improper integration. The last step is accomplished numerically when we are interested in the sought function at finite values of variables. It should be noted that asymptotics of the Padé approximants and of the Padé-Borel summation are limited to integer numbers [1], and more complicated cases of rational or irrational critical indices can not be covered at all. Thus, the standard approach has certain weak points, especially when a direct calculation of critical indices is required.
Application of the self-similar approximants brings some distinct advantages over the standard Padé-Borel summation [3][4][5][6]. In particular, factor and root approximants are defined uniquely for the asymptotic conditions. Most importantly, simple and transparent expressions for the critical indices and amplitudes arise just by following the definition faithfully. In [1], we suggested replacing the whole Padé table with uniquely defined factor approximants, free of a non-uniqueness of Padé tables. It turned out that already from such an approximant applied to the transformed series one can find the critical index, even without having to restore the original asymptotic behavior. However, to find the critical amplitude, one ought to perform an inverse Borel-transform and restore the original asymptotic behavior.
In addition, the self-similar Borel-Leroy summation by means of factor and root approximants leads to an exact reconstruction of the zero-dimensional partition function [1], while the Padé-based techniques are not able to do so. To understand the reasons why factor approximants are more powerful than Padé approximants, we can turn the physical case of a D-dimensional hydrogen atom with slowly convergent energy expansion into powers of 1/D [7,8]. In this case and for similar simple dependencies with power-laws, factor (as well as root) approximants [9] give exact reconstruction of the sought energy. In the case of a zero-dimensional field theory discussed in [1], the Borel-Leroy transform of the exact solution has the form of such a simple dependence with power-law, analogous to the above mentioned energy, and inverse transformation gives exact reconstruction of the zero-dimensional partition function.
Of course, one can approach the class of problems with power-law behavior at large variables from different angles. Recently, the Borel summation was applied also to the hypergeometric functions/approximants [10,11], leading after integration to the so-called hypergeometric-Meijer approximants [12,13]. Such techniques are interesting by all means but are rather difficult to implement in practice. Most notably, non-uniqueness of the solutions comes into play already at the initial stage of asserting asymptotic equivalence. Application to realistic problems thus requires a fitting procedure [14,15], and most results appear only in numerical form. Very honest and direct comments by Mera et al. [12] should be highly commended without disparaging the further attempts to utilize their technique. However, much simpler methods should not be abandoned altogether, since they are much more transparent and could be always rather easily used for estimates before trying any of the time-and effort-consuming techniques [16].
Our choice throughout the current paper of the known self-similar iterated root approximants [9] allows for analytical calculation of critical exponents and amplitudes, with explicit introduction of control parameters and their subsequent evaluation from the general conditions of non-perturbative nature. All that allows the researcher to keep the calculations under control. Again, we can cover rational or irrational critical indices, while also retaining uniqueness in determination of the approximants' coefficients. Summation methods of the current paper are also modular, allowing for introducing the desired modifications quite routinely. There are three modules involved. One is for the transformation per se; another is for the approximants and corresponding amplitudes; and the third is for the optimal selection of the parameters.
Iterated roots are much simpler than factor approximants to define from the asymptotic expansion at small variables, and the formulae can be written in rather high orders. There is no distinction now between odd and even number of terms, as is for the factors. Correspondingly, by using roots, we can avoid an uneven performance of factor approximants for an odd and even number of terms. Resummation procedure leading to the self-similar iterated roots consists of optimizing the flow of approximations viewed as a dynamic system in approximation space. Such flow evolves symmetrically, to satisfy the property of self-similarity. Moreover, simple and symmetric algorithms may be preferred by Nature by itself [17]. See also review [18] on relations between symmetry, simplicity and asymptotic methods.
The control parameter introduced by the Borel-Leroy transformation could be determined from optimization procedure, by imposing either the minimal sensitivity or the minimal difference conditions formulated for the analytically expressed critical amplitudes or critical indices [19,20]. When we apply the Borel-Leroy summation (see, e,g., refs. [3,21,22] and references therein), or yet different Mittag-Leffler summation [5,6,[22][23][24], to the initial truncated series, we again expect that the problem for such a transformed series will be simpler than the original one, and could be tackled by means of the selfsimilar approximants at hand. Applying such transformation also allows for bypassing the problem of the gamma-function poles appearing for the Borel summation, allowing now to calculate critical properties in some important physical cases, such as a shift of Bose-condensation temperature, or even in generic problems of critical indices.
Our method of optimal Mittag-Leffler summation complements the standard Mittag-Leffler summation in two principal ways. Instead of the typically employed Padé approximants, we use iterated root approximants. In order to determine the parameter of Mittag-Leffler transformation, we systematically apply the optimization conditions. The former innovation allows us to express the critical properties in analytical form. The latter gives a systematic way of finding the critical parameters. When the two are considered together, we arrive at the direct analytical method for calculation of critical amplitudes and indices from asymptotic series.

Preliminaries
Generally speaking, we are interested in a real non-negative function f (x) of a real variable x. In addition, we are in possession only of the truncated asymptotic expansion of which can be divergent for finite values of the variable x, and where a 0 > 0. For simplicity, we may set a 0 = 1 sometimes. The most difficult region for approximating is that of the large variable, where approximations are expected to be at their worst performance. We will try to solve the most difficult case of the large-variable behaviour of the function, restricting it to the following asymptotic form: The constant B will be called the critical amplitude and the power β stands for the critical exponent. After the self-similar renormalization (see, e.g., [9,25]), application to the truncated series (2), we arrive to a self-similar approximant f * k (x). Typically, and as is assumed in the current paper, for many vital physical problems, the general expressions for the expansion coefficients are not known. Thus, convergence of the sequence of approximants cannot be determined rigorously. We can discuss only numerical convergence, observing it within the finite sequences of numerical results. It is always preferred, for exceedingly obvious reasons, to apply various methods for finding large-variable exponents. There are several different ways to implement the self-similar renormalization. All of them can be found in the book [25]. The particular realization is dictated by the problem specifics. We opt for such variants that would lead to the possibility of taking explicitly the limit x → ∞. When such a limit is available in analytical form, we could find explicitly the related approximation for the critical amplitude and critical exponent. The task of finding details of critical behavior is far from trivial. In fact, the problem is central for critical phenomena, ubiquitous in Nature. Critical phenomena are of great interest in seismology, phase transitions, polymers, non-equilibrium transport phenomena in structured media, field theory, quantum mechanics, statistics, finance and so on.
Instead of a single critical index, one can generate a spectrum of critical indices dependent on conditions of numerical experiments. When the position of the critical point is known, the problem of critical indices can be reduced to calculation only of the critical indices at infinity. In turn, such a problem can be mapped to an equivalent problem of finding only critical amplitudes at infinity. The latter problem arises most naturally when the value of the exponent β is available from other sources. Then, we need to find only the critical amplitude B. Such problems, even considered on their own, are of a long-standing interest in applied mathematics [26]. However, the tables can be turned, so that, by solving the problem only for critical amplitude, we would also develop the method of calculation for the critical exponent β as explained in the review [9].
Let us denote the critical amplitude following from the approximants at the intermediate stage of the transformed series as C(u) [1]. After the inverse Borel-Leroy transformation, the sought amplitude for the original problem B(u) can be expressed as for the Borel-Leroy summation. The formula simply follows from the more general result from [1] by fixing the critical index to a constant β.
In the current paper, we are going to apply also the Mittag-Leffler transformations [22,23], and find that for the Mittag-Leffler summation. The two formulas for critical amplitudes are central for our study, since they bridge analytically the sought critical amplitude for original coefficients with the amplitude for the transformed coefficients. Only the multipliers involving gamma-function are needed to return from the transformed quantities to the original. The greatest advantage of possessing formulas is that we can perform optimization on formulas as well. Such a convenience arises due to symmetry and simplicity requirements involved in derivation of the selfsimilar approximants.
The parameter u is inherent to the Borel-Leroy and Mittag-Leffler transforms, and it has to be determined from optimization conditions, selected to mimic the conditions of convergence for infinite series, but adapted to the case of truncated series. Such conditions are grounded in some rather intuitive notions, such as a requirement that the numerical solution stop being dependent on the approximation number (minimal difference), or else stop being dependent on the artificially introduced parameter (minimal derivative), as soon as the numerical convergence of the renormalized sequences is achieved.
To find the control parameter u, it we turn to the complete amplitude optimization [1,19,20], in the form of minimal difference [1], or minimal derivative conditions imposed on amplitude B(u) [19,20]. Alternatively, we may solve at first the optimization problem formulated analogously but for the marginal amplitude C(u), and find the control parameter u. From the solution of marginal problems, we can easily restore the solution of the full original problem.
The examples considered in the paper correspond to the coefficients a n behaving in various ways, covering situations from a rapid growth of the coefficients to a fast decay, as well as intermediate cases. The methods of optimal summation manage to adapt to such behaviors.

Iterated Root Approximants
Self-similar root approximants are important for our study. They follow directly from consideration of the flow of approximations viewed as a dynamic system in approximation space. The flow evolves by following the symmetrical trajectory of the approximation cascade. The corresponding symmetry is called self-similarity. Its application to resummation problems formulated in the space of approximations was pioneered by V.I. Yukalov. More details on root approximants and self-similarity can be found, e.g., in our recent paper [27]. The parameters inherent to root approximants are introduced by means of the algebraic transforms and bootstrapping. When the large-variable asymptotic of the sought function is known, they appear to be defined uniquely [9]. If we try to find these parameters only from the small-variable expansion, then there are multiple solutions to the system of equations for the parameters. However, if some additional condition on the parameters is imposed, such as a requirement that all terms in root approximant contribute to the large-variable amplitude [28], then we can have a unique solution for the coefficients of root approximants.
Conditioning the approximants on correct critical exponent β at infinity, we arrived in [28] to the so-called iterated root approximants, with known powers and unknown amplitudes, Now, all relevant parameters A j could be uniquely defined from the asymptotic equivalence with the truncated series in the small x limit. More up-to-date details on derivation of iterated roots can be found in [27].
In the large-variable limit, the approximant (4) behaves as the power-law with the critical amplitude The parameters of iterated root approximations (4) could be found iteratively, adding one new parameter in each order. They are defined for all k = 1, 2, . . . and the question about "odd" and " even" approximants does not arise here at all, unlike the case of factor approximants [1].
Calculations with iterated roots are really simple, but the price to be paid for such convenience amounts to the impossibility to incorporate into (4) nontrivial sub-critical indices. For that end, we should go back to the original form of root approximants, or resort to the so-called additive approximants [25].
Remarkably, critical exponents can be found by using the very same techniques as developed for critical amplitudes. Again, we deal with a function, with power-law asymptotic behavior f (x) Bx β (x → ∞). The critical exponent can be expressed as the limit where ψ(x) ≡ d dx ln f (x), as shown in [1,9,29,30]. The small-variable expansion for the original function is given by the sum f k (x), and we can find the corresponding smallvariable expression for the diff-log function ψ(x), which can be expanded in powers of x, leading to the truncated series However, for x → ∞, we find that simply where the " critical amplitude" is the sought critical index β, and the "critical index" is fixed to δ ≡ −1, to make the limit (7) existent [1,9]. Without much ado, we can simply apply the technique of iterated roots described above for the critical amplitudes and calculate the critical index β. Bearing that in mind, we are going to discuss below the case of critical amplitude B calculation. Then, it remains simply apply the methods to be developed to calculate also the critical index β. However, the case of δ = −1, or β = −1, depending on the context, appears to be divergent when treated by the self-similar Borel summation technique of [1]. In such cases, some other types of summation should be advanced. They will be considered in the next two sections, in order to be able to calculate critical properties in such important cases.

Self-Similar Borel-Leroy Summation
The generalization of the Borel summation method could be devised through the Borel-Leroy transform [1,22], where u plays the role of a control parameter. It has to be chosen to improve the convergence of the sequence of approximants to be applied to the transformed series. The transformation amounts to the coefficients of the original series a n being diminished by the factor (n + u)!. The latter expression should be understood in terms of the factorial function defined for arbitrary m = n + u.
With u = 0 in the Formula (9), we return to the Borel-transform, while, with u → ∞, the Borel-Leroy summation is supposed to be equivalent to the summation of the original series (2) with chosen approximants. The latter statement does follow from the analysis of the Formula (42) from our previous work [1].
Summing up the truncated series L k (x, u) by means of iterated self-similar root approximants yields the approximant The inverse Borel-Leroy transformation gives the approximation f * k (x) for the sought function f (x). It is relatively straightforward to observe that, at large values of the variable i.e., the self-similar Borel-Leroy transform behaves as a power-law with correct critical index β. The marginal amplitude C k (u) follows from the definition of iterated root. Finally, at large values of x, function (12) is expressed as follows: where the critical amplitude is presented as the marginal amplitude C k (u) multiplied by the correcting gamma-function.
In what follows, we are going to approach the optimization problem following the lines of the papers [19,20], by employing either the minimal difference condition or the minimal derivative condition formulated for critical amplitudes. The solution is expected to correspond to a positive u when the coefficients a n are increasing with approximation number n. Covering the case of a n increasing faster than n! may be problematic for the Bose-Leroy methodology and may be better approached with Mittag-Leffler summation.

Self-Similar Mittag-Leffler Summation
Borel summation method can be generalized in a different manner, through the Mittag-Leffler transform [22,23]. For the truncated series (2), it reads as follows: with arbitrary u playing the role of a control parameter. It has to be chosen in such a way that the convergence of the sequence of approximants could be expected to improve. Compared with the conventional Borel-transform valid for a non negative integer m = un, the Mittag-Leffler transform allows for a more general behavior of the series coefficients, formalized by means of the factorial function (10), valid for arbitrary m = nu. From the Formula (10), it is clear that Mittag-Leffler transformation amounts to the coefficients of the original series a n being divided by the factor (un)!.
With u = 1 in the Formula (17), we return to the Borel-transform, while, with u = 0, the Mittag-Leffler summation degenerates to the summation of the original truncated series (2) with chosen approximants. Summing up the truncated series M k (x, u) by means of iterated self-similar root approximants results in the following approximation: The inverse Mittag-Leffler transformation gives the self-similar approximation for the sought function At large values of the variable, the self-similar Mittag-Leffler transform behaves as with the amplitude Finally, it remains to discover that, at large values of x, function (19) acquires the form of a power-law where the amplitude can be expressed through a simple formula, Formula (10) for the gamma-function was employed again to derive the Formula (23). One can think that actual performance of the Mittag-Leffler summation method described in this section depends on how well it performs at the two rather close limits, u = 0 and u = 1. We expect the solution to the optimization problem [19], to be located at positive u when the coefficients a n are increasing with approximation number n. Negative optimal values could also come into play, in particular when a n are decreasing with n.

Bose-Einstein Condensation Temperature
The Bose-Einstein condensation temperature T 0 of ideal uniform Bose gas in threedimensional space is rather well known and can be found in standard textbooks. However, it is also well known that the ideal Bose gas is unstable below the condensation temperature. Atomic interactions do stabilize the system, while also shifting the transition temperature by the amount ∆T c ≡ T c − T 0 . This shift, at a asymptotically small gas parameter Γ ≡ ρ 1/3 a s , where a s is atomic scattering length, and ρ stands for the gas density, is widely assumed to be proportional to Γ, ∆T c T 0 c 1 Γ (Γ → 0) .
Monte Carlo simulations (see, e.g., Refs. [31,32] and multiple references therein) estimate that c 1 = 1.3 ± 0.05. The coefficient c 1 can be defined [33][34][35] as the strongcoupling limit of a function c 1 (g) with the formally obtained expansion in the small parameter g, with g playing the role of effective coupling. For such a problem, we have β = −1, and Borel summation according to the scheme of [1] is not applicable. However, the problem can be addressed with different summations. Let us first apply the method of Mittag-Leffler summation. In the highest available order, we obtain from the minimal-difference condition, B 4 (u) = B 3 (u) [19], the unique solution u * = 0.538618, B 4 = 1.33967.
The result appears to be in a good agreement with the known results mentioned above.
In the same way, one can find the values of c 1 for the O(1) field theory [34], with the following, formally obtained expansion available, and for the problem with β = −1, we calculate from the minimal-difference condition B 4 (u) = B 3 (u), the unique solution The solution for critical amplitude agrees well with Monte Carlo numerical estimate c 1 = 1.09 ± 0.09 , (see [31] and references therein).
For the O(4) field theory, the following expansion can be found in [34], and for the problem with β = −1, from the minimal-difference condition in the form we find the unique solution The solution agrees well with Monte Carlo numerical estimate c 1 = 1.6 ± 0.1 , as discussed in [9,31].
The method of Borel-Leroy summation could be applied similarly and gives rather reasonable results shown in Table 1. In all cases, the solutions appear to be defined uniquely. Results of calculations in the fourth-order by the two summation methods are compared in Table 1. They appear to be in a good agreement with each other and with Monte Carlo simulations. The problems of condensation temperature correspond to slowly decaying, relatively small numbers of the coefficients a n . Rather surprisingly, the summation methods developed for the series with rapidly growing coefficients in mind work rather well. The available truncated series are quite short but nonetheless informative; the minimal difference condition manages to mimic plausibly the role of missing, high-order coefficients.

Fluctuating Fluid Membrane
The pressure p exerted by a two-dimensional membrane can be calculated from perturbation theory with the wall stiffness g playing the role of a small parameter [36]. Extensive calculations produce the following expansion p k (g) = π 2 8g 2 1 + k ∑ n=1 a n g n , (g → 0), with the coefficients The rigid-wall limit of g → ∞ is of primary interest. The coefficients a n are rapidly decaying and the problem, at first glance, does not seem to be treatable accurately by means of Borel-type summations.
In order to apply the techniques developed above, we consider the quantity which diverges at infinity with critical index β = 2. The sought pressure can be calculated as the critical amplitude at infinity. The problem of calculating p(∞) appears to be rather difficult, judging by the fact that the Padé approximants fail to produce reasonable results [9]. Factor approximants also essentially fail, giving very high estimates [9]. Iterated roots application improves the results. For instance, in the case of the iterated roots, we have p * 6 (∞) = 0.071. Calculations with various Borel-type techniques conducted with all known coefficients a n , do not bring any improvement over the iterated root results. In such a situation, one would like to have more coefficients in the expansion in the hope to improve results. In order to have some proxy for the true coefficients, we add some trial, " very small" coefficients, a 7 = a 8 = 0. Routinely, one constructs the approximant M * 8 (g, u), and finds u by employing the minimal-derivative condition [19] in the form Solving the latter equation, we find negative control parameter u * = −0.253341, with the final Mittag-Leffler summation result for the sought pressure treated as the critical amplitude, p 8 (u * ) = 0.077289.
The quality of the result is the same as of the full, self-similar root approximants considered in [38]. The result can be compared also with the value of 0.0821 from the variational perturbation theory [36], which overestimates the Monte Carlo result by about the same amount.
Borel summation gives a much smaller value of 0.06517 for the pressure, while, for the Borel-Leroy summation, the solution to a minimal derivative condition is absent altogether even in the considered high orders. The solution (p ≈ 0.055) exists though in the sixth order and is vastly inferior. The result p * * 3 (∞) = 0.0792, obtained by the analyticalnumerical method of doubly renormalized iterated roots [9], remains the best to the best of our knowledge.

Unitary Limit
The function f * k (x) can be reconstructed from M * k (x) even without integration, by adapting the method of self-similarly corrected Padé approximants [1,31]. To this end, one can use the following formula: function and concentrate on finding the parameters of the correcting, diagonal Padé approximant P n/n (x). Asymptotic equivalence in the formula (28) is enforced by equating terms of the same order of the small-variable expansion of (28) and of the expansion (2). The case of β = 0 can be most conveniently handled by means of the self-similar factor approximants [1]. Even without a transformation, the approximants, when applied to the initial series (2), allow for a direct definition of characteristic amplitudes by extrapolating the initial series (2) to the form where N k = k/2 , k = 4, 6, . . .. However, let us first apply the Mittag-Leffler transform to the truncated series (1), and apply then the factor approximants so that Remember that the condition on the exponent has to be satisfied in each and every order. "Marginal " amplitude for the transformed series is given by the expression suggested in [1], and the critical amplitude is expressed through the marginal amplitude and correcting gamma-function, In the case of β = 0 to be discussed below, B k (u) = C k (u). The amplitudes could be optimized by means of the minimal derivative condition [19], From the latter equation, we have to find the control parameter u = u * . The technique just discussed above could be applied to the problem of extrapolation of the ground-state energy of a dilute Fermi gas. The truncated series for the energy follows from the perturbation theory [39,40]. They are derived with respect to the effective coupling parameter g ≡ |k F a s | , expressed through a Fermi wave number k F , and atomic scattering length a s . The perturbation theory leads to the expansion up to the fourth order, with the coefficients From the experimental viewpoint, it is very important that the effective coupling parameter could be varied by controlling the scattering length, in a rather broad band. The band includes the case of g → ∞. The latter limit corresponds to the so-called unitary Fermi gas, used to model the low-density neutron matter [41]. The often used Bertsch parameter is defined through the energy, ξ = E(∞) × 10 3 [41]. The energy of a dilute Fermi gas in the limit-case g → ∞ can be calculated by various methods. In particular, the Monte Carlo calculations [41,42] yield E(∞) = 0.1116, while the experimentally best known value is E(∞) = 0.1128, (see [41,43]). More recent Monte Carlo simulations give E(∞) = 0.1164 [44].
We construct the approximant M * 4 (g, u) in the highest available order, and find control parameter u by employing the minimal-derivative condition Solving the latter equation, we find negative control parameter with the final Mittag-Leffler-summation result for the sought critical amplitude, Borel-summation can be accomplished by simply setting u = 1, bringing a much smaller result By adding one more trial coefficient a 5 = 0, we can estimate E 6 (1) = 0.1451. The bounds given by the Borel-summation are quite reasonable, but are rather wide. Using all four terms from the expansion and one condition at infinity, we calculate, from the Formula (28) with the control parameter (34), the final value of the energy The latter estimate is closer to the result of [44]. Of course, one can study the unitary limit also with the Borel-Leroy transform, but, in this case, the factor approximants give vastly inferior results E * 4 (∞) = 0.174. We conclude that various Mittag-Leffler summation results are in better agreement with the best available numerical and experimental results than other methods presented in the review [9]. The expressions of the type of (28) can be used for interpolation as described in [27], to construct a ground state energy dependence of low-density neutron matter as a function of interaction strength, possibly generalizing the known formula of [45]. Overall, the optimal Mittag-Leffler summation applied to the series (33), with the coefficients behaving non-monotonously, first decaying and then increasing in absolute values, allows for a rather nimble adaptation to such a situation.

Critical Exponent: Anharmonic Models
We are going to implement below the idea sketched in the end of Section 3. The main target is the critical index at infinity. The Mittag-Leffler transform has to be applied now to the truncated series (8), with the formal result where arbitrary u plays the role of a control parameter. Summing up the latter series by means of iterated self-similar root approximants yields Accomplishing the inverse Mittag-Leffler transformation gives the approximation for the sought function At large values of the variable, the self-similar Mittag-Leffler transform behaves as with the amplitude Finally, at large values of x, the function (19) acquires the following form: where the "amplitude" from (39) is the sought critical index To obtain the final Formula (41), we simply have to adapt the Formula (23). Following equivalent steps, one can also accomplish the Borel-Leroy summation in application to critical indices. Omitting here self-evident intermediate steps, we arrive to the core expression for the critical index: The critical indices β k (u) could be optimized by means, e.g., of the minimal difference condition [19], with k = 2, 3, . . .. From the latter equation, one expects to find the optimal parameter u = u * . The quantum model of one-dimensional quartic oscillator serves mostly for evaluation of the new methods' applicability to various quantum field-theoretical problems [6]. With the Hamiltonian given in dimensionless units as follows: where x ∈ (−∞, ∞) and g ≥ 0, one can calculate the ground-state energy E(g) for the corresponding quantum-mechanical problem. The expansion of E(g) in powers of the coupling constant g results in a divergent series of type (2). Explicit value of the coefficients a n can be found in [46,47]. The leading term of the strong-coupling behavior also has the form a power-law Let us perform the diff-log transformation, find coefficients d n in the formula (8), and proceed step-by-step to the minimal difference condition in the form where the highest available order of truncation where we managed to perform calculations analytically is k = 10.
In the case of Mittag-Leffler summation, we find a unique solution to the minimal difference problem, so that u * = 0.557476, β 10 = 0.330762, while, for the Borel-Leroy summation, we obtained u * = 0.945554, β 10 = 0.324663.
The two estimates (especially the former) are in a good agreement with exact critical exponents. The former estimate is also the best compared with the whole family of estimates obtained in [1].
Let us also bring here more details on the self-similar Mittag-Leffler summation. With increasing approximation order, it leads to the following results: Already in the low orders quite good accuracy is achieved. The results for some "odd" approximations are not shown here, since there are no solutions to the minimal-difference condition in such cases. Numerical convergence appears to be monotonous in the case of self-similar Borel-Leroy summation, brought up below for the higher-orders, Let us consider another popular touchstone, the partition function of the zero-dimensional anharmonic model with the non-negative coupling parameter g. Expanding the integral in powers of g leads to the divergent series (2) with the coefficients a n = (−1) n √ π n! Γ 2n + 1 2 . The strong-coupling behavior of the integral is known as well, with the leading term shown below: We follow the same procedure as in previous example, and write down the minimal difference condition β 10 (u) = β 9 (u).
The two estimates are in a good agreement with exact critical exponent β = −1/4. The two estimates are also the best compared with the whole family of estimates obtained in [1]. In the latter case, there is a small caveat of two more solutions with slightly inferior index, and only the best solution is brought up.
With increasing approximation order, the self-similar Mittag-Leffler summation gives good results already in low orders Let us also bring here more details on the self-similar Borel-Leroy summation dependence on the order of approximation. The summation leads to the following approximants: The results for the approximations with odd k are not shown here, since there are multiple solutions to the minimal-difference condition in such cases. Note that the control parameters are found to be positive for the two touchstone examples with the fast growing coefficients.

Critical Exponent: Expansion Factor of Polymer Chain
Let us proceed with a concrete example of a certain theoretical interest. The excluded volume effect in a polymer chain has been one of the central problems in the field of polymer solutions. The excluded volume effect emerges due to the effective repulsion between segments of the polymer chain, resulting in the enlargement of the chain. The expansion factor α is quantified by the ratio of the mean square end-to-end distance of the chain to its unperturbed value. The expansion factor of polymers gives a parsimonious example of self-organized criticality. The power-laws emerge without an external tuning, in the absence of a or phase transition and for a pure geometrical quantity. However, its statistical mechanics could be developed as well for a zero-component order parameter.
A perturbation theory for the expansion factor leads to a series in a single dimensionless interaction parameter g, quantifying the expansion factor α(g) [48,49]. For the three-dimensional case, as g → 0, the expansion factor can be presented as the truncated series of the type (2), with the coefficients The strong-coupling behavior of the expansion factor as g → ∞, has been found numerically in the paper [50], α(g) ∝ g 0.3508 , while a slightly lower result α(g) ∝ g 0.3504 was obtained in [51]. Let us perform the diff-log transformation, find coefficients d n in the (8), and proceed to the minimal difference condition in the highest available order of truncation k = 5.
In the case of Mittag-Leffler summation, we found two very close positive solutions to the minimal difference problem and took their average for the critical index, while there half-difference serves to measure a margin of error, so that β 5 = 0.350364 ± 0.000777139, while for the Borel-Leroy summation we obtained The two estimates (especially the former) are in a good agreement with numerical estimates quoted above, and smaller than the number 0.3544 quoted in [49]. They are also significantly better than our previous estimates obtained by various methods in [1]. In Table 2, we present the critical index dependencies on the approximation order for the two types of summation. They are obtained similarly to the cases discussed above. Consider also the two-dimensional polymer coil [48], with the expansion factor The critical index at infinity is considered to be known exactly, with β = 1/2, and where the critical amplitude at infinity B is not known exactly [52,53]. The amplitude was estimated in the book [25] by several methods, with a best guess of B ≈ 1.
We apply for optimization the minimal difference condition, The Mittag-Leffler summation brings the following results Correspondingly, the Borel-Leroy summation brings The Borel summation can be performed just by setting u = 0 in the formulas for Borel-Leroy technique, and it gives for the amplitude a very close result, B 4 = 0.969559. All results appear to be only slightly smaller than unity.

Critical Exponent: Susceptibility of the Ising Model
Consider a three-dimensional Ising model for spins S z j = ±1/2, characterized by the HamiltonianĤ defined on a simple cubic lattice, with the ferromagnetic interaction of nearest neighbors with the strength J. The dimensionless interaction parameter g ≡ J k B T is defined as the inverse of the temperature T.
Let us perform the diff-log transformation as explained above and find the coefficients d n in the formula (8). Then, we proceed step-by-step to the minimal difference condition in the form γ 10 (u) − γ 9 (u) = 0, in the 11th order for the original series.
The two estimates (especially the former) are in a good agreement with numerical estimates quoted above. In Table 3, we also present the critical index γ dependencies on the approximation order for the two types of summation. They are obtained similarly to the cases of highest orders discussed above in more detail. Only the results for rather high orders are presented.

Compressibility Factor of Hard-Spheres Fluids
The thermodynamics of hard-spheres fluids is described by the compressibility factor where P is pressure, ρ is density, T is temperature, and a s is the sphere diameter. The parameter f is called packing fraction. The compressibility factor is believed to diverge as a power-law at a finite critical point as with the parameter f c = 1 and conjectured β = 3 [57,58]. This behavior is infused into various phenomenological equations of state. For low packing fraction, the compressibility factor is represented by the virial expansion and more exact terms, up to the eleventh order inclusively are available [59][60][61]. The expansion is also presented and discussed in the book [25].
In order to reduce the consideration to the same type of problems as treated above, we make the same substitution as in the paper [1], with f c = 1: Then, as x → ∞ ( f → f c − 0), the compressibility factor at the critical point behaves as Now, we can apply the same methods as above, imposing the minimal difference condition, β 10 (u) − β 9 (u) = 0, corresponding to the highest available order of truncation k = 11, in the original truncated series (2). In the case of Mittag-Leffler summation, we found two solutions to the minimal difference problem, u * 1 = 0.285573, β 10,1 = 3.55362, u * 2 = 0.46858, β 10,2 = 2.6812. Thus, we take their average to produce the final estimate of the critical index, while their half-difference serves to measure a margin of error, giving β 10 = 3.11741 ± 0.436213.
The two summations lead to differing conclusions about the value of index β. The former result better agrees with the conjectured value of β = 3.
Note that the Mittag-Leffler summation starts generating the predictions, which could be interpreted in favor of the estimate β ≈ 3 already in the eighth-order of perturbation theory, with the averages β 8 = 2.9406 ± 0.39745, β 9 = 3.15744 ± 0.477721.
In the former case, there are three solutions and variance is taken as a measure of margin of error. One can hope that when and if a few more terms in the expansion become available, the conjecture will be confirmed even more reliably.

Conclusions
Our results for the critical indices are based on a few simple formulas and ideas well worth discussing in conclusion.
First, we perform the diff-log transformation of the original series. After the Mittag-Leffler transformation applied to the diff-log transformed series, we can find the part of a critical index S(u), which follows from the self-similar iterated root approximants calculated for the transformed series. Here, u stands for the control parameter entering the Borel-Leroy or Mittag-Leffler transformations. Expressions for the marginal index S(u) could be written analytically. They depend only on the control parameters.
Second, after inverse Borel-Leroy transformation, the sought critical index β(u) for the original problem can be expressed as properly scaled by means of a gamma-function, expression for S(u) β(u) = S(u) Γ(u).
Third, the control parameters ought to be determined from the rather intuitive conditions of a non-perturbative nature, such as minimal derivative or minimal difference conditions. The conditions could be imposed on the critical indices β(u) per se.
The examples considered in the paper correspond to the coefficients with widely varying behaviors. They cover situations with a rapid growth of the coefficients a n and with a fast decay, as well as intermediate cases. The methods of optimal summation are able to adapt to such a variety of a n . Both methods of summation appear to be useful in calculation of the critical amplitudes and indices by means of the uniformly and systematically applied technique. The Mittag-Leffler summation seems to work uniformly well for a wider variety of examples, while Borel-Leroy summation works very well for some selected cases, such as the critical index for the Ising model susceptibility.
Performance of the Mittag-Leffler summation method depends on its performance at the two close limits, u = 0 and u = 1. The solution to the optimization problem is often located at positive u, and close to the bounds set by the two limits. However, negative optimal values could also come into play, in particular when the truncated series coefficients a n are decreasing with n. This is why the Mittag-Leffler summation covers more physical cases than Borel summation corresponding only to u = 1. Similar discussion pertains to the Borel-Leroy summation.
Further generalizations of the Borel summation could be imagined according to Nachbin [62], with exponential function in the formula (12) being replaced by some other function belonging to an exponential type, which includes the sums of exponentials or trigonometric functions. Even stretched exponentials could be employed in place of the exponential function, although they do not belong to the exponential type. One can also envisage applying the transforms iteratively or applying various transforms in a sequence [2,63].
This section is not mandatory, but may be added if there are patents resulting from the work reported in this manuscript.