Evaluation of the transfer matrix of a plasma ramp with squared cosine shape via an approximate solution of Mathieu differential equation

The high longitudinal electric fields generated in plasma wakefields are very attractive for a new generation of high gradient plasma based accelerators. On the other hand, the strong transverse fields increase the demand for a proper matching device in order to avoid the spoiling of beam transverse quality. A solution can be provided by the use of a plasma ramp, a region at the plasma injection/extraction with smoothly increasing/decreasing plasma density. The transport of a beam inside a plasma ramp, beside its parameters, depends on the profile of the ramp itself. Establishing the transfer matrix for a plasma ramp represent a very useful tool in order to evaluate the beam evolution in the plasma. In this paper a study of a cosine squared ramp is presented. An approximate solution of the transverse equation of motion is evaluated and exploited to provide a simple transfer matrix for the plasma ramp. The transfer matrix is then employed to demonstrate that this kind of ramp has the effect to minimize the emittance growth due to betatron dephasing. The behavior of a squared cosine plasma ramp will be compared with an experimentally measured plasma ramp profile in order to validate the applicability of the transfer matrix to real cases.


Introduction
Plasma-based devices have been proven to be a valid candidate for the development of a new generation of accelerating machines [1] due to the high accelerating gradients, up to tens of GV m −1 [2][3][4], that they are able to produce.The main actual task of plasma wakefield acceleration is the preservation of the beam quality during the acceleration process both in terms of emittance and energy spread.Recent results have proven that modern plasma modules and schemes are able to accelerate bunches that can be fully characterized [5,6] and exploited in order to pilot a free electron laser [7,8], a device that notoriously requires high quality bunches only.
Due to these reasons, emittance preservation is a main topic for the design of plasma based facilities [9] and several efforts were performed in order to evaluate optimal conditions for beam transport inside plasma [10].Unfortunately, in most cases of interest, the transverse beam matching inside a plasma requires a focusing at the injection from few micrometers up to sub-micrometer scale in order to completely avoid emittance growth.A plasma ramp [11] is usually referred to as a section of plasma where the density varies smoothly.This is a partial definition that ignores the studies performed on low density constant plasma ramps [12].In this case, the increase of the plasma density is totally sharp, but, for all purposes, this can be considered a ramp.A more correct definition is that a plasma ramp is a section of the plasma profile that introduces a smooth variation of transverse beam parameters.It has been widely demonstrated that for several ramp shapes [12][13][14][15] the presence of a plasma ramp is very helpful to reduce emittance growth inside plasma.However, these studies all refer to adiabatic ramps [11][12][13], a kind of ramp where the focusing strength varies slowly compared to the oscillation wavelength of the bunch.This feature can be mathematically expressed as [12] A = 1 4η 3/2 (z) where A is defined as adiabaticity parameter and η(z) is the plasma focusing strength evolving along the direction of the bunch z.On the other side, analytical and numerical solutions have been derived for a limited number of non-adiabatic ramps with different shapes [12,14].Adiabatic ramps are the most promising in terms of performances, but less appealing in terms of applications.Tapered profiles have been discussed in order to experimentally obtain the smoothly varying shape required for adiabatic focusing [16,17].From a practical point of view, the introduction of a long adiabatic ramp at the injection and extraction of a plasma channel severely reduces the average accelerating field of the whole plasma section, frustrating the efforts for increasing the accelerating gradient that is the strong point of plasma acceleration itself.Besides, in non confined plasma, as for example from capillary discharge, non-adiabatic ramps are naturally formed [18] without further technological effort.As pointed out from Ariniello et al [12], in most cases the longitudinal shapes of the plasma ramps with an analytical solution are described by discontinuous functions or present a discontinuity of the first longitudinal derivative at the junction point between the ramp and the flat-top channel.This kind of ramp is clearly unphysical, despite most of the derived solutions can still be considered as good approximations to reality.The work proposed in this paper is focused on the study of ramps with squared cosine shape, a class of plasma ramps that never presents a fully adiabatic behavior for realistic ramp lengths.The evaluation of the transfer matrix for this ramp, that is a key byproduct of the analytical integration of the motion equations, requires the solution of Hill's differential equation when the focusing term is a sinusoidal function, an equation also known as Mathieu differential equation (MDE) [19].An approximated solution will be provided for this equation in the case of medium range non-adiabatic ramp.The solutions will be applied in order to prove the stabilizing effect of the ramp in terms of transverse emittance growth.Finally, the transfer matrix will be applied to find the matching of a measured plasma profile in order to prove the possibility to practically obtain a plasma ramp with the characteristics and the behavior of a squared cosine ramp.

