Oscillating mushrooms: adiabatic theory for a non-ergodic system

Can elliptic islands contribute to sustained energy growth as parameters of a Hamiltonian system slowly vary with time? In this paper we show that a mushroom billiard with a periodically oscillating boundary accelerates the particle inside it exponentially fast. We provide an estimate for the rate of acceleration. Our numerical experiments confirms the theory. We suggest that a similar mechanism applies to general systems with mixed phase space.


Introduction
Consider a particle which moves freely inside a bounded domain (a billiard) and reflects elastically from the domain's boundary. We assume that the boundary of the billiard domain changes in time and restores its shape periodically. Elastic collisions with the moving boundary change the particle kinetic energy. If the average energy gain over multiple collisions is positive, the particle accelerates. This process is often called "Fermi acceleration" since it resembles the mechanism of cosmic particles acceleration by reflecting from magnetic mirrors proposed by Fermi in 1949 [12]. In recent years the problem has attracted much attention (see e.g. [27,11,23,31,4,14,32] and references therein).
When the particle accelerates indefinitively, the motion of the boundary eventually becomes slow compared to the particle motion, so the separation of time scales will be our standing assumption. From the physical point of view, it is important that the motion of the billiard boundary is not affected by the collisions with the particle, i.e. one may think acceleration in periodically perturbed ergodic billiards.
However, in [31,13] we found that ensemble-averaged Fermi acceleration can be much faster in situations where the chaotic billiards in the frozen family have, for some intervals of the parameter t, more than one ergodic component. In this case most of the initial conditions lead to an exponential energy growth. In other words, dividing and uniting chaotic ergodic components of the frozen billiards can lead to energy transfer which is much more effective than in the case of a single ergodic component.
This observation is of great importance as Hamltonian dynamical systems are rarely ergodic on every energy level. Moreover, it is widely believed that the majority of Hamiltonian systems have a mixed phase space, a rather complicated mixture of chaotic and regular zones. The latter are formed by Cantor-like families of KAM (Kolmogorov-Arnold-Moser) invariant tori. The ergodic adiabatic theory is not applicable in this case, which opens the possibility of exponential Fermi acceleratetion for the majority of initial conditions.
In this paper we further explore the accelerating effect of the violation of ergodicity: We consider an example of a slow-fast system whose fast subsystem has a chaotic set coexisting with an elliptic island filled by invariant tori. We show that the transitions of the fast variables between the chaotic and elliptic zones of the frozen system caused by the slow changes of the billiard shape play a crucial role in the dynamics. Moreover, the effect of these transitions can be quantified and estimated. We believe our analysis should be applicable to a much wider class of slow-fast Hamiltonian systems whose fast subspace contains coexisting chaotic and elliptic components.
The Anosov-Kasuga invariant can be related to the Boltzmann entropy of the system. This means that adiabatic changes of parameters of an ergodic Hamiltonian system do not lead to a systematic increase of the entropy for a very long time. We show that this conclusion does not hold when the frozen system has a mixed phase space. Transitions between the chaotic set and elliptic island break the Anosov-Kasuga adiabatic conservation law, which leads to a systematic increase of the entropy over each period of the slow perturbation. This also means a steady (exponential in time) increase in the energy.
Below we build a quite detailed analysis of this process for a particular type of a planar billiard with periodically moving boundaries, an oscillating Bunimovich mushroom. The frozen billiard shape is shown in Fig. 1. This shape was invented in [7]; the corresponding billiard has a phase space that is sharply divided into a single elliptic and a single chaotic component, each one of positive measure [7,8]. One can easily find explicit expressions for the volumes of the regular and chaotic components. Using this data we propose analytical expressions for the energy distribution after one cycle of the billiard boundary oscillation and predict the energy growth rate. We show that provided the billiard boundary moves along a non-trivial loop in the space of the billiard parameters, a particle inside the billiard accelerates exponentially fast.
Essentially, we show that the process of energy growth over many cycles of boundary oscillations can be modeled by a geometrical Brownian motion, where the particle energy after each cycle is multiplied by an independent random variable. We obtain a formula which evaluates the expectation of the logarithm of this random factor in terms of the volumes of h w θ r d ′ w Figure 1: Bunimovich mushrooms the chaotic and elliptic components, and show that this expectation is always non-negative (and typically strictly positive, which immediately implies the exponential energy growth). The numerical experiments confirm the predicted growth rate with a quite good precision.
We stress that our motivation goes beyond the billiard realm -we believe that for a Hamiltonian system of arbitrary nature whose chaotic behaviour is "spoiled" by the presence of KAM-islands in the phase space a slow change of parameters will typically lead to the steady increase in energy over each period of the parameter variation. In particular, an example of exponential energy growth for a particle in a polynomial potential with periodically changing coefficients is studied in the forthcoming paper [28].
The paper is ordered as follows. In section 2 we introduce geometrical parameters that govern the shape of the Bunimovich mushroom [7] and find explicit expressions for the dependence of the volumes of the regular and chaotic components on these parameters. Section 3 contains our main theoretical developments. We describe there how phase space volume is transferred between the components when the mushroom boundary slowly oscillates, and how the oscillations change the particle energy. First, we derive an adiabatic theory in the presence of particle flux, enabling us to calculate the change in energy of a particle while it stays in the chaotic or regular component. We then find the probability of capture into the regular component. Then, we calculate the averaged growth rate of energy over a period of the mushroom oscillations, and show that this rate is always non-negative. In fact, it is strictly positive for generic oscillations, thus most particles inside the billiard accelerate exponentially fast. In section 4 we present several numerical experiments that confirm the predicted growth rate, the predicted capture/release process, and the predicted energy distribution after one oscilaltion period with a quite good precision. In section 5 we summarize the work and provide several insights regarding the transfer of energy from slow to fast systems.

