Self-consistent time dependent two dimensional and three dimensional space charge distributions with linear force

094202-1 Finding self-consistent distributions of beam particles interacting with each other via the space charge force is one of the challenges of accelerator physics. Exactly solvable models are used for simulation benchmarks, instability threshold calculations, etc. Since such distributions have been found only in one and two dimensions (Kapchinsky-Vladimirsky distribution), it is not possible to apply them to a general three dimensional motion. This paper shows how to construct new sets of self-consistent distributions, extending even to the three dimensional case.


I. INTRODUCTION
In this paper we deal with particles moving in linear external focusing fields and interacting with each other in free space via space charge.The vacuum chamber shielding is neglected throughout the paper.We define selfconsistent distributions with linear force as distributions in which the linearity of the force is conserved not only when the beam particles are subject to external linear focusing, but also to their mutual space charge interaction.These distributions are solutions of the Vlasov equation, which states that the total time derivative of the distribution is equal to zero.The solutions of this equation are all invariants of the motion-their total derivative with respect to time is equal to zero.The following is also true -any function of the invariants of the motion is also a solution of the Vlasov equation.We call invariants (including distributions that satisfy the Vlasov equation) time dependent when they depend not only on coordinates and momenta, but also directly on time (one example is the Courant-Snyder invariant which is a quadratic form of coordinates and momenta with time dependent coefficients).
Kapchinsky and Vladimirsky obtained the first selfconsistent time dependent space charge distribution (KV distribution) in the two dimensional case [1].It has the unique property that under all linear phase space transformations the space charge force depends linearly on the coordinates.The KV distribution depends only on one invariant function of coordinates and momenta, which corresponds to a single Hamiltonian when the focusing is time independent.It was shown in Appendix A of [2] that an extension of this model to three dimensions is not possible.There is no distribution function in three dimensions that depends only on one invariant (the Hamiltonian in the time independent case) and produces linear space charge forces.
Recently, methods to obtain self-consistent time independent space charge distributions were developed for both cases of linear [3] and nonlinear [4] external focusing and space charge forces.Interesting equations for constant space charge density were analyzed in Ref. [5].The general time dependent case is beyond the reach of modern theoretical mechanics.Even single particle motion under nonlinear time dependent forces is not analytically solvable.Exceptions occur in some special cases of forces with symmetries, as, for example, in the case of round beams with special dependence of the force on time and coordinate (see [6], and references therein).Round or spherical beams yield the first examples of the new selfconsistent time dependent space charge distributions with linear forces, where the dependence is now on two invariants, one of which is the angular momentum.The next section shows how to construct such distributions.In Sec.III, it is shown how to construct distributions depending on more than one invariant in cases without special symmetries.Equations describing the evolution of the beam envelope for these distributions are presented in Sec.IV.The practical aspects of our findings are discussed at the end of the paper.

