LETTER Explaining a changeover from normal to super diffusion in time-dependent billiards To cite this article: Matheus Hansen et al 2018 EPL 121 60003   View the article online for updates and enhancements. Related content Fermi acceleration in time-dependent billiards: theory of the velocity diffusion in conformally breathing fully chaotic billiards Benjamin Batisti and Marko Robnik - Phase space interpretation of exponential Fermi acceleration Benno Liebchen, Robert Büchner, Christoph Petri et al. - Separation of particles leading either to decay or unlimited growth of energy in a driven stadium-like billiard André L P Livorati, Matheus S Palmero, Carl P Dettmann et al. - This content was downloaded from IP address 186.217.236.60 on 22/05/2019 at 20:07 https://doi.org/10.1209/0295-5075/121/60003 http://iopscience.iop.org/article/10.1088/1751-8113/44/36/365101 http://iopscience.iop.org/article/10.1088/1751-8113/44/36/365101 http://iopscience.iop.org/article/10.1088/1751-8113/44/36/365101 http://iopscience.iop.org/article/10.1088/1367-2630/13/9/093039 http://iopscience.iop.org/article/10.1088/1367-2630/13/9/093039 http://iopscience.iop.org/article/10.1088/1751-8113/47/36/365101 http://iopscience.iop.org/article/10.1088/1751-8113/47/36/365101 http://iopscience.iop.org/article/10.1088/1751-8113/47/36/365101 March 2018 EPL, 121 (2018) 60003 www.epljournal.org doi: 10.1209/0295-5075/121/60003 Explaining a changeover from normal to super diffusion in time-dependent billiards Matheus Hansen 1 , David Ciro 2 , Iberê L. Caldas 1 and Edson D. Leonel 3 1 Instituto de F́ısica, Universidade de São Paulo - São Paulo - CEP 05314-970 SP, Brazil 2 Instituto de F́ısica, Universidade Federal do Paraná - Curitiba - CEP 81531-990 PR, Brazil 3 Departamento de F́ısica, UNESP - Rio Claro - CEP 13506-900 SP, Brazil received 29 January 2018; accepted in final form 25 April 2018 published online 17 May 2018 PACS 05.45.-a – Nonlinear dynamics and chaos PACS 05.45.Pq – Numerical simulations of chaotic systems PACS 05.40.Fb – Random walks and Levy flights Abstract – The changeover from normal to super diffusion in time-dependent billiards is explained analytically. The unlimited energy growth for an ensemble of bouncing particles in time-dependent billiards is obtained by means of a two-dimensional mapping of the first and second moments of the speed distribution function. We prove that, for low initial speeds the average speed of the ensemble grows with exponent ∼1/2 of the number of collisions with the boundary, therefore exhibiting normal diffusion. Eventually, this regime changes to a faster growth characterized by an exponent ∼1 corresponding to super diffusion. For larger initial energies, the temporary symmetry in the diffusion of speeds explains an initial plateau of the average speed. Copyright c© EPLA, 2018 As coined by Enrico Fermi [1] Fermi acceleration (FA) is a phenomenon where an ensemble of classical and non- interacting particles acquires energy from repeated elastic collisions with a rigid and time-varying boundary. It is typically observed in billiards [2–4] whose boundaries are moving in time [5–9]. If the motion of the boundary is random and the initial particle speeds are small (compared to the moving boundary), the average speed grows at a rate proportional to n1/2, with n denoting the number of collisions. If the initial speeds are large, a transient plateau in their average is observed in a plot V vs. n. This is due to the symmetric diffusion of the speeds to values above and below the initial value [10]. This simmetry is, however, broken when the distribution of speeds reaches a lower bound forcing the mean speed to exhibit a sustained growth. For deterministic oscillations of the border, the scenario is different. Breathing oscillations preserve the shape but not the area of the billiard. It is known that the average speed evolves in a sub-diffusive manner with exponent of the order 1/6 [11,12]. On the other hand, for oscillations preserving the area but not the shape of the billiard there are two regimes of growth. For short times the diffusion of speeds is normal, with growth exponent ∼1/2, then, after a sufficiently large number of collisions, it enters to a super diffusion regime, with exponent ∼1 [13]. This changeover is, so far, not yet explained and our contribution in this letter is to fill up this gap in the theory. This is achieved by studying the momenta of the speed distribution function, noticing that the dynamical angular/time variables have an inhomogeneous distribution in phase space. The results presented in this letter are illustrated by a time-dependent oval-billiard [3] whose phase space is mixed when the boundary is static. The boundary of the billiard is written as Rb(θ, t) = 1 + ε [1 + a cos(t)] cos(pθ), where Rb is the radius of the boundary in polar coordi- nates, θ is the polar angle, ε controls the circle deforma- tion, p > 0 is an integer number determining the shape of the boundary1, t is the time and a is the amplitude of the boundary oscillation. Figure 1 shows a typical scenario of the boundary and three subsequent collisions illustrating the dynamics. The dynamics of the particle is given in terms of a dis- crete, nonlinear four-dimensional mapping A: R4 → R4, that transforms the dynamical variables at collision n to their new values at collision n + 1, (θn+1, αn+1, Vn+1, tn+1) = A(θn, αn, Vn, tn), where Vn = |�Vn| denotes the magnitude of the particle velocity after collision n, and αn corresponds to the angle between the particle 1Non-integer numbers produce open billiard leading to escape of particle through hole on the border. 60003-p1 Matheus Hansen et al. Fig. 1: (Color online) Three consecutive collisions of a parti- cle in a deterministic time time-dependent billiard (red trajec- tory). Illustration of the reflection angle range for a billiard with random motion in the boundary (blue). The relevant an- gles and velocity are shown for the n-th collision. trajectory and the tangent line at the n-th collision with the boundary at the polar angle θn and collision time tn (see fig. 1). Given that each particle moves along a straight line between collisions and with constant speed, the radial position of the particle is given by Rp(t) =√ X2(t) + Y 2(t), where X(t) and Y (t) are the rectangu- lar coordinates of the particle at the time t. The angu- lar position θn+1 is obtained by the numerical solution of Rp(θn+1, tn+1) = Rb(θn+1, tn+1). The time at colli- sion n + 1 is given by tn+1 = tn + √ ΔX2 + ΔY 2/| �Vn|, where ΔX = Xp(θn+1, tn+1) − X(θn, tn) and ΔY = Yp(θn+1, tn+1) − Y (θn, tn). The conservation of momen- tum law is applied in a non-inertial frame where the boundary is at rest. The reflection laws at the instant of collision are �V ′ n+1 · �Tn+1 = ξ �V ′ n · �Tn+1, and �V ′ n+1 · �Nn+1 = −κ �V ′ n · �Nn+1, where the unit tangent and nor- mal vectors are �Tn+1 = cos(φn+1)̂i + sin(φn+1)ĵ, �Nn+1 = − sin(φn+1)̂i + cos(φn+1)ĵ, φ = arctan(Y ′(t)/X ′(t)) with Y ′(t) = dY/dθ, X ′(t) = dX/dθ and ξ, κ ∈ [0, 1] are the restitution coefficients for the tangent and normal direc- tions. The term �V ′ corresponds the velocity of the particle measured in the non-inertial reference frame. The tangen- tial and normal components of the velocity after collision n + 1 are �Vn+1 · �Tn+1 = ξ�Vn · �Tn+1 + (1 − ξ)�Vb · �Tn+1, (1) �Vn+1 · �Nn+1 = −κ�Vn · �Nn+1 + (1 + κ)�Vb · �Nn+1, (2) where �Vb(tn+1) = dRb(t) dt |tn+1 [cos(θn+1)̂i + sin(θn+1)ĵ ] is the velocity of the moving boundary at time tn+1. The speed of the particle after collision n+1 is given by Vn+1 =√ (�Vn+1 · �Tn+1)2 + (�Vn+1 · �Nn+1)2, while the angle αn+1 is written as αn+1 = arctan[ �Vn+1· �Nn+1 �Vn+1·�Tn+1 ]. Given an initial condition (θn, αn, Vn, tn), the dynamical properties of the system can be obtained through appli- cation of the previous equations. We are interested in the behavior of the average speed, obtained from two different kinds of average, namely V = 1 M ∑M i=1 1 n+1 ∑n j=0 |�V |i,j , where the first summation is made over an ensemble of Fig. 2: (Color online) Plot of V vs. n for different initial speeds, as labeled in the figure. Three different regimes are clear in the figure. For large initial speeds, a plateau dominates the dynamics for short n. After a first crossover, the average speed starts to grow as a power law diffusing the speed as a normal diffusion with slope of 0.481(9). Soon after, there is a second crossover where the normal diffusion is replaced by a super diffusion with slope 0.962(6). The right panels show portions of a single realization in the (V, t)-space where A identifies normal diffusion and B super diffusion. The parameters used are ε = 0.08, a = 0.5 and p = 3. different initial conditions randomly distributed in t ∈ [0, 2π], α ∈ [0, π] and θ ∈ [0, 2π] for a fixed initial speed while the second summation corresponds to an average made over the orbit, hence in time. We considered M = 5000 and n = 107 collisions. Figure 2 shows the behavior of V vs. n for different initial speeds. Three different kinds of behavior can be seen from the figure. If the initial speed is large, the curve of average speed exhibits a plateau. The size of the plateau depends on the value of the initial speed (see ref. [10] for a discussion of such kind of behavior in a two-dimensional mapping). A higher initial speed, leads to a longer plateau. A first crossover is observed chang- ing the behavior of the constant regime to the regime of growth with a slope of growth typical of normal diffusion. A numerical fitting gives a slope 0.481(9). The regime of normal diffusion then reaches a second crossover pass- ing to a faster regime of growth named as super diffusion, with slope 0.962(6). The panels on the right-hand side of fig. 2, give the plot of a single realization in the (V, t)- space where A corresponds to normal diffusion and B to super diffusion. We emphasize when the perturbation on the boundary is random, the dynamics in the (V, t)-space is similar to what is observed in A with a constant slope of growth about 1/2 therefore characterized by normal diffu- sion. When the restitution coefficients ξ, κ �= 1, inelastic collisions occur leading to a different scenario [14], where the energy growth is interrupted by the violation of Liou- ville’s theorem and attractors emerge in the phase space. Provided that all particles start with the same speed at a random location and direction they can gain or lose en- ergy upon collisions with the moving boundary and follow very different paths. At a given collision of the full en- semble with the wall, there is roughly the same number of particles gaining and losing speed, then, the mean speed will not change much after the collision [14]. 60003-p2 Explaining a changeover from normal to super diffusion in time-dependent billiards Fig. 3: (Color online) Plot of the evolution of the speed dis- tribution function ρn(V ) for an ensemble of 5000 particles af- ter different number of collisions. The parameters used are ε = 0.08, a = 0.5 and p = 3. However, this process, can not continue indefinitely, for there is a lower bound of the speed, i.e., zero, and no upper bound. Particles reaching the lower bound can only gain speed upon collision with the wall. When a sufficient number of particles reaches the lower bound, the symme- try of the speed diffusion is broken and we have a crossover between the constant regime and the growth regime. Figure 3 shows the evolution of the speed distribution function ρn(V ), for an ensemble of particles as the number of collisions increases. Notice that, after ten collisions the distribution has not reached the lower bound and the mean speed is essentially the initial value of V0. This changes after the lower bound is reached, leading the mean to grow as observed after 100 and 1000 collisions. The speed of a particle after collision with the wall Ṽ , can be written in the form Ṽ (α, θ, t, V ) = V + ζ(α, θ, t, V ), (3) where V is the speed just before the collision. From this, the mean speed of an ensemble of particles just after the n-th collision takes the form V n+1 = V n + δV n, (4) where V n = ∫ ∞ 0 ∫ 2π 0 ∫ 2π 0 ∫ π 0 V Fn(α, θ, t, V )dαdθdtdV, (5) δV n = ∫ ∞ 0 ∫ 2π 0 ∫ 2π 0 ∫ π 0 ζ(α, θ, t, V )Fn(α, θ, t, V ) ×dαdθdtdV, (6) and Fn(α, θ, t, V ), is the phase space distribution function just before collision n. In the case of interest the distribu- tion function can be factored in the form Fn(α, θ, t, V ) = F (θ, α)ρV (t)ρn(V ), (7) Fig. 4: (Color online) Plot of the numerical distribution func- tions ρθ(θ) = R F̃ (θ, α)dα and ρα(α) = R F̃ (θ, α)dθ for deter- ministic and non-deterministic (random) billiards at various amplitude of oscillations. The distributions are mostly deter- mined by the geometry of the static billiard (a = 0), and are weakly modified by the wall oscillation. The inhomogeneous nature of these distributions are due to the presence of low pe- riod saddles in the phase space of the static oval billiard, where individual orbits spend longer times. The control parameters are ε = 0.08 and p = 3. where F (α, θ) is the α − θ distribution, ρV (t) is the colli- sion time distribution and ρn(V ) is the speed distribution function. The F (α, θ) distribution is mainly determined by the billiard geometry and is only weakly modified by the wall oscillation, consequently it can be regarded as in- dependent of the index n, the collision time distribution ρV (t) depends on the speed but not on the index n and the speed distribution function ρn(V ) depends on both the speed and the index. To understand the evolution of the mean speed of an ensemble we do not require to determine the evolution of the global distribution function. Inserting eq. (7) in eq. (5), and defining the partial mean U(V ) = ∫ 2π 0 ∫ 2π 0 ∫ π 0 ζ(α, θ, t, V )F (θ, α)ρV (t)dαdθdt, (8) we obtain a compact expression for the change in the mean speed δV n = ∫ ∞ 0 ρn(V )U(V )dV. (9) Now, consider a second-order expansion of the partial mean U(V ) around the mean speed V n, of the distribution ρn(V ) U(V ) ≈ U(V n) + U ′(V n)(V −V n) + 1 2 U ′′(V n)(V − V n)2. (10) 60003-p3 Matheus Hansen et al. Replacing this expression of U(V ) in eq. (9) we obtain a second-order approximation of the change in the average speed δV n = U(V n) + 1 2 U ′′(V n)(V 2 n − Vn 2 ). (11) The approximation in eq. (10) becomes poor as we move far from the distribution mean. However, the integrand in eq. (9) contributes less where the Taylor expansion of U(V ) is not accurate, because the speed distribution ρn(V ) drops for large and small values of V . However, this argument, becomes weaker when the speed distribu- tion evolves in time to become wider, as observed in fig. 3. Consequently, the precision of this approximation depends on the particular form of U(V ). In the appendix of this manuscript, it is shown that, for thermal speed distributions [15] (similar to those in fig. 3), the relative error of eq. (11) is smaller than an upper bound independent of the mean speed. It is also shown that the upper bound is small for functions U(V ) falling slower than 1/V or growing slower than V 3. Which is satisfied for the functions involved in the following calculations. Replacing eq. (11) in eq. (4) we obtain a second-order approximation of the n + 1 mean speed V n+1 = V n + U(V n) + 1 2 U ′′(V n)(V 2 n − Vn 2 ), (12) which depends on the mean of the quadratic speed at col- lision n. An equation for the evolution of the quadratic mean is also needed. For the quadratic speed after colli- sion Ṽ 2, the collision rule can be also written in the form Ṽ 2(α, θ, t, V ) = V 2 + ψ(α, θ, t, V ), (13) where V is again the speed before collision. Reproducing the same arguments used for the mean speed we obtain V 2 n+1 = V 2 n + W (V n) + 1 2 W ′′(V n)(V 2 n − Vn 2 ), (14) where W (V )= ∫ 2π 0 ∫ 2π 0 ∫ π 0 ψ(α, θ, t, V )F (θ, α)ρV (t)dαdθdt. (15) As will be shown later the inhomogeneous nature of F (θ, α) is a fundamental aspect of super diffusion. It can be related to the presence of low period saddles in the static billiard [16], where the collisions occur more of- ten [17], leading to an increase in the distribution value. Figure 4 we show line-integrated profiles of F (α, θ) for both deterministic and random oscillations [14] of the bil- liard boundary. In the limit of high speeds or small boundary oscillations it can be shown that ζ(α, θ, t, V ) ≈ ψ(α, θ, t, V )/2V , which upon integration in {θ, α, t} leads to U(V ) ≈ W (V )/2V . (16) Fig. 5: (Color online) Plot of the numerical collision time distribution functions ρV (t) for deterministic ((a), (c)) and non-deterministic (random) ((b), (d)) boundary oscillations at various amplitudes and two different initial speeds. For a low initial speed, both cases of the billiard initially have a similar behavior for the collision time distribution ((a), (b)) leading the systems to exhibit the normal diffusion. However, when the initial speed is high, the deterministic billiard exhibits a correlation between the subsequent collisions and consequently an inhomogeneous collision time distribution (c), while for the random billiard such correlation is not observed (d). The cor- responding parameters are ε = 0.08 and p = 3. This results in an approximated form for the two- dimensional mapping of the mean speed and the mean quadratic speed V n+1 = V n + 1 2 WnV 2 n /Vn 3 + 1 2 ( 1 2 VnW ′′ n − W ′ n ) × ( V 2 n /Vn 2 − 1 ) , V 2 n+1 = V 2 n + Wn + 1 2 W ′′ n (V 2 n − Vn 2 ), (17) where Wn = W (Vn),W ′ n = W ′(Vn) and W ′′ n = W ′′(Vn). This mapping is general and its behavior depends on the particular form of the function W (V ) for the system under consideration. For instance, for the time-dependent oval billiard one can show that ψ(α, θ, t) = 4(aε)2 cos2(pθ) sin2(t) − 4aεV cos(pθ) sin(α) sin(t), (18) and the collision time distribution can be approximated by (see fig. 5(a), (c)) ρV (t) = 1 2π [1 − aεχ(V ) sin(t)] , (19) where χ(V ) is a slowly changing function of V that van- ishes for V = 0 and saturates at χ∗ for large V . This distri- bution develops due to the correlation between subsequent 60003-p4 Explaining a changeover from normal to super diffusion in time-dependent billiards collisions for higher speeds for which the time between collisions is small compared to the wall oscillation period. Expectedly, the harmonic part is removed when the wall oscillations are random (see fig. 5(b), (d)). Inserting ψ(α, θ, t) and ρV (t) in eq. (15) and defining the following constants: η1 = (aε)2 ∫ π 0 ∫ 2π 0 F (α, θ) sin(α) cos(pθ)dθdα, (20) η2 = (aε)2 ∫ π 0 ∫ 2π 0 F (α, θ) cos2(pθ)dθdα, (21) we obtain W (V ) = 2η2 + 2η1χ(V )V, (22) which inserted in the mappings for the average speed and the average quadratic speed (17), results in two coupled difference equations Vn+1 − Vn = η1χn + η2V 2 n /Vn 3 , (23) V 2 n+1 − V 2 n = 2η2 + 2η1χnVn, (24) where χn = χ(Vn), and it was used that χ(V ) changes slowly with the speed. More specifically the system sat- isfies χ′(Vn)Vn � χ(Vn). These coupled difference equa- tions can be solved asymptotically for small and large n to give us a picture of the different diffusion regimes ex- hibited by the system during its evolution. If the ensemble of particles begins with small speeds, i.e., of the order aε, we can use χn → 0, and the system can be integrated by taking the continuous limit fn − fn−1 ≈ df(n)/dn, which results in V 2 s (n) = V 2 0 + 2η2n, (25) where the subscript s indicates small n. This solution can be replaced in the mean speed difference equation to obtain an approximated solution for the mean speed valid for few collisions Vs(n) = ( V 2 0 + 2η2n )1/2 . (26) Notice that this solution emerged after assuming an ho- mogeneous phase distribution, i.e., χn → 0, which is also appropriate when the collision phase with the wall is ran- dom. As the mean speed of the ensemble grows by nor- mal diffusion, the time between consecutive collisions is reduced and the collision time distribution becomes inho- mogeneous due to the coherent acceleration and deceler- ation of individual particles (cf. fig. 2, panel B). Leading naturally to a situation where χn �= 0. In the new regime, the value of χn saturates to χ∗ when large speeds are achieved. Then, provided that V 2 n and Vn 2 are of the same order, we have that V 2 n /Vn 3 → 0 for large Vn. The last statement is valid for unimodal distributions, like the Maxwell-Boltzmann one [15] for which V 2 n /Vn 3 = 4/(πV ), but is also valid for delocal- ized situations, like the uniform distribution for which V 2 n /Vn 3 = 4/(3V ). This guarantees that the ratio V 2 n /Vn 3 goes to zero as the mean velocity grows regardless of how localized is the speed distribution. From this, in the deterministic high velocity limit, eq. (23) becomes Vn+1 − Vn ≈ η1χ ∗, (27) which can be integrated to obtain the evolution rule for the high speeds regime Vl(n) = V0 + η1χ ∗n. (28) Here, the subscript l indicates large n. For regular values of η1, η2 and χ∗, the difference Vl(n) − V0 is small com- pared to Vs(n)−V0 for small n, while the opposite happens for large n. Consequently, we can combine eq. (26) and eq. (28) to obtain a single approximate solution valid for all stages of the ensemble evolution V (n) = (V 2 0 + 2η2n)1/2 + η1χ ∗n. (29) An interesting feature of this solution is that η1 vanishes if F (α, θ) is homogeneous, so that, a deterministic billiard without low period saddles in phase space will only ex- hibit normal diffusion because of the uniform distribution of particles in the phase space. The combined solution eq. (29), corresponds to the continuous lines in fig. 2 in excellent agreement with the numerical simulations for all the ensembles considered. As a short summary, in this letter we have shown that an inhomogeneous particle distribution function on the phase space of the static billiard leads to the development of anomalous diffusion regimes in time-dependent situa- tions, and for the particular case of oval billiards, explains the transitions from normal to super diffusion. These re- sults were obtained by identifying the asymptotic dynam- ics of a two-dimensional mapping for the first and second moments of the speed distribution. The presented treat- ment, is sufficiently general to study other anomalous dif- fusion regimes, diverse billiard shapes and more general mappings. ∗ ∗ ∗ MH thanks CAPES for financial support. DC, ILC and EDL thank CNPq (433671/2016-5, 300632/2010-0 and 303707/2015-1) and ILC and EDL thank the São Paulo Research Foundation (FAPESP) (2011/19296-1 and 2017/14414-2) for their financial support. Appendix In this section we are interested in providing a more compeling physical argument for the validity of the ap- proximations used in eq. (11) and eq. (14). First of all consider a large collection of non-interacting particles con- fined in a two-dimensional region in equilibrium with a thermal bath. 60003-p5 Matheus Hansen et al. The probability of a particle having a velocity �V = (Vx, Vy) is proportional to the Boltzmann Factor exp[− m 2kBT (V 2 x +V 2 y )], from which, the probability of hav- ing a speed V , with arbitrary direction is given by ρMB(V , V ) = π 2V 2 V exp [ −π 4 ( V V )2 ] . (A.1) Here, the mean speed V = √ πkBT/2m was used to char- acterize the speed distribution instead of the temperature. This distribution has very similar features to those ob- served in the numerical simulations (fig. 3), but corre- sponds to a system in thermodynamic equilibrium, which is not our case. However we can consider that our system is approximately described as a sequence of equilibrium configurations, where each collision of the ensemble with the moving boundary delivers a bit of energy to the en- semble, changing its temperature (or mean speed in our case). Now, let us consider a given function of the speed with the form Uλ(V ) = V λ, with λ some real number. The mean value of this function under ρMB can be obtained by direct integration and has the form Uλ(V ) = ∫ ∞ 0 ρMB(V , V )V λdV = λ2λ−1π−λ/2Γ(λ/2)V λ = cλV λ . (A.2) Now, let us approximate the value of Uλ with the second- order expansion Ũλ, obtained by expanding Uλ(V ) about the mean of the distribution, which dispenses the integra- tion step, i.e., Ũλ(V ) = Uλ(V ) + 1 2 U ′′ λ (V )(V 2 − V 2 ). (A.3) Here, we can use (V 2 − V 2 ) = (4/π − 1)V 2 , valid for the Maxwell-Boltzmann distribution and U ′′ λ (V ) = λ(λ − 1) V λ−2, so that we have the approximated solution Ũλ(V ) = ( 1 + 1 2π λ(λ − 1)(4 − π) ) V λ = c′λV λ . (A.4) Notice that eq. (A.2) and eq. (A.4) have the same depen- dence in V , then, the relative error of the approximation ελ = |U(V ) − Ũ(V )|/|U(V )| is independent of the mean speed of the distribution. Namely, ελ takes the form ελ = |cλ − c′λ| |cλ| , (A.5) where the coefficients cλ and c′λ are independent of V . Figure 6, we show the relative error of Ũλ for −1 < λ < 3, under a Maxwell-Boltzmann distribution. The difference is below 1% for −0.2 < λ < 2.4, and grows up to ∼10% at λ = −1 and λ = 3. This result suggests that, for speed distributions with thermal characteristics, the second-order approximation Fig. 6: Relative error of eUλ(V ) for −1 < λ < 3. Outside this domain the relative error grows above 10%. presented in eq. (11) is appropriate if the function U(V ) grows or decays slowly with the speed V . For the context of this work, the previous analysis can be extended to functions of the form U(V ) = N∑ i=0 aiUλi (V ) = N∑ i=0 aiV λi , (A.6) where {λi}N i=1 is a set of real numbers. For this function we have the exact and approximated means U(V ) = N∑ i=0 aicλi V λi , Ũ(V ) = N∑ i=0 aic ′ λi V λi , (A.7) with cλi and c′λi given by eq. (A.2) and eq. (A.4). Here, the Cauchy-Schwarz inequality ensures that |U(V ) − Ũ(V )| ≤ ∣∣∣∣∣ N∑ i=0 cλi − c′λi cλi ∣∣∣∣∣ ∣∣∣∣∣∣ N∑ j=0 ajcλj V λj ∣∣∣∣∣∣ . (A.8) Dividing both sides by |U(V )| and using the triangle in- equality in the r.h.s., we obtain an upper bound for the approximation error of the function Ũ(V ), ε[Ũ ] = |U(V ) − Ũ(V )| |U(V )| ≤ N∑ i=0 ελi , (A.9) where the ελi are as given in eq. (A.5) and deppend on λi as depicted in fig. 6. Consequently, the upper bound of the approximation error of Ũ(V ) is independent of the mean speed V , which guarantees the stability of the ap- proximated mapping in eq. (17) for speed distributions with thermal characteristics. REFERENCES [1] Fermi E., Phys. Rev., 75 (1949) 1169. [2] Chernov N. and Markarian R., Chaotic Billiards (American Mathematical Society) 2006. [3] Berry M. V., Eur. J. Phys., 2 (1981) 91. [4] Sinai Ya. G., Russ. Math. Surv., 25 (1970) 137. [5] Lenz F., Diakonos F. K. and Schmelcher P., Phys. Rev. Lett., 100 (2008) 014103. 60003-p6 Explaining a changeover from normal to super diffusion in time-dependent billiards [6] Leonel E. D. and Bunimovich L. A., Phys. Rev. Lett., 104 (2010) 224101. [7] Loskutov A., Ryabov A. B. and Akinshin L. G., J. Phys. A: Math. Gen., 33 (2000) 7973. [8] Gelfreich V. and Turaev D., J. Phys. A: Math. Theor., 41 (2008) 212003. [9] Gelfreich V., Rom-Kedar V. and Turaev D., Chaos, 22 (2012) 033116. [10] Oliveira D. F. M., Silva M. R. and Leonel E. D., Physica A, 436 (2015) 909. [11] Leonel E. D., Oliveira D. F. M. and Loskutov A., Chaos, 19 (2009) 033142. [12] Batistić B. and Robnik M., J. Phys. A, 44 (2011) 365101. [13] Batistić B., Phys. Rev. E., 90 (2014) 032909. [14] Leonel E. D., Galia M. V. C., Barreiro L. A. and Oliveira D. F. M., Phys. Rev. E, 94 (2016) 062211. [15] Reif F., Fundamentals of Statistical and Thermal Physics (Waveland Press) 2009. [16] Hansen M., de Carvalho R. E. and Leonel E. D., Phys. Lett. A, 380 (2016) 3634. [17] Hansen M., da Costa D. R., Caldas I. L. and Leonel E. D., Chaos, Solitons Fractals, 106 (2018) 355. 60003-p7