Frozen mushroom
A Bunimovich mushroom consists of a semi-disk and a stem [7] (see Figure 1). A particle moves freely inside the mushroom D and reflects elastically from its boundary. The particle can go from the cap to the stem and back through a hole. The dynamics of the particle is defined by the Hamiltonian where p = (p x , p y ) and x = (x, y) are respectively the momentum and position of the particle. On the boundary of D the particle velocity is reflected according to the elastic low: the incidence angle is equal to the reflection angle and the absolute value of the velocity is not changed. The phase space of the Bunimovich mushroom consists of two invariant components of positive Lebesgue measure: a regular component (an elliptic island), filled by invariant tori, and a chaotic component, which is ergodic and mixing [7,8]. The regular component consists of points trapped in the mushroom cap.
Let r and w denote the radius of the cap and the semi-width of the hole respectively (see Figure 1 (left)). When the particle moves inside the cap, the incidence angle ϕ with the semi-circular part of the boundary remains constant. So if at some collision | sin ϕ| ≥ ν = w r , then the particle cannot reach the hole and consequently never leaves the cap. Along these trajectories the absolute value of the angular momentum, |xp y − yp x |, is preserved, and consequently the dynamics is integrable. The chaotic component is the complement to the integrable one [8].
In this paper we will consider a tilted Bunimovich mushroom in the domain D = D stem ∪ D cap defined by D stem = |x| ≤ w + y tan θ and − h ≤ y ≤ 0 , D cap = x 2 + y 2 ≤ r and y ≥ 0 .
We assume that w ≤ r and |θ| < π 2 . For simplicity of presentation we assume that for positive θ the cone extends up to y = h, namely, that w ≥ h tan θ. At θ = 0 we obtain the original Bunimovich mushroom (which has a family of parabolic orbits formed by horizontal oscillations inside the stem, whereas at nonzero θ the family of parabolic orbits is destroyed). It can be shown by methods of [6,8] that the tilted mushroom has, for every θ, the same property as the non-tilted one: the complement to the set of points whose orbits are forever trapped in the cap is ergodic and mixing.
Consider a set D 0 ⊂ D and let A(D 0 ) denote its area in the configuration space. Let V (D 0 ) be the phase space volume of the set defined by H(x, y, p x , p y ) = 1 2 and (x, y) ∈ D 0 . Obviously, the volume of all points in the phase space such that H(x, y, p x , p y ) = E and (x, y) ∈ D 0 is equal to 2π √ 2EA(D 0 ), so V (D 0 ) = 2πA(D 0 ). In particular, we find that For the future analysis we will need to find volumes of the integrable and chaotic components. The elliptic island (integrable component) occupies only a part of the cap. Let V ell be the phase space volume of the elliptic island at the energy level H = 1 2 . We claim that V ell = δ(ν)V cap (r) (2) where ν = w r and Indeed, consider a point (x, y) = (ρ cos ψ, ρ sin ψ) inside the cap with a velocity vector (p x , p y ). Then define the angle φ by (p . The absolute value of angular momentum is equal to vρ| sin φ|. Notice that at the cap boundary the angle φ coincides with the incidence angle ϕ. The integrable component consists of trajectories that never cross the stem-cap boundary {ψ = 0, ρ < w} ∪ {ψ = π, ρ < w}, where w is the half width of the hole. Since the angular momentum is conserved, and it is smaller than wv on this boundary, the integrable component at v = 1 is defined by the inequalities: Introducing s 1 = ρ cos ϕ, s 2 = ρ sin ϕ, we get so the volume is given by the following integral Dividing by V cap and using (1), we obtain (2). The chaotic component is the complement of the regular one, therefore its volume is given by

Adiabatic oscillations and capture into an island
Now suppose that the parameters of the mushroom change slowly with time. The velocity of the particle is no longer preserved and we study the time evolution of the particle's energy. The classical adiabatic theory describes dynamics for a slowly changing integrable system [2,3,25], and the ergodic adiabatic theory studies the case when the system is ergodic at every energy level for every frozen moment of time [1,21]. Neither of these theories is directly applicable when the phase space of the frozen system is a mixture of integrable and chaotic components as the particle can transfer from one type of motion to another one due to the changes in the system.

Dynamics in an oscillating circular billiard
While the particle stays in the cap, its dynamics can be described by a circular billiard of the same radius. We assume that the radius is a given smooth function of time, r(t), and the center of the circle does not move. Let the particle hit the boundary at a point P . Let ϕ be the impact angle, i.e., the angle between the pre-collision velocity and the external normal to the boundary at the point of collision. Let v = v sin ϕ denote the component of the particle velocity parallel to the boundary, and v ⊥ = v cos ϕ be the normal component of the velocity. The elastic reflection , this gives the standard elastic law in the coordinate frame that moves with the boundary). Thus, for the outgoing angleφ we have tanφ = sin ϕ We assume that the particle moves much faster than the boundary, i.e. u v. The particle speed after the collision is given bȳ The time interval ∆t to the next collision is Lv −1 where L = |P Q| is the distance to the next collision point Q. Since this interval is O(v −1 )-small, the change in the circle radius is also O(v −1 ), so L is O(v −1 )-close to the value it takes in the static circular billiard of radius r(t), i.e. to the length of the chord that makes the angle ϕ to the radius. This gives us The new value of the radius at the moment of the next collision is Now, by considering the triangle P QO where O is the center of the circle, we find where ϕ is the impact angle at the collision point Q. By (6),(7), we obtain It follows immediately that the impact angle ϕ stays approximately constant over time required for at least O(v) collisions, i.e. it is an adiabatic invariant.
Since the circular billiard keeps being rotationally symmetric even when its radius oscillates, the angular momentum is preserved: v r sin ϕ = vr sin ϕ.
Taking into account equation (8) It follows immediately that the product v(t)r(t) stays approximately constant over time required for O(v) collisions. Taking the square, we see that In terms of the billiard flow, the adiabatic invariance of the impact angle ϕ can be expressed as the adiabatic invariance of the angleφ defined as Indeed, at the moments of collisionφ coincides with ϕ, and it stays approximately constant between impacts. In the oscillating mushroom the above equations are valid for segments of trajectories located entirely inside the mushroom cap.

