From Asymptotic Series to Self-Similar Approximants

The review presents the development of an approach of constructing approximate solutions to complicated physics problems, starting from asymptotic series, through optimized perturbation theory, to self-similar approximation theory. The close interrelation of underlying ideas of these theories is emphasized. Applications of the developed approach are illustrated by typical examples demonstrating that it combines simplicity with good accuracy.


Introduction
The standard way of treating realistic physical problems, described by complicated equations, relies on approximate solutions of the latter, since the occurrence of exact solutions is rather an exception. The most often used method is a kind of perturbation theory based on expansions in powers of some small parameters. This way encounters two typical obstacles: the absence of small parameters and divergence of resulting perturbative series. To overcome these difficulties, different methods of constructing approximate solutions have been suggested.
In this review, we demonstrate how, starting from asymptotic series, there appear general ideas of improving the series convergence and how these ideas lead to the development of powerful methods of optimized perturbation theory and self-similar approximation theory.

Asymptotic Expansions
Let us be interested in finding a real function f (x) of a real variable x. A generalization to complex-valued functions and variables can be straightforwardly done by considering several real functions and variables. The case of a real function and variable is less cumbersome and allows for the easier explanation of the main ideas. Suppose that the function f (x) is a solution of very complicated equations that cannot be solved exactly and allow only for finding an approximate solution for the asymptotically small variable x → 0 in the form There can happen the following cases.
(i) Expansion over a small variable: where the prefactor f 0 (x) is a given function. The expansion is asymptotic in the sense of Poincaré [1,2], since a n+1 x n+1 a n x n → 0 (x → 0) , with a n assumed to be nonzero. (ii) Expansion over a small function: when the function ϕ(x) tends to zero as x → 0 so that a n+1 ϕ n+1 (x) a n ϕ n (x) → 0 (x → 0) .
Here the value of interest corresponds to the limit ε = 1, while the series is treated as asymptotic with respect to ε → 0, hence a n+1 (x)ε n+1 a n (x)ε n → 0 (ε → 0) .
The introduction of dummy parameters is often used in perturbation theory, for instance in the Euler summation method, Nörlund method, and in the Abel method [9]. Dummy parameters appear when one considers a physical system characterized by a Hamiltonian (or Lagrangian) H, while starting the consideration with an approximate Hamiltonian H 0 , so that one has Then perturbation theory with respect to H − H 0 yields a series in powers of ε. Different iteration procedures also can be treated as expansions in powers of dummy parameters. Sometimes perturbation theory with respect to a dummy parameter is termed nonperturbative, keeping in mind that it is not a perturbation theory with respect to some other physical parameter, say a coupling parameter. Of course this misuse of terminology is confusing, mathematically incorrect, and linguistically awkward. Therefore it is mathematically correct to call perturbation theory with respect to any parameter perturbation theory.

Padé Approximants
The method of Padé approximants sums the series f k (x) = k n=0 a n x n (7) by means of rational fractions with the coefficients b n and c n expressed through a n from the requirement of coincidence of the asymptotic expansions As is evident from their structure, the Padé approximants provide the best approximation for rational functions. However, in general they have several deficiencies. First of all, they are not uniquely defined, in the sense that for a series of order k there are C 2 k + 2 different Padé approximants P M/N , with M + N = k, where and there is no uniquely defined general prescription of which of them to choose. Often, one takes the diagonal approximants P N/N , with 2N = k. However, these are not necessarily the most accurate [11]. Second, there is the annoying problem of the appearance of spurious poles. Third, when the sought function, at small x behaves as in expansion (7), but at large x it may have the power-law behavior x β that should be predicted from the extrapolation of the smallvariable expansion, then this extrapolation to a large variable x ≫ 1 cannot in principle be done if β is not known or irrational. Let us stress that here we keep in mind the extrapolation problem from the knowledge of only the small-variable expansion and the absence of knowledge on the behavior of the sought function at large x. This case should not be confused with the interpolation problem employing the method of two-point Padé approximants, when both expansions at small as well as at large variables are available [12].
Finally, the convergence of the Padé approximants is not a simple problem [11,13], especially when one looks for a summation of a series representing a function that is not known. In the latter case, one tries to observe what is called apparent numerical convergence which may be absent.
As an example of a problem that is not Padé summable [14,15], it is possible to mention the series arising in perturbation theory for the eigenvalues of the Hamiltonian where x ∈ (−∞, ∞), g > 0, and m ≥ 8.

Borel Summation
The series (7) can be Borel summed by representing it as the Laplace integral of the Borel transform B k (t) ≡ k n=0 a n n! t n .
This procedure is regular, since if series (7) converges, then f k (x) = k n=0 a n x n = k n=0 a n n! x n ∞ 0 e −t t n dt = Conditions of Borel summability are given by the Watson theorem [9], according to which a series (7) is Borel summable if it represents a function analytic in a region and in that region the coefficients satisfy the inequality |a n | ≤ C n n! for all orders n. The problem in this method arises because the sought function is usually unknown, hence its analytic properties also are not known, and the behavior of the coefficients a n for large orders n is rarely available. When the initial series is convergent, its Borel transform is also convergent and the integration and the summation in the above formula can be interchanged. However, when the initial series is divergent, the interchange of the integration and summation is not allowed. One has, first, to realize a resummation of the Borel transform and after this to perform the integration.
There are series that cannot be Borel summed. As an example, we can mention a model of a disordered quenched system [16] with the Hamiltonian in which ϕ ∈ (−∞, ∞), ξ ∈ (−∞, ∞), and g > 0, so that the free energy, as a function of the coupling parameter, is where the statistical sum reads as By analytic means and by direct computation of 200 terms in the perturbation expansion for the free energy, it is shown [16] that the series is not Borel summable, since the resulting terms do not converge to any limit. Sometimes the apparent numerical convergence can be achieved by using the Padé approximation for the Borel transform under the Laplace integral, which is termed the Padé-Borel summation.

Optimized Perturbation Theory
The mentioned methods of constructing approximate solutions tell us that there are three main ways that could improve the convergence of the resulting series. These are: (i) the choice of an appropriate initial approximation; (ii) change of variables, and (iii) series transformation. However the pivotal question arises: How to optimize these choices.
The idea of optimizing system performance comes from optimal control theory for dynamical systems [17]. Similarly to dynamical systems, the optimization in perturbation theory implies the introduction of control functions in order to achieve series convergence, as was advanced in Refs. [18][19][20] and employed for describing anharmonic crystals [19][20][21][22][23][24] and the theory of melting [25]. Perturbation theory, complimented by control functions governing the series convergence is called optimized perturbation theory.
The introduction of control functions means the reorganization of a divergent series into a convergent one. Formally, this can be represented as the operation converting an initial series into a new one containing control functions u k (x). Then the optimized approximants are The optimization conditions define the control functions in such a way that to make the new series {F k (x, u k (x))} convergent, because of which this method is named optimized perturbation theory. The general approach to formulating optimization conditions is expounded in the review articles [26,27], and some particular methods are discussed in Refs. [28][29][30]. Control functions can be implanted into perturbation theory in different ways. The main methods are described below.