II. TWO DIMENSIONAL ROUND AND THREE DIMENSIONAL SPHERICAL SELF-CONSISTENT DISTRIBUTIONS A. First case: two dimensional round beam
To begin, consider the simple model of an azimuthally symmetric distribution function in an azimuthally symmetric channel with periodic external focusing.Any function that depends only on invariants of the motion satisfies the Vlasov equation for distribution functions (see the discussion of this matter in the Introduction).In our analysis, we use the normalized radial coordinate R r= p and its derivative R 0 with respect to the betatron phase, where r is the conventional radius and denotes the accelerator beta function for round (e.g., solenoidal) focusing.In this case, the betatron phase , which is differentially related to the longitudinal coordinate s by d ds=, serves as the new ''time'' variable.the angle corresponding to radius R. In these coordinates, the Hamiltonian of the motion is time independent.The time dependence of the distribution is related to the time dependent beta function-the equations for it are given in Sec.IV.Consider two dimensional distributions which depend on two invariants: the angular momentum M and the Hamiltonian H, where M XY 0 ÿ YX 0 R 2 0 and H R 02 =2 R 2 =2 M 2 =2R 2 .Our goal is to find all distribution functions, fH; M, that give rise to linear forces.The two dimensional interparticle force (here due to space charge only) in Cartesian coordinates is where ÿ is the region of nonzero beam density, is the linear density of the beam, r 0 is the classical radius of the beam particle, r and r are the relativistic factors, R R ÿ R R i X ÿ X i ; Y ÿ Y i , and the integration is over the variables with subscript i.Here we assume the integral of the distribution function over normalized variables is equal to unity.It is more convenient to express the force in radial coordinates via normalized variables.The transformation X i ; X 0 i ; Y i ; Y 0 i !R i ; R 0 i ; M i ; i from rectangular to radial normalized coordinates is canonical, and therefore the Jacobian is equal to 1.The force in radial coordinates has only a radial component and, after integration over the angle i , it is equal to where is the integration space (yet to be determined) over angular momentum and radial angle.In order for the force to be linear, the inner integral must be proportional to the radius, R i .We substitute for the variable R 0 i in terms of where a, b, c, and d are yet to be determined.To find c and d, we introduce the beam radius R b (all R < R b ) and fix ), the energy H i limit is determined by two stop points (velocity equal to zero).The lower stop point is the value of , and the upper stop point is the value at the beam radius R b , d , the particle never reaches R i and this integration space should be excluded; if , the particle will cross the maximal beam radius and this case should be excluded also -the distribution function is equal to zero due to the constraint R i < R b .The relation between the limits This automatically determines the limits for M i : a 0 and b R b R i .Therefore, the distribution yields a linear force if it satisfies the relation: where is an arbitrary constant.This integral equation has infinitely many solutions.Before giving a variety of solutions, we present one of the simplest distributions to show the basic principle.Equation (4) will be satisfied if the inner integral takes the form, independent of This type of equation is called Abel's integral equation (see [7]) and it has the following solution: Consider a more general case: Then, to obtain a linear force, the function g must satisfy the relation: Intuitively, this equation must have an infinite number of solutions, because the function under the integral depends on two independent variables while the integral equation itself imposes only one restriction.We look for a set of solutions in the simplest form where the integrand depends algebraically on two independent variables, R 2 i M 2 i =R 2 i and M i , multiplied by some simple functions of M i .The idea is to find variables such that, upon substitution, the expression where g 0 is an arbitrary function.It is easy to check by substituting so that the relation ( 8) is automatically satisfied.There are infinitely many choices of g 0 that satisfy but one must take care in choosing.For example, if g 0 is a linear function of its argument, the integral diverges.On the contrary, g 0 t / 1=t yields a finite value for the integral.Therefore, any function g 0 that yields a converging integral provides a good solution of Eq. ( 9).It is not clear how to obtain all possible solutions of Eq. ( 8).We leave this question for future work.
The solution for f from Eq. ( 7) is (see the solution for a similar problem in, e.g., [7])