Capturing into the cap
Suppose that a fast particle moves inside the non-autonomous mushroom. While the particle moves inside the mushroom cap, its dynamics follows the laws of the circular billiard but the sign of its angular momentum is reflected each time the particle hits the bottom of the cap. The hole of the mushroom cap is characterised by the dimensionless parameter Suppose the particle goes through the hole at a moment t = t 0 . At this moment it is located at a point (x, 0) and has velocity (p x , p y ). Let sin ϕ 0 = | sinφ(t 0 )| be the absolute value of the adiabatic invariant defined by equation (9). Since |x| ≤ w(t) and |p y | ≤ √ 2E we conclude that sin ϕ 0 ≤ ν(t).
Depending on the initial conditions, the particle may enter the cap through the hole at t = t 0 and have multiple consecutive collisions inside the cap. While the particle moves inside the cap, | sinφ(t)| is adiabatically invariant. So if the function ν decreases below sin ϕ 0 while the particle remains inside the cap, the particle cannot reach the hole and leave the cap till ν returns back to sin ϕ 0 . We see that the particle can be captured into the cap when ν decreases. If the particle is captured at t = t in , then it is released back into the chaotic zone around the moment t = t out when ν(t out ) = ν(t in ) for the first time. Thus, we can introduce the release function t r (t) so that t out = t r (t in ).