Transverse matching condition of a plasma channel at the plateau
The transverse matching condition for a bunch injected into a focusing channel can be mathematically expressed in several ways, which are equivalent to each other due to the relationship that is intrinsic to the envelope equation [20] where k 2 ext (z) is the normalized focusing strength and β is the Twiss β-function.The envelope equation is derived by means of the Courant-Snyder equation in combination with the equations of motion and is valid assuming the paraxial approximation, a low energy spread and constant emittance of the bunch.The β-function expressed in equation ( 2) is derived from the equations of motions and describes the ability of a transport line to focus a particle bunch, given the required boundary conditions.We neglected the subscript x, y since this equation is valid for an arbitrary choice of transverse directions and in our work we will assume cylindrical symmetry.If the focusing force is constant along the channel, the matching condition is a trivial quasi-stationary solution of equation ( 2) It is easily verified that, assuming any two of these equations, the third one follows immediately from equation (2).The focusing term in our case of interest is given by the varying plasma focusing force.Its description is not a trivial task since plasma wave is formed as a result of the interaction of the driving pulse, either a driving bunch or a laser pulse, with the background bulk.The effects of the plasma ionization, the formation of the plasma wave and the ion motion will be neglected.
From now on, it will be assumed that the bunch transported inside plasma is totally contained in a plasma bubble in nonlinear regime.The focusing force, under these assumptions, can be represented by the ion column model [10].The plasma density n p (z, r) is generally not constant in either longitudinal or transverse direction, but for the purposes of this work, the transverse dependency of the plasma density will be neglected as well.The normalized focusing strength of the ion channel is where γ is the bunch average Lorentz factor and k p = [e 2 n p (z)/ϵ 0 m e c 2 ] 1/2 is the plasma wave number.The matching conditions for a plasma plateau with a constant density is then that for typical plasma and beam parameters is of the order of 1 µm or smaller.An input plasma ramp is a device where the density n p (z) is not constant, but rises from 0 up to the nominal value n 0 of the channel.In order to simplify any further evaluation, it will be assumed that the energy variation of the bunch inside the ramp itself is negligible.In order to verify this assumption, the longitudinal field inside the ramp must be estimated.The average accelerating gradient acting on a ramp can be evaluated from the shape of the ramp, having an equation that links the accelerating gradient to the plasma density.The most natural connection would be given by the wavebreaking limit [21] , that describes the maximum accelerating gradient that can be generated in a plasma wave without destroying its oscillating behavior.In a linear ramp arising from 0 to n 0 in a length L, the average field would be given by 2/3 of the wavebreaking limit evaluated at the end of the ramp.Thus, the assumption of negligible acceleration is certainly valid as long as where m e is the electron mass expressed in electronvolt.However, the usage of the wavebreaking limit is an extreme feature, since the plasma accelerators rarely reach fields that are even comparable.It is safe to state that the average accelerating gradient is usually at least one order of magnitude lower.The overall complexity of this aspect suggests that it is best verified with numerical simulations, afterwards.From now on, the work takes into account bunches with an energy of several hundreds of MeVs, a plasma density around 10 16 cm −3 and ramps with a length of few centimeters, that have been tested to be well below this limit.Anyway, the assumption allows from equation ( 4) to state that n p and k 2 ext are equal, except for a constant, so the focusing strength grows from 0 up to a nominal value as well.We can write k 2 ext (z) = η(z) where we assume that η(z) is a continuous and differentiable function in an arbitrary domain since the behavior of realistic plasma ramps require these assumptions.Further, from equation ( 5), one can also write η(0) = 1/β 2 0 .With a proper configuration, a plasma ramp can help to relax the matching conditions at the entrance of the plasma.We will assume that the beam is traveling with a negative momentum, namely that beam starts from z 0 > 0 at the beginning of the ramp, traveling up to z = 0 where the ramp ends and start traveling inside the channel located at z < 0 for an undefined length.The reverse transfer matrix of the channel, thus the behavior of a beam moving from z = 0 from left to right, will be computed.Based on our previous considerations, a matched beam at z = 0 will have the following Twiss vector If C(z) and S(z) are the even and odd solutions to Hill's equation inside the plasma ramp respectively the transfer matrix for the Twiss parameters [20] is where, for sake of clarity, we omitted the dependency of the functions on z.So we have a new set of equations for the Twiss functions at any point in the ramp In the assumption of a continuous and differentiable η(z), the solutions to Hill's equation are continuous, differentiable twice with continuous second derivative.