B. Second case: three dimensional spherical distribution
We now consider the case of spherically symmetric distributions in three dimensions.In order to treat only the simple electrostatic force, we carry out the analysis in the frame of reference comoving with the beam.Under the assumption of a periodic spherically symmetric external focusing force in the comoving frame, a single set of Courant-Snyder parameters can be obtained.Then, following a procedure analogous to that above for the two dimensional case, we transform to normalized spherical coordinates R, , and .In terms of these coordinates, with the phase angle d ds= as independent variable, the single particle Hamiltonian can be written , where R 0 dR=d P R , P R 2 d=d , and P R 2 sin 2 d=d are the canonical momenta conjugate to R, , and , respectively.Hamilton's equations then show that P is a constant of the motion.For simplicity, we take this constant to be zero, thus giving H fR 02 P 2 =R 2 R 2 g=2.In this case, we now find that P is an invariant, and we consider the class of three dimensional distributions that depend on the two invariants H and P .In order to emphasize the similarity between this case and the two dimensional problem, and because P is an angular momentum, we now change notation slightly by defining M P .As before, our goal is to find all distribution functions, fH; M, that give rise to linear space charge forces.For distributions of this form, the force has only a radial component and, after integration over the angles and , and over P (the distribution is proportional to a delta function in this momentum), this force is equal to where N is the total number of particles in the bunch, is the angular momentum and radial angle integration space, and the integral of the distribution function over R i and is equal to 1.By comparing Eqs. ( 2) and ( 11), we see that the two and three dimensional results are essentially identical except for the 1=R [Eq.( 2)] versus 1=R 2 [Eq.(11)] dependence.Therefore, the analysis of the two dimensional problem carries over with slight modification, and the condition for force linearity given by Eq. ( 4) becomes To obtain the simplest solution, similar to that given by Eq. ( 6), one can again assume that the integrand does not depend on R i .Then the dependence on P is obtained automatically: The simplest solution, analogous to that given by Eq. ( 6), for this Abel's integral equation for fH; P is then Other solutions can be obtained for Eq. ( 12).For example, starting with Eq. ( 7), we now require One infinite set of solutions is similar to Eq. ( 9): with arbitrary g 0 that yields a converging integral over the angular momentum.Because we started with Eq. ( 7) here, the solution for fH i ; M i is given once again by Eq. ( 10).
In general, the distributions presented above do not produce a linear force when the focusing is not azimuthally or spherically symmetric.In the course of simulations we tried to transport one of these round distributions through a transfer line consisting of conventional quadrupoles.To our surprise we found that this particular distribution always remains elliptical with constant density.The force linearity was preserved even in transfers through nonazimuthally symmetric structures.This was the first hint that there exist special outstanding distributions with properties like those of the KV distributionany linear transformation of these distributions in the accelerator structure preserves their space charge force linearity.Moreover, contrary to the KV distribution, they can be extended to three and more dimensions, which makes them realizable in practice.These special distributions are described in the next section.

III. SELF-CONSISTENT DISTRIBUTIONS PRESERVING LINEAR FORCE UNDER ANY LINEAR TRANSFORMATION
We start with the special case of a two dimensional azimuthally symmetric distribution to show how to construct distributions that preserve the force linearity under any linear transformation.We use again the normalized coordinates X; Y and their derivatives X 0 , Y 0 with respect to the betatron phase.We assume azimuthally symmetric focusing and equal horizontal and vertical tunes.The change in beta function due to the space charge force is assumed to be included.Therefore, in normalized variables we have the betatron frequency equal to 1, and coordinates and their derivatives have the same dimension.Consider particles that have the following relations between their coordinates and angles: These relations describe a curve in four dimensional space and produce a circle when projected onto any of the following pairs of variables: X-Y, X-X 0 , Y-Y 0 , or X 0 -Y 0 .In particular, the first relation shows that this curve, projected into X-Y space, represents a circle with the radius R. Now assume we have a continuous set of such curves up to some finite beam radius R b and that the particle density, relative to R, on each curve is proportional to R. This leads to a constant density in X-Y space and, because of the linear relations in Eq. ( 17), the density is constant and equal in each of the four planes of circular projection.Since the density is constant in this case, a linear force is automatically produced.Figure 1 shows the distribution in X-Y space.As is clear from the picture, the bunch is just a rigid rotating disk.
We now prove that this distribution has an outstanding property-the linearity of its space charge force is preserved by any accelerator linear transport. 1.
Geometric argument: the four dimensional elliptical curve given by Eq. ( 17) retains its elliptical form under any linear transformation (one dimensional analog-a Courant-Snyder ellipse always remains an ellipse under all linear transformations).In particular, the number of particles in each ellipse is proportional to its area and the projection to the X-Y plane implies a uniform density and thus a linear force.