Adiabatic theory in the presence of particle flux
The classical ergodic adiabatic theory [1,21,5] relies on the analysis of evolution of volumes in the phase space. Indeed, the classical ergodic adiabatic theory uses two observations which can be presented in a bit oversimplified form in the following way. First, the ergodicity implies that time averages can be replaced by space averages and, consequently, the time evolution of the energy is the same for the majority of initial conditions starting on a given energy level. Second, if the evolution of the energy depended on initial energy only, the dynamics would map an energy level into another energy level. Since the Hamiltonian flow is volume preserving, the volume under the energy level would stay constant. In a two-dimensional billiard this volume is proportional to EV , which is indeed the ergodic adiabatic invariant.
If the system is not ergodic on energy levels, there is a possibility of the phase space flux to or from an ergodic component of the frozen system. Then volume preservation argument is no longer valid and the adiabatic theory requires a modification.
In the non-autonomous mushroom billiard, the flux between the regular and chaotic zones is governed by the parameter ν(t) defined by (10). Ifν < 0, phase volume "leakes" from the chaotic zone to the regular one, whereasν > 0 leads to the opposite effect.
We analyse this situation by considering short time intervals, [t, t + dt], at which the mushroom's shape does not change noticeably, yet, the particle experiences a large number of collisions with the billiard boundary. Suppose that during this time interval the radius of the hole w has changed by dw, the length of the stem h has changed by dh and the radius of the cap r has changed by dr.
In this discussion we need to distinguish two cases depending on the sign of which determines the direction of the particle flux between the integrable and chaotic components. In particular, if rdw < wdr, the particle in the integrable component cannot leave the cap but the particle in the chaotic zone can be captured into the integrable component. Let us consider the case of capture, dν < 0, in more details. We conjecture that if the particle is sufficiently fast and the changes in the billiard's shape are sufficiently small then, from the statistical point of view, the distribution of the energies at t + dt depends only on the initial and final billiard shapes, i.e., it depends only on the values of dw, dh, dr and dθ, and is independent of the particular form of the evolution of the billiard's shape in the intermediate moments of time. Now, in order to separate the process of capturing into the cap from the adiabatic evolution of the energy, we represent the change of the billiard's shape as a composition of two steps. At the first step, we allow the mushroom to take the intermediate shape shown on Figure 1 (right): we make the hole in the cap slightly narrower by inserting two straight-line segments into the hole that symmetrically extend the cup bottom line.
We use the notation d z and d z to label changes of a variable z during the first and, respectively, the second of these steps (z ∈ {r, w, h, θ}). So, over the time interval dt we get dz = d z + d z.
The shape of the billiard at the end of Stage 1 is uniquely defined by the following requirements: On Stage 1 the particle energy is not changed but all possible transitions between the chaotic and integrable component take place during this step, i.e., d E = 0 and d ν = dν < 0. On Stage 2 the particle cannot move from the integrable to the chaotic component or vice versa. This condition is achieved by the requirement d ν = 0. All changes in the energy are on this step: dE = d E.
On Stage 1 all walls remain static (hence d h = 0, d θ = 0, d r = 0) except for the bottom of the cap, which extends to cover a part of the hole as shown on Figure 1 (right). By the end of this stage the radius of the hole is changed by which is negative since dν < 0 by assumption, and the parameter ν takes its final value: We note that the total phase space volume does not change during this stage but a part of the phase volume is transferred from the chaotic to the regular component: The last equality is a consequence of equation (2) and the requirement d r = 0. It is important to note that on Stage 1 the energy remains constant. Indeed, as the straight-line segments that are inserted into the hole slide along themselves, and the other parts of the billiard boundary do not move during the Stage 1, the normal velocity of the boundary is non-zero only at two points, the end points of these two segments. Therefore, the particle energy can change only if it hits the boundary exactly at one of these points at some moment of time, but this is a probability zero event.
On Stage 2 the billiard is slowly deformed from its intermediate shape to the final one in such a way that the parameter ν = w r remains constant. Namely, w and r are changed simultaneously in such a way that d r = dr and d ν = 0. Then equation (12) implies that d w + d w attains its correct final value dw. At the same time h, θ are changed to ensure d h = dh, d θ = dθ. It follows that on Stage 2 the particle is trapped either in the regular or in the chaotic zone and cannot transfer from one component to the other due to the conservation of impact angles in the circular part.
If the particle is in the regular zone, its dynamics is described by the circular billiard and its energy changes according to the adiabatic law d (EV cap ) = 0 as derived in Section 3.1.
Since both E and V cap are constant on Stage 1, and we get dE = d E and dV cap = d V cap and, consequently, Since the flux is absent the chaotic zone remains invariant for the non-autonomous billiard, and we assume that the standard ergodic adiabatic theory can be applied to the restriction of our dynamical system onto its invariant subset. So if the particle is inside the chaotic zone, then the absence of flux allows us to claim that d (EV cha ) = 0. Using equation (13) we get Since d E = dE, we conclude that for particles in the chaotic zone where V cha is given by (5). Note that in contrast with the classical ergodic adiabatic theory, the right-hand side of this equality contains an additional term which takes care of the phase volume flux from the chaotic zone. One can check that the same equations describe the evolution of the energy in the case when dν > 0 (flux from the integrable to chaotic component). The equations can be derived in a similar way but Stages 1 and 2 are to be swapped.