Initial Approximation
Each perturbation theory or iterative procedure starts with an initial approximation. It is possible to accept as an initial approximation not a fixed form but an expression allowing for variations. For concreteness, assume we are considering a problem characterized by a Hamiltonian H containing a coupling parameter g. Looking for the eigenvalues of the Hamiltonian, using perturbation theory with respect to the coupling, we come to a divergent series As an initial approximating Hamiltonian, we can take a form H 0 (u) containing trial parameters. For brevity, we write here one parameter u. Then we define the Hamiltonian To find the eigenvalues of the Hamiltonian, we can resort to perturbation theory in powers of the dummy parameter ε, yielding Setting ε = 1 and defining control functions u k (x) from optimization conditions results in the optimized approximants The explicit way of defining control functions and particular examples will be described in the following sections.

Change of Variables
Control functions can be implanted through the change of variables. Suppose we consider a series Accomplishing the change of the variable we come to the functions f k (x(z, u)). Expanding the latter in powers of the new variable z, up to the order k, gives In terms of the initial variable, this implies Defining control functions u k (x) yields the optimized approximants (13). When the variable x varies between zero and infinity, sometimes it is convenient to resort to the change of variables mapping the interval [0, ∞) to the interval (−∞, 1], passing to a variable y, where u > 0 and ω > 0 are control parameters [31,32]. The inverse change of variables is The series (18) becomes Expanding this in powers of y, we obtain Defining control functions u k = u k (x) and ω k = ω k (x) gives the optimized approximants Other changes of variables can be found in review [27].

Sequence Transformations
Control functions can also be implanted by transforming the terms of the given series by means of some transformation,T Defining control functions u k (x) gives F k (x, u k (x)). Accomplishing the inverse transformation results in the optimized approximant As an example, let us consider the fractal transform [26] that we shall need in what follows, For this transform, the scaling relation is valid: The scaling power s plays the role of a control parameter through which control functions s k (x) can be introduced [33][34][35].

Statistical Physics
In the problems of statistical physics, before calculating observable quantities, one has to find probabilistic characteristics of the system. This can be either probability distributions, or correlation functions, or Green functions. So, first, one needs to develop a procedure for finding approximations for these characteristics, and then to calculate the related approximations for observable quantities. Here we exemplify this procedure for the case of a system described by means of Green functions [18][19][20][21][22][23]. Let us consider Green functions for a quantum statistical system with particle interactions measured by a coupling parameter g. The single-particle Green function (propagator) satisfies the Dyson equation that can be schematically represented as where G 0 is an approximate propagator and Σ(G) is self-energy [36,37]. Usually, one takes for the initial approximation G 0 the propagator of noninteracting (free) particles, whose self-energy is zero. Then, iterating the Dyson equation, one gets the relation which is a series in powers of the coupling parameter g. Respectively the sequence of the approximate propagators {G k (g)} can be used for calculating observable quantities that are given by a series in powers of g. This is an asymptotic series with respect to the coupling parameter g → 0, which as a rule is divergent for any finite g. Instead, it is possible to take for the initial approximation an approximate propagator G 0 (u) containing a control parameter u. This parameter can, for instance, enter through an external potential [38] corresponding to the self-energy Σ 0 . Then the Dyson equation reads as Iterating this equation [39] yields the approximations for the propagator This iterative procedure is equivalent to the expansion in powers of a dummy parameter. Being dependent on the control parameter u, the propagators G k (g, u) generate the observable quantities A k (g, u) also depending on this parameter. Defining control functions u k (g) results in the optimized approximants for observable quantities.

Optimization Conditions
The above sections explain how to incorporate control parameters into the sequence of approximants that, after defining control functions, become optimized approximants. Now it is necessary to provide a recipe for defining control functions. By their meaning, control functions have to govern the convergence of the sequence of approximants. The Cauchy criterion tells us that a sequence {F k (x, u k )} converges if and only if, for any ε > 0, there exists a number k ε such that for all k > k ε and p > 0.
In optimal control theory [17], control functions are defined as the minimizers of a cost functional. Considering the convergence of a sequence, it is natural to introduce the convergence cost functional [26] in which the Cauchy difference is defined, To minimize the convergence cost functional implies the minimization of the Cauchy difference with respect to control functions, for all k ≥ 0 and p ≥ 0. In order to derive from this condition explicit equations for control functions, we need to accomplish some rearrangements. If the Cauchy difference is small, this means that it is possible to assume that u k+p is close to u k and F k+p is close to F k . Then we can expand the first term of the Cauchy difference in the Taylor series with respect to u k+p in the vicinity of u k , which gives Let us treat F k+p as a function of the discrete variable p, which allows us to expand this function in the discrete Taylor series where a finite difference of m-th order is As examples of finite differences, we can mention Thus the first term in the Cauchy difference can be represented as Keeping in the right-hand side of representation (44) a finite number of terms results in the explicit optimization conditions. The zero order is not sufficient for obtaining optimization conditions, since in this order hence the Cauchy difference is automatically zero, In the first order, we have which gives the Cauchy difference The minimization of the latter with respect to control functions implies Minimizing the first part in the right-hand side of expression (47), we get the minimal-difference condition min for the control functions u k = u k (x). The ultimate form of this condition is the equality The minimization of the second part of the right-hand side of expression (47) leads to the minimal-derivative condition The minimum of condition (50) is made zero by setting When this equation has no solution for the control function u k , it is straightforward to either set or to look for the minimum of the derivative In this way, control functions are defined by one of the above optimization conditions. It is admissible to consider higher orders of expression (44) obtaining higher orders of optimization conditions [27].
Control functions can also be defined if some additional information on the sought function f (x) is available. For instance, when the asymptotic behavior of f (x), as x → x 0 , is known, where then the control functions u k (x) can be defined from the asymptotic condition

Thermodynamic Potential
As an illustration of using the optimized perturbation theory, let us consider the thermodynamic potential of the so-called zero-dimensional anharmonic oscillator model with the statistical sum and the Hamiltonian Taking for the initial approximation the quadratic Hamiltonian in which ω is a control parameter, we define where the perturbation term is Employing perturbation theory with respect to the dummy parameter ε, and setting ε = 1, leads to the sequence of the approximants Control functions for the approximations of odd orders are found from the minimal derivative condition ∂F k (g, ω k ) ∂ω k = 0 (k = 1, 3, . . .) .
Thus we obtain the optimized approximants Their accuracy can be characterized by the maximal percentage error comparing the optimized approximants with the exact expression (56). These maximal errors are As we see, with just a few terms, we get quite good accuracy, while the bare perturbation theory in powers of the coupling parameter g is divergent. Details can be found in review [27]. This simple model allows for explicitly studying the convergence of the sequence of the optimized approximants. It has been proved [40,41] that this sequence converges for both ways of defining control functions, either from the minimal derivative or minimal-difference conditions.