R b
Algebraic proof: Consider an arbitrary nonsingular linear transformation, U; U 0 ; V; V 0 MX; X 0 ; Y; Y 0 with inverse X; X 0 ; Y; Y 0 mU; U 0 ; V; V 0 , so that mM Mm I.The nonsingularity requirement eliminates all singular linear transformations (see footnote 1 ), a set of measure zero.The linear relations in Eq. ( 17) can then be expressed in terms of U; U 0 ; V; V 0 as follows: Unless the determinant of the coefficients of U 0 and V 0 in Eq. ( 18) is zero (another set of transformations with measure zero), these equations can be solved for U 0 and V 0 as linear equations in terms of U and Substitution of these expressions for U 0 and V 0 into X; X 0 ; Y; Y 0 mU; U 0 ; V; V 0 yields a further linear equation X; X 0 ; Y; Y 0 sU; V, where U 0 and V 0 have now been eliminated.This result allows us to express the quadratic equation in Eq. ( 17) as follows: Equation ( 19) shows that, unless the coefficient of U2 or that of V 2 is zero (another set of transformations having measure zero), the U-V space transformation of Eq. ( 17) is an ellipse.If we can show that the density remains a constant in these ellipses, then the transformed force is linear.The quadratic elliptical form of Eq. ( 19) together with the linear dependence of the density with respect to R suggests that this is true.A further illustration of this is the fact that, for linear transformations, each element of area in X-Y space transforms by a constant multiplier to an area in U-V space.We see this as follows: Because of the linear relationships in Eq. ( 17), we can write Then the element of area in U-V space is dU ^dV detTdX ^dY.Thus, each element of area in the X-Y space is mapped to an area in U-V space by a constant multiplier, detT.Therefore, unless detT 0, areas are mapped by a constant multiplier and densities remain uniform.
A similar approach works for the three dimensional case, but rather than two linear relations as in Eq. ( 17), one needs to introduce three linear relations and one quadratic relation to produce a surface in six dimensional phase space.Then this surface is populated proportional to the value of R 2 .Now we generalize to a much broader class of distributions.Since we deal with arbitrary linear focusing below, the coordinates and angles are not normalized from here.
Consider the case of n dimensional space, or 2n dimensional phase space. 2 All distributions we consider have a form, written for n dimensions where X X; X X 0 are the n dimensional vector coordinates and momenta, g is a function (to be found later), H is the Hamiltonian if the case is time dependent -in the general case it is an arbitrary quadratic positively defined function of initial conditions.H b is an upper bound of the H, and in the product of m delta functions we have m n and the matrix E 0 has rank not less than m to avoid degeneracy.In the following work we refer to this as the fn; mg case.The bounding surface H H b defines the bounding ellipsoid for this problem.Appendix A shows that constant elliptical density yields a linear force and that such a force remains linear under almost all linear transformations.
The fn; mg case here is related to n dimensional motion.The distribution function in this case depends on the quadratic form of initial coordinates and angles, which is an invariant of the motion, and on m linear combinations of initial coordinates and angles.Finding the appropriate distribution of particles, i.e., finding the function g in Eq. (20), is sometimes impossible.For example, it is shown in the Appendix of [2] that physical distributions depending only on the Hamiltonian do not exist in the three dimensional case, thus the f3; 0g case has no solution.In the fn; ng cases, the distributions are products of n delta functions of linear expressions independent in the momentum space variables and the function gH ÿ H b is just the identity, so that the integration over momentum space gives a constant coordinate space density trivially.We list now solutions for the one, two, and three dimensional cases.In all cases, we normalize f in terms of the constant real space density, .We note that the value of this density changes under transformations of the distribution.We take the function H [see (20)] for the one, two, and three dimensional cases to be H and d z are positive constants and H b is the limiting value of H.In all cases, even including fn; ng in which H does not explicitly appear in f, the condition H H b constrains the coordinates to the ellipsoidal volume so defined.Therefore, in all cases below, the coordinates and angles are constrained to lie within some boundary; the distribution function outside this boundary is equal to zero.
The f1; 1g case is very simple.It represents a rotating segment in the 2D phase space.If the line density of the segment is constant, its projection to 1D coordinate space gives constant density and linearity of the force.The distribution function is where e is an arbitrary constant.For the f1; 0g case, f is a solution of Abel's equation (see, e.g., [2], p. 13) that can be written The f2; 2g case is described at the beginning of the Section.Again, we state the obvious solution for the distribution function: The distribution for the f2; 1g case is given by f gHY 0 ÿ e 0 yx X 0 ÿ e yx X ÿ e yy Y, where g is the function to be determined and e 0 yx , e yx , and e yy are arbitrary constants.For the sake of simplicity, we present the method of finding function gH when e 0 yx 1 and e yx e yy 0. The integration over angles has to be equal to the constant density as long as the Hamiltonian condition is satisfied.
This is again Abel's equation [see Eq. ( 5)] with the solution When e 0 yx , e yx , and e yy are arbitrary constants, the same method yields the distribution function: The f2; 0g case is given by the famous KV distribution: We do not present the details of its derivation since it has been a topic of many investigations (see, e.g., [2]).The f3; 3g case follows the simple fn; ng delta-function form: where the e ij are arbitrary constants.One can see that the form of the distributions is the same for the f1; 1g, f2; 2g, and f3; 3g cases and is valid for any fn; ng case, because the integration over angles is automatically constant and, since the coordinates are restricted within the elliptical boundary, the distribution in coordinate space is a uniformly filled ellipsoid.
In the f3; 2g case, the procedure to obtain the distribution function is similar to that for the f2; 1g case.Namely, after integration over two angles, the two delta functions in the distribution f3; 2g reduce the integral equation again to Abel's form.The final result is The f3; 1g case is similar to the f2; 0g distribution.Namely, one delta function eliminates one angle variable.The integration over the two remaining angles has to yield a constant density.But this is exactly the KV case.Therefore, except for one delta-function factor, the distribution coincides with the KV, or f2; 0g case.It is x e 02 zy d 0 y e 02 zx q 2H b ÿ HZ 0 ÿ e 0 zx X 0 ÿ e 0 zy Y 0 ÿ e zx X ÿ e zy Y ÿ e zz Z; (28) where the e ij and e 0 ij are the arbitrary constants.For the f3; 0g case, as shown in the Appendix of Ref. [2], the distribution does not exist.