Probability of capture
In order to describe the acceleration mechanism induced by capture and release by the integrable region we consider an ensemble of non-interacting particles inside the billiard. Then we can discuss the probabilities for a particle to be captured.
For an ensemble of non-interacting particles uniformly distributed inside the chaotic zone we can estimate the probability of capture into the elliptic island due to the decreasing of the quotient ν(t) = w(t)/r(t). Suppose thatν ≤ 0 on a time interval [t 0 , t 1 ] and we start with n cha (t 0 ) particles with initial conditions uniformly distributed inside the chaotic zone. We assume that the distribution of the particles in the chaotic zone remains uniform for all times. Let t ∈ [t 0 , t 1 ). We consider the transfer of particles from the chaotic zone during the time interval [t, t + dt] following the two-steps approximation described in the previous section. On Stage 1, a volume of size is transferred from the chaotic zone to the regular component, and on Stage 2 the number of particles in the chaotic zone remains unchanged. Since the particles are uniformly distributed we get Let p cha (t) = n cha (t)/n cha (t 0 ) be the probability of being in the chaotic zone at the time t.
On the time interval [t 0 , t 1 ], the function p cha satisfies the relation Integrating equation (16) we get for t ∈ [t 0 , t 1 ]. We stress that these equations are valid only during the period of capture when no particles are released from the elliptic zone.

Energy growth rate over a cycle
Suppose that the parameters of the mushroom are periodic functions of time and let T be the period. During this cycle a particle from the chaotic zone may be captured into the island when ν(t) decreases. If the particle is captured at a time t, then it is released at the time t r (t) defined by (11). Note that according to this definition, if ν(t) is increasing at some t, then t r (t) = t. We also note that ν(t) = ν(t r (t)). Since all particles captured between t and t r (t) are released back into the chaotic zone by the time t r (t), we also get n cha (t) = n cha (t r (t)) for all t. We assume that at t = 0 the parameter ν(t) takes its maximum, so all the particles which were captured inside the chaotic zone during the cycle are releases back by the beginning of the next cycle. Let us define a compression factor which describes the change in the relative size of the chaotic zone during the time of capture: (14) and (15) imply that