Symmetric ramps
The most simple case occurs when the ramp function η(z) is an even function.A linear differential equation of nth order is symmetric as long as the even terms, u k=2j (z) and t(z), and the odd terms, u k=2j+1 (z), have opposite parities.Recalling equation (7) one recognizes that t(z) = 0 and u 2 (z) = 1 are both even functions while u 1 = 0 is an odd function.If we assume that u 0 = η(z) is an even function too, equation ( 7) is a symmetric differential equation of the second order that admits an even, C(z), and an odd, S(z), solution.These functions can be normalized so that C(0) = 1 and S ′ (0) = 1.The boundary conditions at z = 0 can be summarized as follows We recall that the derivative of an even function is odd and vice versa.Since C and S are assumed to be differentiable, we can also write the following boundaries It is easy to verify that substituting equation ( 12) into equation ( 9), the conditions equation ( 11) are automatically respected, thus for an even function the only requirement for boundary is re-normalization of the solutions.

Squared cosine plasma ramp
The symmetric choice of the form for η(z) is the following that is a continuous, differentiable function with continuous first derivative, for which it is also possible to find a set of Twiss parameters at injection that guarantee the matching at z = 0.For this kind of ramp, the adiabaticity parameter from equation ( 1) can be evaluated as In figure 1 the value of A is shown as a function of z/L from 0 (plateau) to 1 (end of the ramp) for different values of L/β 0 .
The adiabaticity parameter presents a divergence at z = L and goes to 0 at z = 0, meaning that for any value of L this kind of ramp is non-adiabatic.For high values of L/β 0 , the values of z/L where A ≪ 1 increase, meaning that non-adiabatic behavior is concentrated in the first part of the ramp, getting adiabatic when the bunch approaches to the plateau.As a consequence, adiabatic approach is not possible and C and S shall be evaluated from Hill's equation.The Hill's equation for the ramp can be written as that, by the use of trigonometric identities and substitution of variables, can be turned into the following equation where v = π z/2L.Equation ( 16) can be identified as an MDE [22], whose general form is Assuming periodic boundary conditions, analytical solutions of MDE exist only for given couples of a and q.This does not occur in the description of a plasma ramp since the boundaries are given in equation ( 12).The general solutions to MDE are known as Mathieu even C(a, q, v) and odd S(a, q, v) functions.
As stated, the Mathieu functions can be normalized in order to respect the boundary conditions, so we can write the general solutions for the Hill's equations in the ramp as where S ′ is the Mathieu odd function derivative.The solution equation ( 18) respects all the criteria from the boundary conditions.The normalized solution can be expressed in the compact form where ξ = z/L and R = L/πβ 0 is the ratio between the length of the ramp and the betatron wavelength at the plateau.By inserting equation (19) in equation ( 8) one can evaluate the evolution of Twiss parameters inside the ramp as It is important to notice equations ( 19) and ( 20) only depend on parameter R, describing an entire class of ramps.The matching behavior of the ramps is exactly the same as long the R parameter is the same.In figure 2 the evolution of the Twiss β-function and α-function is shown, evaluated by use of equations ( 19) and (20) compared to the numerical integration of envelope equation and equation (20) and a numerical simulation performed with the hybrid kinetic-fluid code Architect [23].The considered bunch energy is E 0 = 500 MeV, the peak plasma density is n 0 = 10 16 cm −3 and the ramp length is 10 mm.As can be clearly seen, equation ( 20) is substantially equivalent to numerical integration of envelope equation and they are both in good agreement with simulation results.