Eigenvalue Problem
Another typical example is the calculation of the eigenvalues of Schrödinger operators, defined by the eigenproblem Hψ n = E n ψ n .
Let us consider a one-dimensional anharmonic oscillator with the Hamiltonian in which x ∈ (−∞, ∞) and g > 0.
For the initial approximation, we take the harmonic oscillator model with a control parameter ω. Following the approach, we define where Employing the Rayleigh-Schrödinger perturbation theory with respect to the dummy parameter ε, we obtain the spectrum E kn (g, ω), where k enumerates the approximation order and n = 0, 1, 2, . . . is the quantum number labeling the states. The zero-order eigenvalue is For odd orders, control functions can be found from the optimization condition For even orders, the above equation does not possess real-valued solutions, because of which we set Using optimized perturbation theory results in the eigenvalues Comparing these with the numerically found eigenvalues E n (g) [42], we define the percentage errors Then we can find the maximal error of the k-th order approximation which gives for the ground state are Again we observe good accuracy and numerical convergence. Recall that the bare perturbation theory in powers of the anharmonicity parameter g diverges for any finite g. The convergence of the sequence of the optimized approximants can be proved analytically [43,44]. More details can be found in Ref. [27].

Nonlinear Schrödinger Equation
The method can be applied to strongly nonlinear systems. Let us illustrate this by considering the eigenvalue problem with the nonlinear Hamiltonian Here N is the number of trapped atoms, the potential is an external potential trapping atoms whose interactions are measured by the parameter where a s is a scattering length. This problem is typical for trapped atoms in Bose-Einstein condensed state [45][46][47][48][49][50][51][52].
The trap anisotropy is characterized by the trap aspect ratio It is convenient to introduce the dimensionless coupling parameter Measuring energy in units of ω ⊥ and lengths in units of l ⊥ , we can pass to dimensionless units and write the nonlinear Hamiltonian as with a dimensionless wave function ψ.
Applying optimized perturbation theory for the nonlinear Hamiltonian [53,54], we take for the initial approximation the oscillator Hamiltonian in which u and v are control parameters. The zero-order spectrum is given by the expression with the radial quantum number n = 0, 1, 2, . . ., azimuthal quantum number m = 0, ±1, ±2, . . ., and the axial quantum number j = 0, 1, 2, . . .. The related wave functions are the Laguerre-Hermite modes. The system Hamiltonian takes the form where the perturbation term is Perturbation theory with respect to the dummy parameter ε gives the energy levels E yielding u k = u k (g) and v k = v k (g). Applications to trapped atoms are discussed in Refs. [53,54].

Hamiltonian Envelopes
When choosing for the initial approximation a Hamiltonian, one confronts the problem of combining two conditions often contradicting each other. From one side, the initial approximation has to possess the properties imitating the studied problem. From the other side, it has to be exactly solvable, providing tools for the explicit calculation of the terms of perturbation theory. If the studied Hamiltonian and the Hamiltonian of the initial approximation are too much different, perturbation theory, even being optimized, may be poorly convergent. In such a case, it is possible to invoke the method of Hamiltonian envelopes [27,55].

General Idea
Suppose we take as an initial approximation a Hamiltonian H 0 that, however, is very different from the considered Hamiltonian H. The difficulty is that the set of exactly solvable problems is very limited, so that sometimes it is impossible to find another Hamiltonian that would be close to the studied form H and at the same time solvable. In that case, we can proceed as follows.
Notice that, if a Hamiltonian H 0 defines the eigenproblem then a function h(H 0 ) satisfies the eigenproblem enjoying the same eigenfunctions. The function h(H) can be called the Hamiltonian envelope [27,55]. Note that, because of the property (92), h(H 0 ) can be any real function.
Accepting h(H 0 ) as an initial Hamiltonian, we obtain the system Hamiltonian with the perturbation term If we find a function h(H 0 ) that better imitates the studied system than the bare H 0 , then the convergence of the sequence of approximations can be improved. The general idea in looking for the function h(H 0 ) is as follows. Let the system Hamiltonian be And let the eigenproblem for a Hamiltonian enjoys exact solutions, although poorly approximating the given system. Looking for the function h(H 0 ), we keep in mind that the most influence on the behavior of wave functions is produced by the region, where the system potential V (r) displays singular behavior tending to ±∞. Suppose this happens at the point r s . Then the function h(H 0 ) has to be chosen such that that is the function h(H 0 ) needs to possess the same type of singularity as the potential of the studied system. Below we illustrate how this choice is done for concrete examples.

Power-Law Potentials
Let us consider the Hamiltonian with a power-law potential in which x ∈ (−∞, ∞), ω 0 > 0, A > 0, and ν > 0. To pass to dimensionless units, we scale the energy and length quantities as The dimensionless coupling parameter is In what follows, in order not to complicate notation, we omit the bars above dimensionless quantities. In dimensionless units, we get the Hamiltonian In order to return to the dimensional form, it is sufficient to make the substitution Taking for H 0 the Hamiltonian we compare the potentials As is evident, the singular point here is x s = ∞. To satisfy condition (97) for ν < 2, we have to take while for ν > 2, we need to accept since now In that way, the Hamiltonian envelope is given by the function

Inverse Power-Law Potentials
The radial Hamiltonian with an inverse power-law potential has the form in which r ≥ 0, l = 0, 1, 2, . . ., A > 0, and ν > 0. Again we can introduce the dimensionless quantities and the dimensionless coupling parameter where ω is arbitrary. Since ω is arbitrary, it can be chosen such that the coupling parameter be unity, In dimensionless units the Hamiltonian becomes This reminds us the Coulomb problem with the Hamiltonian Here u is a control parameter. Comparing the potentials we see that to satisfy condition (97) we have to take the envelope function as as far as Then the Hamiltonian envelope reads as

Logarithmic Potential
As one more example, let us take the radial Hamiltonian of arbitrary dimensionality with the logarithmic potential where r > 0, B > 0, b > 0, and the effective radial quantum number is Again, we need to work with dimensionless quantities, defining and the dimensionless coupling parameter Then, for the simplicity of notation, we omit the bars over the letters and get the dimensionless Hamiltonian Accepting at the starting step the oscillator Hamiltonian we have to compare the potentials Now the singular points are r s = 0 and r s = ∞. This dictates the choice of the envelope function Some explicit calculations can be found in Refs. [27,55]. Optimized perturbation theory, whose main points are expounded above, has been applied to a great variety of problems in statistical physics, condensed matter physics, chemical physics, quantum field theory, etc, as is reviewed in Ref. [27].

Optimized Expansions: Summary
As is explained above, the main idea of optimized perturbation theory is the introduction of control parameters that generate order-dependent control functions controlling the convergence of the sequence of optimized approximants. Control functions can be incorporated in the perturbation theory in three main ways: by choosing an initial approximation containing control parameters, by making a change of variables and resorting to a reexpansion trick, or by accomplishing a transformation of the given perturbation sequence. Control functions are defined by optimization conditions. Of course, there are different variants of implanting control functions and choosing the appropriate variables. In some cases, control functions u k (x) can become control parameters u k , since constants are just a particular example of functions.
Below we summarize the main ideas shedding light on the common points for choosing control functions, the variables for expansions, on the convergence of the sequence of optimized approximants, and on the examples when control functions can be reduced to control parameters. Also, we shall compare several methods of optimization. To make the discussion transparent, we shall illustrate the ideas on the example of a partition function for a zero-dimensional ϕ 4 field theory and on the model of one-dimensional anharmonic oscillator.