Let S(t) = log E(t) and S(t) = E[S(t)]. Then equations
The total number of particles is n 0 = n cha + n ell . Consequently, We can rewrite this equation in the following form Integrating over a complete cycle and taking into account that p cha = n cha /n 0 (all particles are in the chaotic zone at the beginning) we get The first term may be integrated by parts: where we grouped together contributions from capture-release pares taking into account that p cha (t) = p cha (t r (t)) and δ(t) = δ(t r (t)). For the second integral we get )dp cha where we used equation (16) which is valid during the capture process. Defining p ell = 1−p cha we get the formula Since p ell (t) is a non-decreasing function of time during the capture and log g < g − 1 for any g > 0, g = 1, we conclude that the energy growth rate m 1 is non-negative. Moreover, it is strictly positive if g(t) ≡ 1. As p ell increases during the capture process, we have shown that for any periodic cycle which is nontrivial, namely, for any cycle having a non-trivial interval of capture with g(t) ≡ 1 on this interval. Since ν(t) = ν(t r (t)) and δ is a function of ν only, we get δ(ν(t)) = δ(ν(t r (t))), i.e., at the moments of capture and release, the regular zone takes the same proportion of the cap's phase space volume. We conclude that g(t) = 1 if and only if V cap (t)/V stem (t) = V cap (t r )/V stem (t r ) (see (18), (5)).
This observation can be restated in the following way. A billiard cycle defines a closed curve V ell (t) where t ∈ [0, T ]. If this curve encloses a non-empty interior, the cycle is nontrivial. On the other hand, any cycle which reduces to a one dimensional curve is trivial. In particular, if one changes only one-parameter of the billiard, the above mechanism of exponential acceleration is not created.
Next, we examine the evolution of the energy distribution after one complete cycle of the billiard boundary. Suppose that a particle is captured at t = t in and let E 1 (t in ) be its energy at the end of the cycle. In the adiabatic approximation its energy can be obtained by integrating equation (19): where the second line uses the definition of the compression factor (equation (18)) and the fact that T 0 d(log V cha ) = 0. Since the probability of capture at t ∈ [t in , t in + dt in ] is equal to −ṗ cha (t in )dt in where p cha is given by equation (17), we obtain the probability distribution of the energy after the cycle in an implicit form. Later in this paper we will use these implicit equations to reconstruct the distribution of the energy for specific examples of the mushroom cycles.
On the other hand if the particle is not captured over the cycle its energy is defined in the adiabatic approximation (see (15)): This equation implies that, if the particle stays in the chaotic zone over the whole cycle, its energy at the end of the cycle does not need to be equal to the initial energy. This conclusion is in a strong contrast with the ergodic case, where the adiabatic theory predicts that the energy returns close to its initial value at the end of a cycle. So these changes are determined solely by the corrections terms in (15), which take into account the phase flux due to the non-ergodicity of the frozen billiards. We check numerically the prediction of equation (22) for several examples of billiard cycles in the next section. Any closed curve in the space of parameters defines two billiard cycles which corresponds to two different direction of motion along the curve. The values of the integrals (22) for these two cycles have the same absolute value but opposite signs. So the energy of the noncaptured particles may increase or decrease after the completion of the cycle, but in both cases the energy averaged over all initial conditions increases.

Examples of billiard cycles
The theory developed in the previous sections relies upon several assumptions, which are difficult to prove analytically. In particular, we assumed that the ergodic averaging theory is applicable to the chaotic component in a system which does not satisfy some of the assumptions of the original averaging theory. The most dramatic violation here is that there are transitions between the ergodic components -we propose that formula (15) replaces the usual adiabatic law. Additionally, the billiard is not a smooth system (see discussion in [14]). Finally, we also assumed here that the distribution of the particle in the chaotic zone of the breathing billiard is close to the stationary uniform one. So we carry out numerical tests to check the correctness of the theoretical predictions for the energy growth rate and for the distribution of the energy after a cycle of the billiard boundary.