Approximated transfer matrix
By setting ξ = 1 it is possible to describe the evolution of the bunch from the beginning to the end of the ramp, thus treating the whole ramp as a single focusing element.From equations ( 8) and ( 19) one can write the normalized transfer matrix for the Twiss function for an injection ramp and for an extraction ramp where c = c(R) and s = s(R).This kind of matrix is defined for the transport of Twiss parameters that are normalized with respect to β 0 , namely where β * = β/β 0 and γ * = γβ 0 .The subscript index 1 refers to the bunch parameters at the end of the ramp while the index 2 refers to the parameters in the plateau.Since Mathieu functions are transcendental, it is preferable to find an approximated form that allows to handle equations ( 21) and ( 22) in a more simple way.From figure 3, one can notice that the normalized Mathieu transfer elements have a sinusoidal behavior of period π, together with a slow evolution of the envelope.The envelope of c, s increases with R while the counterpart for c ′ , s ′ decreases.The function c has a phase delay of approximatively π/2 with respect to s, as well as c ′ with respect to s ′ .In the same way the Mathieu transfer elements present a phase delay with respect to their derivatives.Further, from Liouville's theorem we have that the transfer matrix is unimodular [20], namely Taking into account all these properties, one can assume the following approximated form for these equations:  Inserting equation ( 25) into condition equation ( 24) leads to cos(2ω) = cos 2 ω − sin 2 ω, meaning that Liouville's theorem is satisfied by equation ( 25) for any value of ω.In order to evaluate this phase delay, an useful information can be retrieved by plotting the quantity (c 2 + s 2 )(c ′2 + s ′2 ), shown in figure 4. As one can notice, this function saturates quickly at 2 for values of R > 5, whereby the consistency of the approximation holds if ω is such that In a similar fashion, it is possible to fit the values of a and b taking into account the following equation where we are considering the average between the two functions √ c 2 + s 2 and 2/(c ′ 2 + s ′ 2 ) since taking into account only one term of the sum gives lightly different fit results.The fit is shown in figure 5 and the resulting parameters are a = 1/4 and b = √ 3. The final result for the approximation is A comparison between equation ( 19) evaluated at ξ = 1 and equation ( 27) as functions of R is shown in figure 6 and the corresponding difference between the two expressions is shown in

Stabilizing effect of plasma ramp
The presence of ramp, as previously stated, helps to relax the matching conditions for the β-function.If a beam does not meet perfect matching conditions, equation ( 21) can be applied in order to retrieve β * = β(0)/β 0 and α = α(0).If these values are not exactly 1 and 0, the envelope will oscillate inside the plateau, causing an emittance growth due to betatron dephasing [24].The maximum expected emittance growth due to betatron dephasing can be evaluated through the following equation Figure 8 reports the increase of the stability valley as a function of R if we impose an emittance increase within 20% of the initial value.What one can deduce from this result is that the presence of the ramp dampens the envelope oscillations that are occurring inside the ramp itself, leading to a condition where the mismatching on the plateau is lower in a much wider range of Twiss parameters.This result is consistent and complementary with the claim performed by Dornmair et al [25] that the introduction of a plasma ramp dampens the oscillations of the centroid of the beam due to a transverse injection misalignment of a witness bunch respect to the plasma wake and limits the emittance growth.

Comparison with a realistic ramp
The assumption of a squared cosine shape was motivated by their smoothness properties that makes them good candidates for representing a realistic experimental setup.The ability to effectively control the plasma density, creating the desired shape, is of primary importance for any possible application of the present work.As previously discussed, ramps naturally arise from discharge capillary generated plasma and they can be manipulated by modifying the capillary geometry.With a proper tapering of the capillary tips, it is already possible to obtain a squared cosine plasma ramp.In figure 9 several measurements of the plasma density profile of a capillary with a step tapering are shown.In this particular configuration, the capillary is 3 cm long in total, with the gas entering in a single inlet located at the center.The diameter is 1 mm in the capillary center and becomes 1.3 mm on the edges in a 4 mm length.Plasma density is measured at various times after the discharge occurred and the peak density is in the range 1 × 10 17 cm −3 − 6 × 10 17 cm −3 .This range of densities is too high for the purposes of EuPRAXIA working point [9], but still provides an indication that it is possible to engineer the ramp shape at any density.In figure 10 it is shown the evolution of the normalized Twiss parameters for the ramps shown in figure 9.The plateau density was set at 1 × 10 16 cm −3 and the Twiss parameters at the injection are evaluated by means of equations (27).The evolution for both ramps is very consistent, meaning that small deviation from the cosine shape do not introduce major variations in the beam dynamics.The matching is not perfect, since the outcoming normalized Twiss parameters at the plateau are of the order of α ≈ 0.02 and β/β 0 = 0.85.The expected emittance growth for this kind of bunch, according to equation (28), is below 1.4%; therefore, from all practical purposes, this can be considered an optimal matching.