IV. ENVELOPE EQUATIONS FOR LINEAR FORCE DISTRIBUTIONS
Recent studies [9] show that two dimensional envelope equations give satisfactory descriptions of beam dynamics and conditions for emittance growth.Therefore it is useful to show how to obtain these equations for distributions ( 21) -(28).
To obtain the envelope equations for the round beam, we take the radial beam force from (2), assuming it is linear: where R is the normalized radial coordinate, R R b R 0 dMgR 2 M 2 =R 2 ; M=R, is the linear beam density (in units of inverse length), g is determined by ( 9), R b " p , and " is the beam boundary emittance ".Equation (10) shows how to calculate the distribution function from g.
Substitution of this force into the well-known equation for the square root w of the beta function (see, e.g., [10]), which is linearly proportional to the beam envelope, gives where primes designate differentiation with respect to the longitudinal coordinate s, Ks is the external azimuthally symmetric focusing, r 0 is the particle classical radius, and and are the conventional relativistic parameters.Using the relation a w " p between envelope a and w, the equation (30) can be rewritten for a: For the 3D spherical beam case, the force, expressed in normalized variables, is where R is the normalized radial coordinate, N is the number of particles, R R b R 0 dP gR 2 P 2 =R 2 ; P =R 2 , g is determined by Eq. ( 16), R b " p , and " is the beam boundary emittance ".In this case, Eq. ( 10) again generates the distribution function from g.The envelope equation in this case can be obtained in a manner, similar to the 2D case: Omitting the oversimplified 1D case, we will obtain equations for 2D and 3D cases [Eqs.( 23) -( 28)].Unfortunately, for the two and three dimensional cases the dynamics can no longer be described by a single envelope parameter.For example, the two dimensional ellipsoid also has an orientation angle as a parameter.For this reason we no longer deal with normalized variables.All the coordinates below are assumed to be real lengths and angles.Instead of using the complicated notations of envelopes and their orientations, we present the equations for special parameters of the beam, through which all other beam variables can be expressed.Consider the f2; 2g case with initial distribution given by Eq. ( 23).In the course of linear transformations, this distribution may lose its symmetric form.In general, the coordinates and angles of the bounding elliptical line in 4D phase space can be parametrized in the following way: x a cos b sin ; x 0 c cos d sin ; y e cos f sin ; y 0 g cos h sin ; (34) where a; b; . . .are parameters for the bounding envelope, and is a free parameter running from 0 to 2. The projection of this line onto the XY plane gives an ellipse: We have equations for x and y: x 00 K x sx ÿ e 2 3 mc 2 E x ; y 00 K y sy ÿ e 2 3 mc 2 E y ; (36) where E x and E y are the electric field components (we have ignored for simplicity the external coupling between x and y motion).To find the electric field in the free space from elliptical distribution for (35) we first need to rotate the ellipse to eliminate the coupling term, i.e., to express old coordinates x; y via new coordinates x N ; y N x 094202-7 094202-7 x N cos' y N sin', y y N cos' ÿ x N sin' such that the ellipse relation (35) in the new system does not have a x N y N term.The expression for the angle tangent can be easily found: (37) In the new system the equation for the ellipsoid reads where b 2 cos 2 ' ÿ 2fb ea cos' sin'.The electric fields for the ellipse (38) [see the Appendix A, Eq. (A3)] after the recalculation back from the rotated frame are Substituting ( 34) and ( 39) into (36) and keeping in mind that the angle is fixed for particles and should not be differentiated, we get the equations for the parameters a, b, e, and f by equating separately the factors of sin and cos : where 4r 0 = 2 3 .The equations for the remaining constants c, d, g, and h can be obtained by differentiating the coordinates and equating them to angles: The f2; 1g case with distribution of Eq. ( 24) is more complicated, but the envelope equation can be obtained in the same manner.Since the boundary in phase space of the distribution is a surface in this case, we need to parametrize the coordinates with two angles and , calculate the force via these parameters, and substitute the result into Eq.(36).We omit the details in favor of a three dimensional example.We consider the f3; 3g case, with distribution given by Eq. ( 26), which can be parametrized with two angles 0 < , < 2. We consider a diagonalized ellipsoid with no coupling between the horizontal, vertical, and longitudinal directions: x a cos cos ; y b sin cos ; z c sin : ( The potential for N uniformly distributed particles inside this ellipsoid is given in Appendix A. Substituting the potential and the coordinates into the equations of motion yields where 3Nr 0 =2 2 3 .The cases f3; 2g and f3; 1g can be obtained in the same manner, but are not shown here.