Fixed cap
Our first example corresponds to the following protocol for the time dependence of the mushroom parameters: the radius of the gap and the length of the stem follow the straight lines which connect the points on the plane (w, h) while θ is fixed and r = 1 (see Figure 2). For definiteness we assume that w 0 < w 1 and h 0 < h 1 . The particle can be trapped inside the mushroom cap during the first stage of the process, when the hole shrinks, and then it is released during the third stage. Since the radius of the cap is fixed, the particle's energy stays constant while the particle remains in the cap.
The processes of capture and release are determined by the width of the gap: if the particle is captured into the cap at t = t in then it is released at t = t out such that w(t in ) = w(t out ). So we can rewrite equations (20) and (18) E log where where the volumes are defined by (1) and (2). In the derivation of the last equality we take into account that V cap remains constant and V ell = δ(w)V cap . Then p nc = p(w 0 ) is the probability of avoiding the capture completely. If the particle avoids capture, equation (22) implies that the energy after a complete cycle is given by If a particle is captured at the moment when the hole has size w then its energy at the end of the cycle is described by (21), which takes the form These equations can be used to construct a distribution for log E 1 E 0 at the end of cycle. We note that they involve integrals of explicitly defined functions.
Another example is obtained when we reverse the protocol (23) by following the same path in the space of parameters in the clockwise direction: The distribution of the energy and the acceleration rate are described by the same equations but h 0 and h 1 are swapped.
In the first series of experiment we follow the two above described protocols by moving the billiard boundaries with a piece-wise constant acceleration. We generate N = 10 5 initial points uniformly distributed inside the billiard. All initial conditions have the same energy E 0 = 10 9 and a random direction of the initial velocity. Then we follow numerically the trajectory for each of the initial conditions during the time required for one complete cycle of the billiard boundary.
The distributions of the final energy E 1 is described by histograms which represent relative frequency density for log E 1 E 0 . The histograms are shown on Figure 3 where the dashed red line represents the theoretical predictions for the energy distribution and the vertical red line marks the position of the theoretically predicted energy for non-captured particles. In parallel, we mark the capture and release time for each of the initial condition. In this way we test the accuracy of the prediction for the particle flux into the cap and, in particular, for p nc .
We see that the prediction for the energy growth rate is in good agreement with the numerically obtained average growth rate:  Figure 3 consist of two components which correspond to captured and noncaptured particles. The distribution of the energy for the captured particles is apparently in good agreement with the theory. The distribution for the non-captured particles looks like Gaussian and its width scales as E 1/4 0 . This behaviour has been observed for system where the ergodic averaging theory is applied [5,14]. The centre of the distribution is reasonably close to the value predicted by equation (25). We see that the proposed correction to the adiabatic theory (in particular, equation (22)) is realized. We also note that numerical experiments show that the relative frequency of escaping the capture is in excellent agreement with the theoretical prediction given by p nc . Notice that changing the loop direction leads to changing the role of heating in the stem and cap -clockwise motion means that the heating occurs in the cap and the cooling in the stem whereas anticlockwise motion reverses their role. Yet, as predicted, the overall averaged growth rate of energy is positive in both cases.

Example with moving cap
In the second set of examples we change the billiard parameters in the following way: In this cycle all parameters of the billiard (except the slope θ) are changed simultaneously. The capture-release process is determined by ν(t) = w(t) r(t) . Since ν(t) = ν(2π − t) and it is monotonically decreasing on (0, π) there is a simple relation between the time of capture and the time of release: We find the compression factor from equation (18) and the energy growth rate from (20). Then the distribution of the energy after a complete cycle is found from equations (21), (22) and (17). For the purpose of plotting the distributions we evaluated these integrals numerically using the Simpson rule. We compare the theoretical predictions and the results of the numerical computations. We set the amplitudes to the following values The slope of the billiard stem is kept equal to tan θ = 0.1111 during all experiments of the present section. The protocol is illustrated in Figure 4. The right hand side plot shows (clockwise direction). Parameters are r 0 = h 0 = 1, a = 0.5, b = −0.5, c = 0.8. On the right: E 1 (t in ) -theoretical prediction, see equation (21). the predicted energy after a full cycle for a particle captured at time t in , see equation (21). The positivity of energy gain for most of the captured particles corresponds to the right hump of the distribution shown in Fig 5. The negative value for t in = π demonstrate that here the non-captured particles have, on average, a decrease in their energy, indeed, here log E nc 1 E 0 = −0.422465. We select N = 25 000 uniformly distributed initial conditions inside the billiard, which have the same initial energy and randomly chosen initial directions of velocity. The distribution of the energy after one cycle is shown on Figure 5 for two selected values of the initial energy E 0 . The dashed line represents the theoretical distribution. We see that the numerical data are rather close to the theoretical prediction and, as it should be expected, the agreement is better for the higher initial energy. The distribution has two modes. The left mode corresponds to the particles which are not captured into the cap during the cycle. The position of the left mode is close to the theoretical prediction given by log The theoretical prediction for the growth rate is m 1 = 0.122768. During the numerical experiment with E 0 = 10 7 we obtain the average growth rate to be m * 1 = 0.1213 ± 0.0033 (±σ N ) which is in excellent agreement with the theory. We also trace the capture and release times for each of the initial conditions. The distribution of the capture times is illustrated by a histogram shown on Figure 6. The theoretical prediction obtained from equation (17) is plotted using the dashed line. In this experiment 38.78% of particles is not captured in the cap, which is in a good agreement with the probability of non-capture being p nc = 38.474.
We conclude that the theory is overall in good agreement with the data from the numerical experiments, including both the distributions of energy and of time of capture.
We see that after a single cycle the distribution of log E 1 E 0 is quite far from being Gaussian.   The central limit theorem suggest that if the increments of the logarithm of energy are not correlated over consecutive cycles, then the distribution of 1 n log En E 0 should be close to the normal one centered around m 1 for large values of n. Figure 7 represent the distributions for n = 10 and n = 30. It is clearly seen that while the central part of the distribution is close to the predicted Gaussian shape, the tails deviate noticeably from the normal law.
Finally, we stress that the non-ergodicity of the billiard plays the central role in the creation of the exponential acceleration. In order to illustrate this difference we consider a billiard cycle with the same parameters as above but setting c = 0. So while the radius of the cap and the length of the stem are the same as in the previous experiments, the frozen billiard table remains chaotic for each time. Then the ergodic adiabatic theory predicts that the energy should come to its initial value for the majority of the initial conditions. This is confirmed by a numerical experiment: the distribution of the energy after one cycle is shown in Figure 8. We see that the final energy is distributed in a Guassian-like way with quite small standard deviation around zero: |m * 1 | < 3 · 10 −5 and therefore we observe no exponential acceleration on average.