Expansion over Dummy Parameters
The standard and often used scheme of optimized perturbation theory is based on the incorporation of control functions through initial approximations, as is mentioned in Sec. 4.1. Suppose we deal with a Hamiltonian H(g) containing a physical parameter g, say coupling parameter. When the problem cannot be solved exactly, one takes a trial Hamiltomian H 0 (u) containing control parameters denoted through u. One introduces the Hamiltonian in which ε is a dummy parameter. One calculates the quantity of interest F k (g, u, ε) by means of perturbation theory in powers of the dummy parameter ε, after which sends this parameter to one, ε → 1.
Employing one of the optimization conditions discussed in Sec. 6, one finds the control functions u k (g). The most often used optimization conditions are the minimal-difference condition and the minimal-derivative condition Substituting the found control functions This scheme of optimized perturbation theory was suggested and employed in Refs. [18][19][20][21][22][23] and in numerous following publications, as can be inferred from the review works [26][27][28][29][30]. As is evident, the same scheme can be used dealing with Lagrangians or action functionals. Instead of the notation ε for the dummy parameter, it is admissible to use any other letter, which, as is clear, is of no importance. Sometimes one denotes the dummy parameter as δ and, using the same standard scheme, one calls it delta expansion. However, using a different notation does not compose a different method.

Scaling Relations: Partition Function
The choice of variables for each particular problem is the matter of convenience. Often it is convenient to use the combinations of parameters naturally occurring in the considered case. These combinations can be found from the scaling relations available for the considered problem.
Let us start with the simple, but instructive, case of the integral representing the partition function (or generating functional) of the so-called zero-dimensional ϕ 4 field theory with the Hamiltonian where g > 0.
Invoking the scaling ϕ −→ λϕ leads to the relation By setting λ = g −1/4 yields the equality And setting These relations show that at large coupling constant the expansion is realized over the combination ω 0 / √ g, while at small coupling constant the natural expansion is over g/ω 2 0 .

Scaling Relations: Anharmonic Oscillator
The other typical example frequently treated for demonstrational purposes is the one-dimensional anharmonic oscillator with the Hamiltonian where g > 0. Let the energy levels E(g, ω 0 ) of the Hamiltonian be of interest. By scaling the spatial variable x −→ λx results in the relation Setting while for λ = ω −1/2 0 we get the relation In particular, for the quartic anharmonic oscillator, with p = 4, we have and Again these relations suggest what are the natural variables for expansions over large or small coupling constants.

Optimized Expansion: Partition Function
The standard scheme of the optimized perturbation theory has been applied to the model (129) many times, accepting as an initial Hamiltonian the form in which ω is a control parameter. Then Hamiltonian (124) becomes Note that Hamiltonian (130) transforms into (141) by means of the replacement Following the standard scheme of optimized perturbation theory for the partition function, and using the optimization conditions for defining control functions, it was found [40,41,44] that at large orders the control functions behave as The minimal-difference and minimal-derivative conditions give α = 1.0729855. It was proved [40,41] that this scheme results in the sequence of optimized approximants for the partition function that converges to the exact numerical value. The convergence occurs for any α > α c = 0.9727803.

Optimized Expansion: Anharmonic Oscillator
The one-dimensional quartic anharmonic oscillator with the Hamiltonian (134), where p = 4 and g > 0, also serves as a typical touchstone for testing approximation methods. The initial approximation is characterized by the harmonic oscillator in which ω is a control parameter. The Hamiltonian (124) takes the form As is seen, the transformation from (134) to (145) is realized by the same substitution (142), with the substitution for ω 0 that can be represented as This shows the appearance of the characteristic combination 1 − (ω 0 /ω) 2 that will be used below. Calculating the energy eigenvalues following the standard scheme, one finds [43,44] the control function with α ≈ 1 for both the minimal-difference and minimal-derivative conditions. The convergence of the sequence of optimized approximants to the exact numerical values [42], found from the solution of the Schrödinger equation, takes place for α > α c = 0.9062077.

Order-Dependent Mapping
Sometimes the procedure can be simplified by transforming the initial expansion, say in powers of a coupling constant, into expansions in powers of other parameters. By choosing the appropriate change of variables, it can be possible to reduce the problem to the form where control functions u k (g) are downgraded to control parameters u k . The change of variables depends on the approximation order, because of which it is called the order-dependent mapping [56].

Change of Variables
Let us be given an expansion in powers of a variable g, f k (g) = k n=0 a n g n .
By analyzing the properties of the considered problem, such as its scaling relations and the typical combinations of parameters arising in the process of deriving perturbative series, it is possible to notice that it is convenient to denote some parameter combinations as new variables. Then one introduces the change of variables where is treated as a control parameter. By substituting (149) into (148) gives the function f k (g(z, u)), which has to be expanded in powers of z up to order k, leading to the series The minimal-difference condition defining the control parameters u k . Since, according to (150), the value u k denotes the combination of parameters u k = u k (g, z), hence it determines the control functions z k (g). The pair u k and z k (g), being substituted into (151), results in the optimized approximants Thus, the convenience of the chosen change of variables is in the possibility of dealing at the intermediate step with control parameters instead of control functions that appear at a later stage.

Partition Function
To illustrate the method, let us consider the partition function (129) following the described scheme [56]. From the substitution (142) it is clear that natural combinations of parameters appearing in perturbation theory with respect to the term with ε in the Hamiltonian (141) are and Then the combination of parameters (150) reads as In order to simplify the notation, it is possible to notice that the parameter ω always enters the equations being divided by ω 0 . Therefore, measuring ω in units of ω 0 is equivalent to setting ω 0 → 1. In these units, .
Finding from the minimal-difference condition (152) the control parameter u k and using definition (157) gives the control function Then relation (155) results in the control function Finally, one gets the partition function Z k (z k (g), u k ). This procedure, with the change of variables used above, has been shown [57] to be equivalent to the standard scheme of optimized perturbation theory resulting in optimized approximants Z k (g).

Anharmonic Oscillator
Again using the dimensionless units, as in the previous section, one sets the notations Then the combination (150) becomes Similarly to the previous section, one finds the control parameter u k and from (161) one obtains the control functions z k (g) and ω k (g). The resulting energy levels E k (z k (g), u k ) coincide with the optimized approximants E k (g), as has been proved in [57].

Variational Expansions
The given expansion over the coupling constant (148) can be reexpanded with respect to other variables in several ways. One of the possible reexpansions has been termed variational perturbation theory [31]. Below it is illustrated by the example of the anharmonic oscillator in order to compare this type of a reexpansion with other methods. Let us consider the energy levels of the anharmonic oscillator with the Hamiltonian (134) with p = 4. As is clear from the scaling relations of Sec. 11, the energy can be represented as an expansion We have the identity that is a particular case of the substitution (142) with the control parameter ω and ε = 1. Employing the notation where it is straightforward to rewrite the identity (163) in the form This form is substituted into expansion (162), which then is reexpanded in powers of the new variable g/ω 3 , while keeping u untouched and setting ω 0 to one. The reexpanded series is truncated at order k. Comparing this step with the expansion in Sec. 11, it is evident that this is equivalent to the expansion over the dummy parameter ε. And comparing the expansion over g/ω 3 with the expansion over z in Sec. 12, we see that they are also equivalent. Thus we come to the expansion where Then one substitutes back the expression (165) for u = u(g, ω).
The control function ω k (g) is defined by the minimal derivative condition, or, when the latter does not have real solutions, by the zero second derivative over ω of the energy E k (g, ω). The found control function ω k (g) is substituted into E k (g, ω), thus giving the optimized approximant The equivalence of the above expansion in powers of g/ω 3 to the expansions with respect to the dummy parameter ε, or with respect to the parameter z, becomes evident if we use the notation of the present section and notice that the substitution (146) can be written as This makes it immediately clear that the expansion over g/ω 3 , with keeping u untouched, is identical to the expansion over the dummy parameter ε.