V. PRACTICAL CONSIDERATION
Some of the distributions obtained here can be realized in practice.For example, two dimensional distributions of the f2; 2g kind, shown in Fig. 1, can be constructed in the process of injection into a ring with painting.For bunches much longer than the vacuum chamber radius, the longitudinal space charge force is significantly reduced by the vacuum chamber shielding, and a two dimensional approximation normally works for self-consistent stationary (or periodic) distributions.To carry out the injection leading to distribution (23), it is necessary to satisfy certain conditions, both in the injection painting scheme and in the lattice: (i) The painted horizontal and vertical emittances should be increased linearly in time in order to yield a distribution of constant charge density; (ii) the painted horizontal and vertical phases should differ by PRST-AB 6 SELF-CONSISTENT TIME DEPENDENT TWO . . .

(2003)
094202-8 094202-8 90 to obtain an elliptical distribution in transverse space, and (iii) the vertical and horizontal tunes must be nearly equal to maintain a constant relationship between the horizontal and vertical particle phases.Figure 2 demonstrates these constraints in the creation of a f2; 2g distribution by painting into the SNS lattice [11].The calculations were carried out using the ORBIT code [12].For simplicity, a first order linear approximation was used for single particle transport and the space charge forces were calculated using the 2.5D model.Painting was done for a 1 GeV proton beam with final charge of 3 10 14 protons after 1060 turns, or about 1 ms.The painted distribution was taken to be monoenergetic and uniform in longitudinal phase space, and the rf cavities were turned off to give a coasting beam.The figure shows the particle distributions in transverse phase space at the end of injection for three cases.In all cases the beams were painted with horizontal and vertical emittances " x t " y t " f t=t f , where " f and t f are the final emittances and times, respectively, and the time dependence of the painted emittance was chosen to provide a uniform charge density.This is illustrated by the essentially linear decrease of the emittance curves in the fourth plot of Fig. 2. In the first plot (in red), the beam was painted with linear emittance increase in time and with equal horizontal and vertical tunes Q x Q y 6:23, but the bump sizes and angles were chosen to give equal horizontal and vertical betatron phases.Because of this, horizontal and vertical positions are directly correlated, and a circular distribution cannot be formed.In the second plot (in green), the beam was painted with linear emittance increase in time and with horizontal and vertical betatron phases 90 apart, so that horizontal and vertical positions are inversely correlated.Although this is a necessary condition for painting a round beam, it is not sufficient.In this case the tunes were taken to be Q x 6:28 and Q y 6:18, which causes the horizontal and vertical phases to decohere, resulting in a rectangular beam.In the third plot (in blue), the beam was painted with linear emittance increase in time, the horizontal and vertical betatron phases were painted 90 apart, and the horizontal and vertical tunes were taken to be equal, Q x Q y 6:23.Because of the 90 phase relationship at the injection foil, the equal tunes, and the constant density of the injected distribution, a round beam of type f2; 2g is painted and preserved.The 3D case is more suitable for short linac beams.Spherical distributions (14), if obtained in beam sources, require symmetric focusing.This can be achieved with solenoidal focusing by combining longitudinal kicks with solenoids to provide simultaneous equal focusing for all three dimensions-all three tunes should be always equal.It is an open question how to obtain in beam sources the distributions of Eqs. ( 26) -(28) with arbitrary linear focusing.