Summary and Discussion
We propose a mechanism for achieving an averaged exponential acceleration rate for majority of initial conditions in a slowly varying system in which the fast dynamics has mixed phase space. The mechanism is examined by considering an oscillating Bunimovich mushroom, where analytical predictions are derived and are corroborated numerically. Specifically, after finding explicit expressions for the volumes of the regular and chaotic components, analytical expressions for the energy distribution after one cycle and of the predicted exponential energy growth rates are derived. Numerical experiments support our predictions for both the energy distribution and the growth rate, and support our claim that the violation of ergodicity is essential for getting exponential acceleration. While we examined time periodic oscillations, the described mechanism does not require precise periodicity of the process -it only assumes that the billiard approximately restores its shape and size from time to time.
To obtain these estimates we develop a generalised ergodic adiabatic theory, which takes into account the presence of flux between different ergodic components. We showed that the averaged exponential growth rate attains the following formula where g denotes the compression factor for the phase volumes between the moments of capture into the elliptic island and release (equation (18)) and p ell denotes the probability of the particle being captured in the elliptic component. This formula is derived under quite general conditions and we hope that it may be applicable to other systems with mixed phase space. Averaged exponential acceleration is achieved when the billiard parameters change along a non-trivial loop, for which the compression factor g of equation (29) is not identically one. A non-trivial loop bounds a non-empty interior when projected on a parameter plane in which one axis corresponds to the phase space volume of the chaotic zone and the other axis corresponds to the flux between the chaotic and integrable components. An important conclusion is that the exponential acceleration rate vanishes if the motion of the billiard boundary can be described by periodic oscillations of a single parameter (as is often done in numerical simulations of Fermi acceleration). Similarly, it vanishes if one of these two ingredients -the flux between the ergodic components or the volume change of the ergodic components -is missing. When the exponential rate vanishes, we still expect to observe some slow acceleration, typically quadratic in time (see e.g. [22,27,23,26]).
We conjecture that the described mechanism, of exponential acceleration due to interior flux and volume changes for different ergodic components of systems with mixed phase space that are adiabatically deformed, is quite general. However, in contrast with the Bunimovich mushroom, in generic billiards and in generic smooth systems, the separation of the frozen fast system into continuously depending on parameter ergodic components is more problematic. The existence of chaotic components with positive phase space volume is unknown, and the numerically observed boundary between the seemingly integrable and chaotic components is often "sticky" and fractal. Nonetheless, we may envision that some rough estimates distinguishing the regular from the chaotic components may be derived (e.g. by calculating Lyapunov exponents), from which the compression rate g and the capture probability p ell may be found. Then, formula (29) may formally connect these geometrical features of the frozen system with the energy growth rate in the adiabatically perturbed system. We should note that the influence of sticky or parabolic orbits on the statistics may be non-trivial. In fact, our initial numerical experiments with the classical Bunimovich mushroom, which has a rectangular stem and thus a family of parabolic periodic orbits, showed that the energy distribution in the chaotic zone was different from the one obtained with tilted geometry. The influence of these effects when multiple slow cycles are considered is yet to be explored.
While we do not anticipate that all the analytical predictions provided here may be carried over to the general case literally, we expect that the main principle, of achieving exponential acceleration by changing volumes of ergodic components on a non-trivial loop of parameters, is quite general.