Conclusions
In this paper a wide-ranging study regarding the treatment of plasma ramps in linear optics approximation has been presented.In order to uniquely define the focusing strength of a plasma ramp, the ion column model has been adopted together with the assumption of negligible acceleration in the ramp, deriving a linear dependency between the focusing strength and the plasma density as a function of the longitudinal coordinate.The insertion of the ramp shape into the Hill's equation evidenced that the treatment of plasma ramps with an even functional form is more simple since there is a solid demonstration that it is always possible to find for this kind of ramp a solution that satisfy the matching condition at the plateau.The choice of a ramp with squared cosine shape has been performed since this function is continuous, derivable and even.This kind of ramp is non-adiabatic for every value of length, leading to the necessity of a solution for the differential equation without any adiabatic approximation.The resulting differential equation results to be equivalent to an MDE, that has been widely studied in literature, but none of the cases of study could be reconducted to the framework of the paper itself.The solutions given by the Mathieu functions were fully satisfying, showing an high degree of agreement not only with the numerical solution of envelope equation, but with a numerical simulation of a witness bunch injected in the wake of a driver inside plasma.This provided a validation for the use, in the proposed framework, of the ion column model with negligible deceleration for the treatment of plasma ramp.A further treatment of this case of study led to the choice of treating the ramp as a single focusing element, studying the normalized even and odd solutions of the MDE at the end of the ramp only.An approximation of this kind of solutions was found, leading to the definition of the Mathieu transfer matrix for the plasma ramp.This empirical approximation assumed the character of an analytical approximation after the evaluation of the parameters.The fact that the parameters can be expressed as functions of integer numbers and that this approximation results to be progressively more precise with the increase of R, suggests that could be possible to derive a solid and analytical demonstration of the equation ( 27).Anyhow, this kind of demonstration lies beyond the aim of this work.Beyond the great simplification of the evaluation of the matching conditions in presence of a ramp with the squared cosine shape, a successful use of the approximated equations has been performed, evaluating in a purely analytical way that one of the effects of the plasma ramp is to dampen the betatron oscillations on the plateau even in mismatched cases.Finally, the experimental feasibility of this kind of ramps have been shown, comparing the transverse evolution of a bunch both in a ramp with a squared cosine shape and a segmented profile that was measured experimentally.The convergence of the results has shown that it is possible to design and realize a plasma ramp that matches the behavior of the squared cosine plasma ramp and that it is possible as well to find an optimal matching with the only use of the proposed matching equations.In conclusion, the theoretical approach that has been performed in this paper shows great potential and a wide range of possible application in the design of plasma based accelerators.

Figure 1 .
Figure 1.Evolution of the adiabatic parameter as a function of z/L for different values of L/β 0 .

Figure 2 .
Figure 2. Comparison of the evolution of Twiss β-function and α-function inside a plasma ramp (purple) evaluated by means of numerical simulation (red), numerical integration of envelope equation (yellow) and transfer matrix of the plasma ramp (black).

Figure 3 .
Figure 3. Set of equation (19) evaluated at ξ = 1 as a function of R.

Figure 6 .
Figure 6.Comparison of the numerical solution for Mathieu functions with the approximation from equation (27).

Figure 7 .
Figure 7. Difference between exact and approximated solutions for Mathieu functions.

Figure 8 .
Figure 8. Emittance growth as a function of β inj /β 0 and α inj for different values of R. The presence of a ramp greatly relaxes the matching conditions, broadening the stability region where the emittance increase is lower than 20%.

Figure 9 .
Figure 9. Plasma density measurements in the tapered capillary.Capillary profile is schematized in the background, considering the y-axis as millimeters.On the top it is shown the normalized plasma profile measured at several times after the discharge.After 2 µs the profile stabilizes at a squared cosine like plasma profile.On the bottom it is shown the measurement at 1900 ns delay, together with a squared cosine fit and an approximation of the line made of 20 segments.

Figure 10 .
Figure 10.Twiss functions evolution of a 500 MeV bunch inside the squared cosine plasma ramp and the segmented plasma ramp shown in figure 9.The longitudinal coordinate on the x-axis is normalized respect to the ramp length.The injection conditions of this bunch was evaluated by means of equation (27) and the normalized Twiss parameters at the end of the ramp are −0.02< α < 0.02 and β/β 0 ≈ 0.85.