VI. CONCLUSION
This paper demonstrates that there exist a variety of one, two, and three dimensional self-consistent distributions with linear space charge forces.For some special distributions linearity is preserved under any linear transformation of the phase space.The practical use of these distributions is discussed.
where D X ; D 0 X are positive definite diagonal n n matri-ces and 0 is the zero matrix.The bounding surface H H b defines the bounding ellipsoid for this problem.We now demonstrate two important properties of such distributions: (1) Ellipsoids of constant space charge density in real space give linear forces; and (2) these forces remain linear under almost all linear transformations.It is not necessary that the transformation be symplectic.We note here that if the coordinate system X X; X X 0 undergoes a linear transformation, then the Hamiltonian in terms of the transformed coordinates will be symmetric and positive definite.
First, let us show that the force is linear inside any ellipsoidal distribution of constant charge density in real space.Assume in n dimensions that we have a vector (force) field and a scalar (density) that are related by r Ẽ E .To this force field there corresponds a scalar potential such that Ẽ E ÿr and r 2 ÿ.Assume also that we have an ellipsoidal volume of constant charge density given in some orthonormal coordinate system by the relation , where A is a symmetric matrix.Because A is symmetric, it can be diagonalized using orthogonal matrices, and the equation for the ellipsoid becomes b in the new orthonormal coordinates Y i .In two dimensions, the problem of constant charge density in an ellipse can be solved in terms of elliptic cylindrical coordinates, and the field inside the ellipse is linear.In three dimensions, Smythe [13] has treated the case of a perfectly conducting ellipsoidal shell, and Comunian et al. [14] have extended this work to the calculation of the field and potential of a uniform ellipsoidal charge distribution. 3Although they presented their formulation specifically for three dimensions, it is valid in the n dimensional case.The formula for the electrostatic potential is the following: where and is given by the relation P n i1 Y 2 i =c 2 i 1 for Y i outside the ellipsoid and by 0 for Y i inside the ellipsoid.Given this expression for the potential, the electric field is easily calculated to be It is easy to show that both and Ẽ E are continuous everywhere (including across the boundary of the ellipsoid), that r Ẽ E ÿr 2 inside and 0 outside the ellipsoid, and that the behavior as Ỹ Y ! 1 is correct.Furthermore, inside the ellipsoid, the electric field is clearly linear in the coordinates, thus establishing the first property.It is interesting that, in two dimensions, although the potential diverges logarithmically, the electric field converges and gives the KV result: where is the number of particles per unit length.
Let us now establish the second property, namely, the preservation of the linearity of the force under almost all linear transformations.According to the first property, this will be demonstrated if we can show that the integral of the transformed distribution function over the transformed momentum space is a constant density in the transformed real space.We want to change variables in such a way that the integral of the distribution function over the angles has the initial form.If the initial force is linear, it is linear for the new variables and vice versathis follows from the fact that all our changes of variables preserve the linearity of the force.Below we deal with the distribution fn; mg after it is integrated over m angles.
Because it has m delta-function factors, m angles disappear, and one needs to substitute linear functions (coming from arguments of delta functions) of the remaining angles and coordinates into the quadratic form H.
Starting with the X X; X X 0 coordinate system and the associated distribution function and Hamiltonian, we have the constant beam density: H< H b f X X X; X X 0 dX 0 1 ^ ^dX 0 nÿm const: Now consider the nonsingular linear transformation Then the transformed distribution function can be related to the original distribution by the conservation of the number of particles in an arbitrary phase space volume element, which gives f U Ũ U; Ũ U 0 f X X X; X X 0 detM.In the Ũ U; Ũ U 0 coordinate system, the Hamiltonian is no longer diagonal, but it is symmetric and positive definite: where S U ; S 0 U are symmetric and positive definite.The strategy is to diagonalize this expression for the Hamiltonian in a manner that retains the constancy of U during integration over momentum space.The first step is to remove the matrix C U , thus separating the momen-tum and coordinate subspaces.This can be accomplished by defining new generalized coordinates by Choosing the matrix B to satisfy S 0 U B C T U , which is possible since S 0 U is invertible, the Hamiltonian in terms of the Ṽ V; Ṽ V 0 coordinates becomes where S V S U ÿ B T S 0 U B and S 0 V S 0 U are again symmetric and positive definite.Note that integration over the momentum Ṽ V 0 space for fixed V also corresponds to fixed U, since U V. Because both S V and S 0 V are symmetric matrices, it is possible to diagonalize them separately through orthogonal transformations.Defining a final coordinate system where O; O 0 are the appropriate orthogonal matrices, the Hamiltonian can be expressed in the W W; W W 0 coordinate system as where D W O T S V O and D 0 W O 0T S 0 V O 0 are diagonal positive definite matrices.Again, because of the separation of the coordinate and momentum subspaces, integration over the momentum W W 0 space for fixed W corresponds to fixed V, and hence to fixed U. Furthermore, because the Jacobeans of the linear transformations from Ũ U, Ũ U 0 $ Ṽ V, Ṽ V 0 $ W W, and W W 0 are all unity, we find that f W W W; W W 0 f U Ũ U; Ũ U 0 f X X X; X X 0 detM detMgH b ÿ H, where the matrices F; F 0 are obtained by substituting for X X; X X 0 in terms of in the arguments of the delta functions, E ij X j E 0 ij X 0 F ij W j F 0 ij W 0 j and collecting terms.Because of the nonsingularity of the transformation relating X X; X X 0 and W W; W W 0 , the matrix F has rank m except for a set of transformations of measure zero.Let us now consider the evaluation of the density obtained by integration over momentum space in the Ũ U; Ũ U 0 coordinate system.It is obvious from the above transformations that fixing the point Ũ U and integrating over Ũ U 0 is equivalent to fixing V or W W and integrating over V 0 or W 0 , respectively.

FIG. 1 .
FIG. 1. (Color) Distribution with linear force.Arrows and circles show particle velocities and trajectories, respectively.The filled area shows the region of uniformly distributed particles.

FIG. 2 .
FIG. 2. (Color) Results of painting with linearly increasing emittance, as indicated by linearly decreasing emittance distributions, for (red distribution) equal tunes and same horizontal and vertical betatron phases, (green distribution) unequal tunes and 90 betatron phase separation, and (blue distribution) equal tunes and 90 betatron phase separation.
ÿ e xy YY 0 ÿ e yx X ÿ e yy Y; (23) where e xx , e xy , e yx , and e yy are arbitrary constants.
e xy Y ÿ e xz ZY 0 ÿ e yx X ÿ e yy Y ÿ e yz ZZ 0 ÿ e zx X ÿ e zy Y ÿ e zz Z;

3
All the expressions are given for the open space case without any shielding.