Control Functions and Control Parameters
It is important to remark that it is necessary to be cautious introducing control functions through the change of variables and reexpansion. Strictly speaking, such a change cannot be postulated arbitrarily. When the change of variables is analogous to the procedure of using the substitutions, such as (142), (146) or (169), naturally arising in perturbation theory, as in Sec. 4, then the results of these variants will be close to each other. However, if the change of variables is arbitrary, the results can be not merely inaccurate, but even qualitatively incorrect [27,58]. It is also useful to mention that employing the term control functions, we keep in mind that in particular cases they can happen to become parameters, although order-dependent. Then instead of functions u k (x) we can have parameters u k . There is nothing wrong in this, as far as parameters are a particular example of functions. The reduction of control functions to control parameters can occur in the following cases.
It may happen that in the considered problem there exists such a combination of characteristics that compose the quantities u k depending only on the approximation order but not depending on the variable x. For instance, this happens in the mapping of Sec. 12, where the combinations u k = u k (g, ω k (g) play the role of control parameters. In the case of the partition function, this is the combination (157) and for the anharmonic oscillator, it is the combination (161).
The other example is the existence in the applied optimization of several conditions restricting the choice of control parameters. The typical situation is when the optimization condition consists in the comparison of asymptotic expansions of the sought function and of the approximant. Suppose that, in addition to the small-variable expansion we know the large-variable expansion of the sought function Let us assume that we have found the optimized approximant F k (x, u k (x)), where the control functions u k (x) are defined by one of the optimization conditions of Sec. 6. These conditions provide a uniform approximation of the sought function on the whole interval of its definition. However the resulting approximants F k (x, u k (x)) are not required to give exact coefficients of asymptotic expansions either at small or at large variable x. If we wish that these asymptotic coefficients would exactly coincide with the coefficients of the known asymptotic expansions (170) and (171), then we have to implant additional control parameters and impose additional asymptotic conditions. This can be done by using the method of corrected Padé approximants [27,[59][60][61][62]. To this end, we define the optimized approximant as where is a diagonal Padé approximant, whose coefficients c n and d n , playing the role of control parameters, are prescribed by the accuracy-through-order procedure, so that the asymptotic expansions of (172) would coincide with the given asymptotic expansions of the sought function at small x, and at large x, The number of the parameters in the Padé approximant is such that to satisfy the imposed asymptotic conditions (174) and (175).

Self-Similar Approximation Theory
As has been emphasized above, the idea of introducing control functions for the purpose of governing the convergence of a sequence stems from the optimal control theory, where one introduces control functions in order to regulate the trajectory of a dynamical system, for instance so that to force the trajectory to converge to a desired point. The analogy between perturbation theory and the theory of dynamical systems has been strengthened even more in the self-similar approximation theory [26,27,[63][64][65][66][67]. The idea of this theory is to consider the transfer from one approximation to another as the motion on the manifold of approximants, where the approximation order plays the role of discrete time. Suppose, after implanting control functions, as explained in Sec.4, we have the sequence of approximants F k (x, u k ). Recall that the control functions can be defined in different ways, as has been discussed above. Therefore we, actually, have the manifold of approximants associated with different control functions, This will be called the approximation manifold. Generally, it could be possible to define a space of approximants. However the term approximation space is used in mathematics in a different sense [68]. So, we shall deal with the approximation manifold. The transfer from an approximant F k to another approximant F k+p can be understood as the motion with respect to the discrete time, whose role is played by the approximation order k. The sequence of approximants F k (x, u k ) with a fixed choice of control functions u k = u k (x) defines a trajectory on the approximation manifold (176). Let us fix the rheonomic constraint defining the expansion function x k (f ). Recall that in the theory of dynamical systems a rheonomic constraint is that whose constraint equations explicitly contain or are dependent upon time. In our case, time is the approximation order k. The inverse constraint equation is Let us introduce the endomorphism by the definition acting as This endomorphism and the approximants are connected by the equality The set of endomorphisms forms a dynamical system in discrete time with the initial condition y 0 (f ) = f .
By this construction, the sequence of endomorphisms {y k (f )}, forming the dynamical system trajectory, is bijective to the sequence of approximants {F k (x, u k (x))}. Since control functions, by default, make the sequence of approximants F k (x, u k (x)) convergent, this means that there exists a limit And as far as the sequence of approximants is bijective to the trajectory of the dynamical system, there should exist the limit y * (f ) = lim k→∞ y k (f ) .
This limit, being the final point of the trajectory, implies that it is a fixed point, for which Thus to find the limit of an approximation sequence is equivalent to determining the fixed point of the dynamical system trajectory. We may notice that for large p, the self-similar relation holds: which follows from conditions (185) and (186). As far as in the real situations it is usually impossible to reach the limit of infinite approximation order, we assume the validity of the selfsimilar relation for finite approximation orders: This relation implies the semi-group property y k · y p = y k+p , y 0 = 1 .
The dynamical system in discrete time (182) with the above semi-group property is called cascade (semicascade). The theory of such dynamical systems is well developed [69,70]. In our case, this is an approximation cascade [27]. Since, as is said above, in realistic situations we are able to deal only with finite approximation orders, we can find not an exact fixed point y * (f ), but an approximate fixed point y * k (f ). The corresponding approximate limit of the considered sequence is If the form F k (x, u k ) is obtained by means of a transformation like in (27), then the resulting self-similar approximant reads as

Embedding Cascade into Flow
Usually, it is more convenient to deal with dynamical systems in continuous time than with systems in discrete time. For this purpose, it is possible to embed the approximation cascade into an approximation flow, which is denoted as and implies that the endomorphism in continuous time enjoys the same group property as the endomorphism in discrete time, that the flow trajectory passes through all points of the cascade trajectory, and starts from the same initial point, The self-similar relation (194) can be represented as the Lie equation in which v(y(t, f )) is a velocity field. Integrating the latter equation yields the evolution integral where t k is the time required for reaching the fixed point y * k = y * k (f ) from the approximant y k = y k (f ). Using relations (181) and (190), this can be rewritten as where F k = F k (x, u k (x)) and F * k = F * k (x, u k (x)). The velocity field can be represented resorting to the Euler discretization This is equivalent to the form in which We may notice that the velocity field is directly connected with the Cauchy difference (39), As is explained in Sec. 6, the Cauchy difference of zero order equals zero, hence in that order the velocity is zero, and F * k = F k . The Cauchy difference of first order is nontrivial, being given by expression (46). In this order, the velocity field becomes The smaller the velocity, the faster the fixed point is reached. Therefore control functions should be defined so that to make the velocity field minimal: Thus we return to the optimization conditions of optimized perturbation theory, discussed in Sec. 6. Opting for the optimization condition

Stability Conditions
The sequence {y k (f )} defines the trajectory of the approximation cascade that is a type of a dynamical system. The motion of dynamical systems can be stable or unstable. The stability of motion for the approximation cascade can be characterized [27,67,71] similarly to the stability of other dynamical systems [69,72,73]. Dealing with real problems, one usually considers finite steps k. Therefore the motion stability can be defined only locally. The local stability at the k-th step is described by the local map multiplier The motion at the step k, starting from an initial point f , is stable when The maximal map multiplier defines the global stability with respect to f , provided that The maximum is taken over all admissible values of f . The image of the map multiplier (207) on the manifold of the variable x is The motion at the k-th step at the point x is stable if Respectively, the motion is globally stable with respect to the domain of x when the maximal map multiplier The map multiplier at the fixed point y * k (f ) is The fixed point is locally stable when and it is globally stable with respect to f if the maximal multiplier satisfies the inequality µ * k < 1 .
The above conditions of stability can be rewritten in terms of the local Lyapunov exponents The motion at the k-th step is stable provided the Laypunov exponents are negative. The occurrence of local stability implies that the calculational procedure should be numerically convergent at the considered steps. Thus, even not knowing the exact solution of the problem and being unable to reach the limit of k → ∞, we can be sure that the local numerical convergence for finite k is present.

Free Energy
In order to demonstrate that the self-similar approximation theory improves the results of optimized perturbation theory, it is instructive to consider the same problem of calculating the free energy (thermodynamic potential) of the model discussed in Sec. 7, with the statistical sum (57). Following Sec. 7, we accept the initial Hamiltonian (59) and define Hamiltonian (60). Expanding the free energy (220) in powers of the dummy parameter ε, we have the sequence of approximants (62). The control functions ω k (g) are defined by the optimization conditions (63) and (64), which give ω k (g) = 1 2 1 + 1 + 12s k g where The rheonomic constraint (177) takes the form From here, we find the expansion function The endomorphism (180) reads as with the coefficients A kn given in Refs. [27,71,74,75], and where The cascade velocity (200) becomes Taking the evolution integral (199), with t k = 1, we come to the self-similar approximants The accuracy of the approximations is described by the percentage errors where f (g) is the exact numerical value of expression (220). Here we have The map multipliers (207) are The coupling parameter g pertains to the domain [0, ∞), Then f ∈ [0, ∞), and α(f ) ∈ [0.1). The maximal map multiplier (209) is found to satisfy the stability condition (210).

Fractal Transform
As is explained in Sec. 4, control functions can be incorporated into a perturbative sequence either through initial conditions, or by means of the change of variables, or by a sequence transformation.
In the above example of Sec. 14, we have considered the implantation of control functions into initial conditions. Now we shall study another way, when control functions are incorporated through a sequence transformation. Let us consider an asymptotic series in which f 0 (x) is a given function. Actually, it is sufficient to deal with the series To return to the case of series (230), we just need to make the substitution Following the spirit of self-similarity, we can remember that the latter is usually connected with the power-law scaling and fractal structures [76][77][78]. Therefore, it looks natural to introduce control functions through a fractal transform [79], say of the type [26,27,[33][34][35] The inverse transformation is With the series (231), we have a n x n+s .
As is mentioned in Sec. 4, the scaling relation (30) is valid. The scaling exponent s plays the role of a control parameter.
In line with the self-similar approximation theory, we define the rheonomic constraint yielding the expansion function The dynamic endomorphism becomes And the cascade velocity is What now remains is to consider the evolution integral.

Self-Similar Root Approximants
The differential equation (197) can be rewritten in the integral form Substituting here the cascade velocity (239) gives the relation where Accomplishing the inverse transformation (234) leads to the equation The explicit form of the latter is the recurrent relation Using the notation and iterating this relation k − 1 times results in the self-similar root approximant This approximant is convenient for the problem of interpolation, where one can meet different situations.
(i) The k coefficients a n of the asymptotic expansion (231) up to the k-th order are known and the exponent β of the large-variable behavior of the sought function is available, where although the amplitude B is not known. Then, setting the control functions s j = s, from Eq. (244), we have and the root approximant (245) becomes For large variables x, the latter behaves as with the amplitude and exponent β k = km k .
Equating β k to the known exponent β, we find the root exponent All parameters A n can be found from the comparison of the initial series (231) with the smallvariable expansion of the root approximant (248), which is called the accuracy-trough-order procedure. Knowing all A n , we obtain the large-variable amplitude B k .
(ii) The k coefficients a n of the asymptotic expansion (231) up to the k-th order are available and the amplitude B of the large-variable behavior of the sought function is known, but the large-variable exponent β is not known. Then the parameters A n again are defined through the accuracy-through-order procedure (253). Equating the amplitudes B k and B results in the exponent (iii) The k coefficients a n of the asymptotic expansion (231) are known and the large-variable behavior (246) is available, with both the amplitude B and exponent β known. Then, as earlier, the parameters A n are defined from the accuracy-through-order procedure and the exponent m k is given by Eq. (252). The amplitude B k can be found in two ways, from expression (250) and equating B k and B. The difference between the resulting values defines the accuracy of the approximant.
Then considering the root approximant (245) for large x → ∞, and comparing this expansion with the asymptotic form (255) we find all parameters A n expressed through the coefficients b n , and the large-variable internal exponents are while the external exponent is It is important to mention that the external exponent m k can be defined even without knowing the large-variable behavior of the sought function. This can be done by treating m k as a control function defined by an optimization condition from Sec. 6. This method has been suggested in Ref. [33].
Notice that when it is more convenient to deal with the series for large variables, it is always possible to use the same methods as described above by transferring the large-variable expansions into small-variable ones by means of the change of the variable z = 1/x.

Self-Similar Nested Approximants
It is possible to notice that the series can be represented as the sequence etc., through up to the last term Applying the self-similar renormalization at each order of the sequence, considering ϕ j as variables, we obtain the renormalized sequence in which b j = t j n j , n j = −s j (j = 1, 2 . . . , k − 1) .

Using the notation
we come to the self-similar nested approximant For large x, this gives with the amplitude B k = A n 1 1 A n 1 n 2 2 A n 1 n 2 n 3 3 . . . A n 1 n 2 n 3 ...n k k (267) and the exponent β k = n 1 + n 1 n 2 + n 1 n 2 n 3 + . . . + n 1 n 2 n 3 . . . n k .
If we change the notation for the external exponent to and keep the internal exponents constant, then the large-variable exponent becomes When the exponent β of the large-variable behavior is known, where then, setting β k = β, gives The parameter m should be defined so that to provide numerical convergence for the sequence {f * k (x)}. For instance, if m = 1, then using the asymptotic form we get In the latter case, the nested approximant (265), with the notation The same form can be obtained by setting in the root approximant (245) all internal exponents n j = 1.
The external exponent m k can also be defined by resorting to the optimization conditions of Sec. 6. Several applications of the nested approximants are given in [85].

Self-Similar Exponential Approximants
When it is expected that the behavior of the sought function is rather exponential, but not of power law, then in the nested approximants of the previous section, we can sent n j → ∞, hence b j → 0 and A j → 0. This results in the self-similar exponential approximants [86] in which C n = a n a n−1 t n (n = 1, 2, . . . , k) .
The parameters t n are to be defined from additional conditions [26,27], so that the sequence of the approximants be convergent. It is often sufficient to set t n = 1/n. This expression appears as follows. By its meaning, t n is the effective time required for reaching a fixed point from the previous step. Accomplishing n steps takes time of order nt n . The minimal time corresponds to one step. Equating nt n and one gives t n = 1/n. Some other ways of defining the control parameters t n are considered in Refs. [26,27,86].

Self-Similar Factor Approximants
By the fundamental theorem of algebra [87], a polynomial of any degree of one real variable over the field of real numbers can be split in a unique way into a product of irreducible first-degree polynomials over the field of complex numbers. This means that series (259) can be represented in the form with the coefficients b j expressed through a n . Applying the self-similar renormalization procedure to each of the factors in turn results in the self-similar factor approximants [88][89][90] where N k = k/2 , k = 2, 4, 6, . . .
The control parameters A j and n j are defined by the accuracy-through-order procedure by equating the like order terms in the expansions f * k (x) and f k (x), In the present case, it is more convenient to compare the corresponding logarithms This leads to the system of equations N k j=1 n j A n j = D n (n = 1, 2, . . . , k) , in which This system of equations enjoys a unique (up to enumeration permutation) solution for all A j and n j when k is even, and when k is odd, one of A j can be set to one [27,91].
At large values of the variable, we have where the amplitude and the large-variable exponent are If the large-variable exponent is known, for instance from scaling arguments, so that then equating β k and β imposes on the exponents of the factor approximant the constraint The self-similar factor approximants have been used for a variety of problems, as can be inferred from Refs. [27,[88][89][90][91][92]].

Self-Similar Combined Approximants
It is possible to combine different types of self-similar approximants as well as these approximants and other kinds of approximations.

Different Types of Approximants
Suppose we are given a small-variable asymptotic expansion which we plan to convert into a self-similar approximation. At the same time, we suspect that the behavior of the sought function is quite different at small and at large variables. In such a case, we can combine different types of self-similar approximants in the following way. We take in series (287) several initial terms, a j x j (n < k) (288) and construct of them a self-similar approximant f * n (x). Then we define the ratio and expand the latter in powers of x as Constructing a self-similar approximant C * k/n (x), we obtain the combined approximant The approximants f * n (x) and C * k/n (x) can be represented by different forms of self-similar approximants. For example, it is possible to define f * n (x) as a root approximant, while C * k/n (x) as a factor or exponential approximant, depending on the expected behavior of the sought function [93].

Self-Similar Padé Approximants
Instead of two different self-similar approximants, it is possible, after constructing a self-similar approximant f * n (x), to transform the remaining part (290) into a Padé approximant P M/N (x), with M + N = k − n, so that The result is the self-similarly corrected Padé approximant, or briefly, the self-similar Padé approximant [59][60][61] The advantage of this type of approximants is that they can correctly take into account irrational behavior of the sought function, described by the self-similar approximant f * n (x), as well as the rational behavior represented by the Padé approximant P M/N (x).
Note that Padé approximants (8) actually are a particular case of the factor approximants (278), where M factors correspond to n j = 1 and N factors, to n j = −1. This is because the Padé approximants can be represented as (1 + C n x) −1 .

Self-Similar Borel Summation
It is possible to combine self-similar approximants with the method of Borel summation. According to this method, for a series (287), one can define [9,94] the Borel-Leroy transform where u is chosen so that to improve convergence. The series (294) can be summed using one of the self-similar approximations, and converting u into a control parameter u k , thus getting B * k (t, u k ). Then the self-similar Borel-Leroy summation yields the approximant The case of the standard Borel summation corresponds to u k = 0. Then the self-similar Borel summation gives In addition to the considered above combinations of different summation methods, one can use other combinations. For example, the combination of exponential approximants and continued fractions has been employed [95].

Self-Similar Data Extrapolation
One often meets the following problem. There exists an ordered dataset {f n : n = 1, 2, . . . , k} labeled by the index n, and one is interested in the possibility of predicting the values f k+p outside this dataset. The theory of self-similar approximants suggests a solution to this problem [27,96]. Let us consider several last datapoints, for instance the last three points How many datapoints one needs to take depends on the particular problem considered. For the explicit illustration of the idea, we take three datapoints. The chosen points can be connected by a polynomial spline, in the present case, by a quadratic spline defined so that From this definition it follows Treating polynomial (299) as an expansion in powers of t makes it straightforward to employ self-similar renormalization, thus, obtaining a self-similar approximant g * (t). For example, resorting to factor approximants, we get with the parameters The approximants g * (t), with t ≥ 2 provide the extrapolation of the initial dataset. The nearest to the dataset extrapolation point can be estimated as This method can also be used for improving the convergence of the sequence of self-similar approximants. Then the role of datapoints f k is played by the self-similar approximants f * k (x). In that case, all parameters a = a(x), b = b(x), c = c(x), as well as A = A(x) and m = m(x) become control functions. This method of data extrapolation has been used for several problems, such as predictions for time series and convergence acceleration [27,96,97].

Self-Similar Diff-Log Approximants
There is a well known method employed in statistical physics called diff-log transformation [98,99]. This transformation for a function f (x) is The inverse transformation, assuming that the function f (x) is normalized so that reads as When we start with an asymptotic expansion the diff-log transformation gives Expanding the latter in powers of x yields with the coefficients b n expressed through a n . This expansion can be summed by one of the selfsimilar methods giving D * k (x). Involving the inverse transformation (305) results in the self-similar diff-log approximants A number of applications of the diff-log transformation can be found in Refs. [61,83,99], where it is shown that the combination of the diff-log transform with self-similar approximants gives essentially more accurate results than the diff-log Padé method.

Critical Behavior
One says that a function f (x) experiences critical behavior at a critical point x c , when this function at that point either tends to zero or to infinity. It is possible to distinguish two cases, when the critical behavior occurs at infinity, and when at a finite critical point. These two cases will be considered below separately.

Critical Point at Infinity
If the critical behavior happens at infinity, the considered function behaves as Then the diff-log transform tends to the form Here B is a critical amplitude, while β is a critical exponent. The critical exponents have a special interest for critical phenomena. If we are able to define a self-similar approximation f * k (x) directly to the studied function f (x), then the critical exponent can be found from the limit Otherwise, it can be obtained from the equivalent form where a self-similar approximation for the diff-log transform D * k (x) is needed. The convenience of using the representation (313) is in the possibility of employing a larger arsenal of different self-similar approximants. Of course, the factor approximants can be involved in both the cases. However, the root and nested approximants require the knowledge of the large-variable exponent of the sought function, which is not always available. On the contrary, the large-variable behavior of the diff-log transform (311) is known. Therefore for constructing a self-similar approximation for the diff-log transform, we can resort to any type of self-similar approximants.
It is necessary to mention that the root and nested approximants can be defined, without knowing the large-variable behavior, by invoking optimization conditions of Sec. 6 prescribing the value of the external exponent m k , as is explained in Ref. [33]. However this method becomes rather cumbersome for high-order approximants.

Finite Critical Point
If the critical point is located at a finite x c that is in the interval (0, ∞), then Here the diff-log transform behaves as Again, the critical exponent can be derived from the limit provided a self-similar approximant f * k (x) is constructed. But it may happen that the other form is more convenient, where a self-similar approximant for the diff-log transform D * k (x) is easier to find. This is because the nearest to zero pole of D * k (x) defines a critical point x c , while the residue (317) yields a critical exponent.
Note that by the change of the variable the problem of a finite critical point can be reduced to the case of critical behavior at infinity. For instance one can use the change of the variable z = x/(x c −x) or any other change of the variable mapping the interval [0, x c ) to [0, ∞). Numerous examples of applying the diff-log transform, accompanied by the use of self-similar approximants, are presented in Refs. [61,83,99], where it is also shown that this method essentially outperforms the diff-log Padé variant.

Non-Power-Law Behavior
In the previous sections we were mainly keeping in mind a kind of power-law behavior of considered functions at large variables. Now it is useful to make some comments on the use of the described approximation methods for other types of behavior. The most often met types of behavior that can occur at large variables are the exponential and logarithmic behavior. Below we show that the developed methods of self-similar approximants can be straightforwardly applied to any type of behavior.

Exponential Behavior
The exponential behavior with respect to time happens in many mathematical models employed for describing the growth of population, mass of biosystems, economic expansion, financial markets, various relaxation phenomena, etc. [100][101][102][103][104][105].
When a sought function at a large variable displays exponential behavior, there are several ways of treating this case. First of all, this kind of behavior can be treated by self-similar exponential approximants of Sec. 18. The other way is to resort to diff-log approximants of Sec. 22 or, simply, to consider the logarithmic transform If the sought function at large variable behaves as then Therefore the function behaves as Keeping in mind the asymptotic series (306), we have which can be expanded in powers of x giving

Logarithmic Behavior
When there is suspicion that the sought function exhibits logarithmic behavior at large variables, it is reasonable to act by analogy with the previous subsection, but now defining the exponential transform For the asymptotic series (306), we have whose expansion in powers of x produces This can be converted into a self-similar approximation E * k (x), so that the final answer becomes As an example, let us consider the function with the logarithmic behavior at large variables, This function has the expansion Defining factor approximants E * k (x), we obtain the approximants (331), whose large-variable behavior is of correct logarithmic form with the amplitudes B k , Calculation of critical exponents is one of the most important problems in the theory of phase transitions. Here we show how the critical exponents can be calculated by using self-similar factor approximants applied to the asymptotic series in powers of the ε = 4 − d, where d is space dimensionality. We shall consider the O(N) ϕ 4 field theory in d = 3. The definition of the critical exponents can be found in reviews [27,117]. One usually derives the so-called epsilon expansions for the exponents η, ν −1 , and ω. The other exponents can be obtained from the scaling relations In three dimensions, one has The number of components N corresponds to different physical systems. Thus N = 0 corresponds to dilute polymer solutions, N = 1, to the Ising universality class, N = 2, to superfluids and the so-called XY magnetic models, N = 3, to the Heisenberg universality class, and N = 4, to some models of quantum field theory. Formally, it is admissible to study arbitrary N. In the case of N = −2, the critical exponents for any d are known exactly: For the limit N → ∞, the exact exponents also are available: The latter for d = 3 reduce to The epsilon expansion results in the series obtained for ε → 0, while at the end we have to set ε = 1. Direct substitution of ε = 1 in the series (348) leads to the values having little to do with real exponents. These series require to define their effective sums, which we accomplish by means of the self-similar factor approximants f * k (ε) = f 0 (ε) Then we set ε = 1 and define the final answer as the half sum of the last two factor approximants f k (1) and f k−1 (1). Let us first illustrate the procedure for the O(1) field theory of the Ising universality class, where there exist the most accurate numerical calculations of the exponents, obtained by Monte Carlo simulations [117][118][119][120][121]. The epsilon expansions for η, ν −1 , and ω can be written [122] as η ≃ 0.0185185ε 2 + 0.01869ε 3 − 0.00832877ε 4 + 0.0256565ε 5 , ν −1 ≃ 2 − 0.333333ε − 0.117284ε 2 + 0.124527ε 3 − 0.30685ε 4 − 0.95124ε 5 , ω ≃ ε − 0.62963ε 2 + 1.61822ε 3 − 5.23514ε 4 + 20.7498ε 5 . (350) If we set here ε = 1, we get senseless values η = 0.0545, ν = 2.4049 and ω = 17.5033. However by means of the self-similar factor approximants we obtain the results shown in Table 3, which are in good agreement with Monte Carlo simulations [117][118][119][120][121].
The use of the self-similar factor approximants can be extended to the calculation of the critical exponents for the arbitrary number of components N of the O(N) symmetric ϕ 4 field theory in d = 3. In the general case, the epsilon expansions [122] read as Summing these series by means of the self-similar factor approximants [123,124], we obtain the exponents presented in Table 4. The found values of the exponents are in good agreement with experimental data as well as with the results of numerical methods, such as Padé-Borel summation and Monte Carlo simulations. It is important to stress that when the exact values of the exponents are known (for N = −2 and N → ∞), the self-similar approximants automatically reproduce these exact data.

Conclusion
In this review, we have presented the basic ideas of the approach allowing for obtaining senseful results from divergent asymptotic series typical of asymptotic perturbation theory. The pivotal points of the approach can be emphasized as follows.
(i) The implantation of control functions in the calculational procedure, treating perturbation theory as optimal control theory. Control functions are defined by optimization conditions so that to control the convergence of the sequence of optimized approximants. The optimization conditions are derived from the Cauchy criterion of sequence convergence. The resulting optimized perturbation theory provides good accuracy even for very short series of just a few terms and Table 1: Coefficients a n of the loop expansion c 1 (x) for the number of components N.   [116] makes it possible to extrapolate the validity of perturbation theory to arbitrary values of variables, including the limit to infinity.
(ii) Reformulation of perturbation theory to the language of dynamical theory, handling the motion from one approximation term to another as the motion in discrete time played by the approximation order. Then the approximation sequence is bijective to the trajectory of the effective dynamical system, and the sequence limit is equivalent to the trajectory fixed point. The motion near the fixed point enjoys the property of functional self-similarity. The approximation dynamical system in discrete time is called cascade. The approximation cascade can be embedded into a dynamical system in continuous time termed approximation flow. The representation in the language of dynamical theory allows us to improve the accuracy of optimized perturbation theory, to study the procedure stability, and to select the best initial approximation.
(iii) Introduction of control functions by means of a fractal transformation of asymptotic series, which results in the derivation of several types of self-similar approximants. These approximants combine the simplicity of their use with good accuracy. They can be employed for the problem of interpolation as well as extrapolation.
The application of the described methods is illustrated by several examples demonstrating the efficiency of the approach.  Table 4: Critical exponents for the O(N)-symmetric ϕ 4 field theory obtained by the summation of ε expansions using self-similar factor approximants.