Advanced fuel fusion, phase space engineering, and structure-preserving geometric algorithms

Non-thermal advanced fuel fusion trades the requirement of a large amount of recirculating tritium in the system for that of large recirculating power. Phase space engineering technologies utilizing externally injected electromagnetic fields can be applied to meet the challenge of maintaining non-thermal particle distributions at a reasonable cost. The physical processes of the phase space engineering are studied from a theoretical and algorithmic perspective. It is emphasized that the operational space of phase space engineering is limited by the underpinning symplectic dynamics of charged particles. The phase space incompressibility according to the Liouville theorem is just one of many constraints, and Gromov's non-squeezing theorem determines the minimum footprints of the charged particles on every conjugate phase space plane. In this sense and level of sophistication, the mathematical abstraction of phase space engineering is symplectic topology. To simulate the processes of phase space engineering, such as the Maxwell demon and electromagnetic energy extraction, and to accurately calculate the minimum footprints of charged particles, recently developed structure-preserving geometric algorithms can be used. The family of algorithms conserves exactly, on discretized spacetime, symplecticity and thus incompressibility, non-squeezability, and symplectic capacities. The algorithms apply to the dynamics of charged particles under the influence of external electromagnetic fields as well as the charged particle-electromagnetic field system governed by the Vlasov-Maxwell equations.


I. INTRODUCTION
Among all light-ion nuclear fusion reactions that are potentially suitable for energy production, deuterium-tritium (D-T) fusion is advantageous in terms of having the largest cross-section and lowest energy for peaked cross-section.D-T fusion reactions can be sustained in a confined, thermalized plasma consisting of deuterium and tritium ions at 10KeV temperature without significant radiation loss.This is desirable since fusion cross-sections are always much smaller than those of Coulomb scatterings between charged particles, and the fusion plasma thermalizes quickly before fusion unless powerful external electromagnetic fields are applied to maintain the non-thermal particle distributions.The cross-sections for advanced fuel fusion reactions [1][2][3][4], such as the deuterium-helium-3 (D-He3) fusion and proton-boron-11 (P-B11) fusion, peak at much higher energies.These fusion reactions are difficult to sustain in thermalized plasmas because of the significantly increased radiation loss of energy at higher temperatures [3,5].From the perspective of reducing transport due to collisions and turbulence, D-T fusion is also favorable because of the lower temperature required.
The disadvantages of D-T fusion are also prominent.First, the fusion energy is carried by 14MeV neutrons, which damage the first wall and other necessary plasma-facing hardware, such as RF antennas.Moreover, magnetic confinement devices such as tokamaks or stellarators will require massive superconducting magnets, which will need to be fully shielded from these neutrons.Much effort has been directed toward developing materials that can resist the energetic fusion neutrons.There is another difficulty.D-T fusion power plants are believed to be clean and have an unlimited supply of fuel.This is certainly true for deuterium.But for tritium, fuel self-sustainability is still an engineering challenge [6,7].Assuming initial tritium inventory is available through other venues, tritium breeding using fusion neutrons and lithium-6 is the currently envisioned solution for tritium self-sustainability.However, this requires a near-perfect tritium recycling and recovery rate, given the tritium burning fraction and tritium breeding ratio currently achievable.For an analogy to draw below, let β denote the tritium burning fraction in one particle confinement time, r the tritium breeding return (1+r is the tritium breeding ratio), and l the tritium recycling and recovery loss rate.For simplicity, we use the same loss rate for the tritium recycling from divertors/first walls and tritium recovery from the breeder blankets.Tritium self-sustainability means β + (1 − β)l < β(1 + r)(1 − l), which, when β and l are small, simplifies to l < βr. (1) For the typical values of β = 1% and r = 10% [6,7], the tritium recycling and recovery loss rate l must be less than 0.1%.This is a stringent requirement imposed by the large amount of recirculating tritium in the system associated with the low tritium burning ratio (see Fig. 1).Even if a near-perfect recycling and recovery rate is achievable, it translates to expansive capital and operational costs.
Tritium self-sustainability is still an engineering challenge.For the tritium burning fraction and breeding return currently achievable, there is a large amount of recirculating tritium in the system, and the required tritium recycling and recovery rate is in the range of 99.9%.
For commercial fusion energy production, advanced fuel fusion using D-He3 or P-B11 has been considered as a possible alternative with certain advantages [1][2][3][8][9][10][11][12][13].In addition to being tritium-free, D-He3 and P-B11 reactions are also almost aneutronic, removing the need to develop neutron-resisting materials.In exchange for these two major advantages, advanced fuel fusion needs to be operated in a non-thermal setting, for the reason discussed above, and the challenge is how to maintain the non-thermal particle distributions and the associated recirculating power at an affordable cost [1,2,8,14,15].In this sense, advanced fuel fusion trades the requirement of a large amount recirculating tritium for that of large recirculating power.It is reasonable to expect that recirculating power is easier to generate and maintain in comparison.To meet the challenge, innovative phase space engineering using electromagnetic fields can be applied [16][17][18][19][20][21].For charged particles in ionized plasmas, injecting high-frequency electromagnetic fields seems to be the only way to manipulate phase space distributions faster than the collisional thermalization, and such methods have been extensively investigated for achieving hot ion modes via α-channeling in D-T fusion [22][23][24][25][26].In the present study, we analyze the physics of phase space engineering from a theoretical and algorithmic perspective.In particular, we emphasize the phase space constraints imposed by the Hamiltonian nature of the dynamics and discuss structure-preserving geometric algorithms that are indispensable in the study and design of phase space engineering techniques.In Sec.II, we discuss the need for phase space engineering for non-thermal advanced fuel fusion.Symplectic dynamics of charged particles, the underpinning physics of phase space engineering, and the associated phase space constraints are analyzed in Sec.III.
Structure-preserving geometric algorithms suitable for phase space engineering are presented in Sec.IV.
For non-thermal advanced fuel fusion, large recirculating power exists in the system to maintain the non-thermal distributions of particles, which does not necessarily entail significant energy loss if enough attention is paid to preserve the power flow.An example is the energy-recovering accelerator technology [38].If the success of D-T fusion energy is only possible with lossless tritium recirculation in the system, non-thermal advanced fuel fusion depends on phase space engineering that can sustain efficient recirculating power.
Let's demonstrate conceptually how an RF field-driven Maxwell demon [22,[39][40][41][42][43][44][45][46] can be utilized to maintain the directed energy of a fusing ion beam.In the co-moving frame of the beam, a Maxwell demon is implemented as a one-way door that allows ions with forward velocity to pass through but reflects ions with backward velocity, see Fig. 2. The one-way door guarded by the demon can also be called a one-way wall.It converts the thermal energy into the directed kinetic energy of the beam.In this process, the entropy of the ion beam is decreased, which is associated with certain energy costs and increased entropy of the entire system including the external hardware generating the RF field.The technical challenge is to minimize the cost within the operational space allowed by physics.Not all phase space manipulations are possible if the Maxwell demon was driven by electromagnetic fields.For example, the phase space volume of the ions needs to be conserved according to the Liouville theorem, which is known as phase-space incompressibility (see Sec. III).When the backward-moving ions are pushed forward by the demon, they have to occupy the same volume in the left part of the phase space as shown in Fig. 3(a).If the operation is confined in the (x, v x ) plane, incompressibility rules out the possibility of forming a cold, focused, energetic beam as shown in Fig. 3(b), which would be ideal for fusion reactions.
One of the advantages of aneutronic fusion is that the fusion energy is released in the form of kinetic energy carried by charged particles, and the fusion energy can be extracted directly and efficiently as electromagnetic energy.However, because the kinetic energy of the fusing ions is much smaller than the fusion energy released, the released energy is largely thermalized and occupies a certain phase space volume that needs to be conserved during the direct energy extraction via electromagnetic channels.For a given distribution of fusion products, e.g., α-particles in the P-B11 fusion, what is the maximum fusion energy that can be extracted electromagnetically [47][48][49][50]?Or equivalently, what is the ground state of the system via electromagnetic channels?Because of the phase-space incompressibility, the ground state energy is not zero.
Gardner studied the ground state of a phase space distribution accessible by all possible volume-preserving maps of phase space [51].Let f 0 (x, v) denote the initial distribution and ϕ(x) an external potential.Gardner showed that the ground state accessible via volumepreserving maps is f 1 (ϕ(x)+mv 2 /2), where f 1 is a function of one scalar argument satisfying: (1) f is monotonically decreasing and (2) for any ξ > 0, the phase space volume of the region where f 1 (ϕ(x) + mv 2 /2) > ξ is equal to the phase space volume of the region where The function f 1 is uniquely determined.Let's call this ground state accessible via volume-preserving maps Gardner ground state.Intuitively, the Gardner ground state is reached by "doing a Lebesgue integral", which can be explained as follows.Suppose a coin purse dropped and the coins inside scattered on the floor.There are two ways to pick up the coins and count your assets.You can pick up all the coins on the first square foot of the floor, then all the coins on the second square foot, and so on.This is the Riemann integral.Alternatively, you can pick up all the quarters first, then all the dimes, then all the nickels, and finally all the pennies.This is the Lebesgue integral, and in Gardner's terminology, ξ = 25¢, 10¢, 5¢, and 1¢.Gardner's construction was later interpreted as Gardner re-stacking algorithm [52,53].The algorithm moves all the quarters to the center of the floor without changing the size of their footprints, then all the dimes to the region next to the quarters, and so on.See Fig. 4.

Lebesgue integral
Gardner re-stacking It turns out that incompressibility, or phase space volume conservation, is just one of many constraints in phase space.Another important constraint is known as non-squeezability, which is a consequence of the symplectic nature of the Hamiltonian dynamics of charged particles under the influence of electromagnetic fields.To design effective phase space engi-neering schemes, it is crucial to understand these dynamic constraints and adopt structurepreserving simulation algorithms that preserve these constraints.As we will see, the Gardner ground state may not be accessible by the Hamiltonian dynamics of charged particles under the influence of externally injected or self-consistently generated electromagnetic fields.

III. SYMPLECTIC DYNAMICS OF CHARGED PARTICLES -THE PHYSICS OF PHASE SPACE ENGINEERING
In this section, we briefly describe the constraints of phase space dynamics as imposed by the law of physics and the relevance to phase space engineering.The dynamics of charged particles under the influence of an external or self-consistently generated electromagnetic field is governed by Hamilton's equation, where H = H(x i , p i , t) is the Hamiltonian function.For a general electromagnetic field on spacetime specified by 4-potential (ϕ(x, t), A(x, t)), Here, we have taken the non-relativistic limit of H for simplified presentation.The dynamics according to Eqs. ( 2)-( 4) is given by a solution map It is easy to verify that the solution map φ t is symplectic, which, by definition, means that the Jacobian matrix of φ t , is symplectic.
In general, a 2n × 2n matrix S is called symplectic if S T JS = J, where which defines an almost complex structure on R 2n , i.e., J : R 2n → R 2n and J 2 = −1.
More fundamentally, the symplecticity of the solution map φ t is a geometric property.The geometric form of the canonical Hamilton's equation in R 2n is where is the canonical symplectic form.In Eq. ( 8), "i" denotes inner product between the 2-form Ω and the vector field ( ẋ, ṗ), and in Eq. ( 9) "d" denotes exterior derivative.
The solution map φ t being symplectic means geometrically that Ω is an invariant of the dynamics, i.e., where φ * t denotes the pullback by φ t .An immediate consequence of Eq. ( 10) is the volumeform conservation [54], where is the volume-form defined by the symplectic 2-form Ω in the 2n dimensional phase space.
Equation ( 11) is the well-known Liouville theorem in the geometric form.Note that the dimension of Ω, as a 2-form on R 2n , is n(2n − 1), and the dimension of Λ 2n is 1.Volume conservation (11), or incompressibility is necessary but not sufficient for symplecticity (10).
One may be content with the volume conservation as stated in the classical Liouville theorem and question the importance of symplecticity.The fact is that from the viewpoint of physics, phase space volume can only be meaningfully defined through the symplectic 2-form Ω. Recall that length, area, and volume in spacetime can only be defined by the metric tensor that is determined by Einstein's equation.However, the spacetime metric tensor does not define a volume-form in phase space, the symplectic 2-form Ω does [55,56].
At the fundamental level, phase space volume conservation is one of many properties of symplectic maps.
In addition to incompressibility (volume-preserving), the non-squeezability of symplectic maps proved by Gromov in 1985 is the most celebrated one [57][58][59][60][61][62][63][64][65].Gromov's non-squeezing theorem states: If there exists a symplectic map φ in R 2n sending the ball B 2n (r) into some cylinder Z 2n j (R), then r ≤ R. (See Fig. 5) Here, This is a strong constraint on symplectic maps considering the fact that B 2n (r = R) ⊂ Z 2n j (R).The ball B 2n (R) is already in the cylinder Z 2n j (R), there is no need to squeeze it.But if we increase the size of the ball by an infinitesimal amount ε, then we can't push the ball B 2n (R + ε) into the cylinder Z 2n j (R) via any symplectic map.In R 2n , the Euclidean metric ⟨ , ⟩, the canonical symplectic form Ω, and the almost complex structure J are compatible, i.e., ⟨a, b⟩ = Ω(a, −Jb).
The distance in Eqs. ( 13) and ( 14) is measured using this metric.In addition, all maps in the present study are assumed to be diffeomorphisms (smooth) unless explicitly stated otherwise.
It is easily shown [64] that the non-squeezing theorem is equivalent to the fact that under a symplectic map φ, the projection or footprint of B 2n (r) on any (x j , p j ) plane is always larger than πr 2 , i.e., Area of Non-squeezability stipulates that the phase space footprint on each (x j , p j ) plane (1 is non-diminishing relative to an initial set that is a ball B 2n (r) in phase space.
Non-squeezability puts strong constraints on phase space engineering.Take the example of the Maxwell demon discussed in Sec.II.Starting from a phase space ball in R 6 , after sending the system through a Maxwell demon in the (x, p x ) plane, we may wish to reduce the footprint of the particle system in the (x, p x ) plane to reduce the energy cost.This would be achievable by squeezing a certain amount of phase space volume across the (y, p x ) or (z, p z ) plane if phase space volume conservation were the only constraint.But, the nonsqueezing theorem prohibits this type of operation.For a given initial distribution f 0 , the corresponding Gardner ground state may not be accessible by symplectic maps.
If the initial set V of particles in phase space is not a ball, what will be the minimum footprint in a conjugate plane (x j , p j ) (1 ≤ j ≤ n) under a symplectic map φ?In this case, we can search for the largest ball that can be squeezed inside V by a symplectic map, where Sym is the set of all symplectic diffeomorphisms (symplectomorphisms) in R 2n .Obviously, C o (V ) is a symplectic invariant and can be used as a measure of the minimum footprint.Similarly, we can also search for the smallest cylinder that V can be squeezed into by a symplectic map, and C 1 (V ) is another symplectic invariant that can be used to measure the minimum footprint.Because there exist symplectic maps to switch (x j , p j ) and (x i , p i ) for any i and j, the definition of C 1 (V ) is independent of j.According to the non-squeezing theorem, Thus, for the purpose of phase space engineering, C 0 (V ) and C 1 (V ) can be used to obtain a good estimate of the minimum footprint on any phase space conjugate plane (x j , p j ).Here, C 0 (V ) and C 1 (V ) are examples of a family of symplectic invariants called symplectic capacities [59,61,[65][66][67][68] that characterize symplectic maps in a similar manner that phase space volume characterizes volume-preserving maps.In fact, there are infinitely many symplectic capacities.It turns out that C 0 (V ) is the smallest and C 1 (V ) the largest.
We now display a proof of the non-squeezing theorem for the case of linear symplectic maps to illustrate the physics related to minimum footprints and symplectic capacities.The proof is a slightly modified version of a known proof given in Refs.[65,68].(The proof of the general non-squeezing theorem is difficult but not inaccessible to physicists.) Linear symplectic maps are the solutions of linear Hamiltonian systems, which are found in many branches of physics.In particular, the dynamics of high-intensity charged particle beams in accelerators, storage rings, and transport/focusing channels is governed by timedependent linear Hamiltonian systems [37,[69][70][71][72][73][74][75][76][77].For non-thermal fusion technologies, the dynamics of high-intensity particle beams is expected to play an important role.The solution of a linear Hamiltonian system is a linear map specified by a 2n × 2n symplectic matrix S, Here, z = (x 1 , x 2 , ..., x n , p 1 , p 2 , ..., p n ) .
For linear symplectic maps, Gromov's non-squeezing theorem reads: If there exist a symplectic matrix S and an integer j ≤ n such that S (B 2n (r)) ⊂ Z 2n j (R), then R ≥ r.Without losing generality, we prove the case of r = 1 and j = 1.The strategy of proof is to show that for some z 0 on the boundary of B 2n (1), the absolute value of the x 1 or p 1 coordinate of Sz 0 is no less than 1.However, by definition, B 2n (1) is open, and we need to select a sequence of points in B 2n (1) to approach z 0 .Using the metric defined in Eq. ( 15 where the second equal sign is because S T is symplectic, and inequality is the Cauchy-Schwartz inequality.Therefore, either |ā| ≥ 1 or b ≥ 1. Assume |ā| ≥ 1 and let The 1), and its p 1 coordinate will be no less than 1 − ε.Thus, the radius of Z 2n j (R) should be no less than 1 − ε, since Sz ε ∈ Z 2n j (R).Let ε → 0 and we have R ≥ 1.This completes the proof.
In the proof, we find on the boundary of B 2n (1) a point z 0 , whose footprint under any map S on the (x 1 , p 1 ) plane sticks out the unit circle.But neither z 0 nor Sz 0 is necessarily on the (x 1 , p 1 ) plane.z can be anywhere on the surface of B 2n (1), depending on S. The fact that the footprint of S (B 2n (1)) sticks out the unit circle on the (x 1 , p 1 ) plane for any S is attributed to the information on the entire B 2n (1), and this situation is exactly the same for every (x j , p j ) plane.In terms of physics and phase space engineering, the minimum footprint of a group of particles on any (x j , p j ) plane is the same as measured by a symplectic capacity, which is a scalar assigned to the entire group of particles.But we have choices as to which symplectic capacity to use.
For the last topic in this section, let's go back to the canonical Hamilton's equation for charged particles specified by Eqs. ( 2)-( 4).The kinetic momentum p is the sum of kinetic momentum m ẋ and unit-normalized vector potential A(x, t) at the particle's location.If we choose to use x and v = ẋ as the phase space coordinates, Hamilton's equation still assumes the form of Eq. ( 2) except that the canonical symplectic form Ω will be replaced by a non-canonical symplectic form ω, i.e., It is straightforward to verify that Eqs. ( 23)-( 25) are equivalent to the familiar equations of motion, If the phase space engineering needs to consider the electromagnetic field self-consistently generated by the fusion plasma itself, the charged particle-electromagnetic field system is an infinite dimensional non-canonical Hamiltonian system specified by the Morrison-Marsden-Weinstein (MMW) bracket [78][79][80][81][82][83], where F and G are functionals of the function space {(E, B, f The bracket {h, g} xv inside the first term on the right hand side is the canonical Poisson bracket for functions h and g of (x, v), Hamilton's equation for the dynamics of a functional F is where Hamiltonian functional is The equations of motion for the fields are obtained by expressing them as functionals indexed by x and v, The Hamiton's equation (28) for F (x, v), E (x) , and B (x) recovers the standard Vlasov-

IV. STRUCTURE-PRESERVING GEOMETRIC ALGORITHMS FOR PHASE SPACE ENGINEERING
While phase space incompressibility has been well-known for centuries, phase space nonsqueezability has only been established for 40 years.The associated dynamic invariants, symplectic capacities, have not been widely recognized in most physics applications.For example, Gardner's ground state analysis for plasmas only considered the phase space incompressibility.When the non-squeezability is included in the analysis, the ground state is expected to be very different, and the ground state energy will increase.From the perspective of fusion energy production, this means less energy can be extracted via conservative electromagnetic channels.In this sense and at this level of sophistication, phase space engineering for non-thermal advanced fuel fusion is abstractly symplectic topology.Admittedly, it is difficult to directly incorporate phase space non-squeezability and other constraints of symplectic topology in analytical studies or numerical simulations of phase space engineering.For instance, the minimum footprint estimates according to Eqs. ( 17) and ( 18) are important design parameters, but it is not clear how they can be accurately calculated.
In this section, we explain in detail the mechanisms by which these algorithms preserve exactly the symplecticity and thus incompressibility, non-squeezability, and symplectic capacities of the Hamiltonian systems.To be concise, we will focus only on the high-order, explicitly solvable, non-canonical symplectic algorithms for the charged particle dynamics governed by Eqs. ( 23)- (25).
Let φ ∆t denote the exact solution map of Eq. ( 23) in a time interval ∆t, Denote by ν ∆t a numerical solution map of Eq. ( 23) in a time interval ∆t, ν ∆t , called one-step map of the algorithm, is an approximation to φ ∆t .We would like ν ∆t to be exactly symplectic as well, i.e., Even though ν ∆t is not exactly the same as φ ∆t in general, we prefer ν ∆t to come from the same family of symplectomorphisms because all important properties of Hamiltonian dynamics, such as incompressibility and non-squeezability, are associated with the symplectomorphisms as a group, instead of with a specific member of the group.An algorithm is called symplectic if its one-step map is symplectic.
For a general non-canonical symplectic system, condition ( 39) is difficult to satisfy unless ν ∆t is the exact solution φ ∆t .In fact, no generic symplectic algorithm exists for general non-canonical systems.But for a specific given non-canonical system, bespoke symplectic algorithms might exist.Equations ( 23)- (25) for charged particle dynamics are such a system.
The bespoke symplectic algorithms are constructed by a Hamiltonian splitting technique.
The symplectic condition (39) is satisfied if ν ∆t is an exact solution.Of course, we can't write down the exact solution corresponding to the ω defined by Eq. ( 24) and H defined by Eq. ( 25).However, observe that for the ω defined by Eq. ( 24), the symplectic condition ( 39) is satisfied exactly by the exact solution φ ∆t generated by any Hamiltonian function.For many simple Hamiltonian functions, the exact solution maps can be written down explicitly, and the compositions of these exact solution maps will also preserve the ω defined by Eq. ( 24).
If we can split the Hamiltonian H defined by Eq. ( 25) into different parts such that the exact solution for each sub-Hamiltonian system can be written down explicitly, then these exact solutions of the sub-systems can be composed to approximate the exact solution of H to arbitrary high orders, and such approximate solutions automatically satisfy the symplectic condition (39).Such a splitting was found.Specifically, we split H into four sub-Hamiltonians [85-87, 90, 91] , Here, (x, y, z) are the Cartesian coordinates of the configuration space, and (v x , v y , v z ) are the corresponding velocity coordinates.
For the sub-system specified by H ϕ , Eq. ( 23) is which, in terms of ẋ and v, is where Note that in Eq. ( 46), the time derivative of the vector potential enters through the noncanonical symplectic form ω defined in Eq. ( 24).The exact solution map of Eq. ( 43) can be explicitly written down, φ ϕ,∆t : For the sub-system specified by H x , Eq. ( 23) is In Eq. ( 50), the magnetic field B(x, t) = ∇ × A(x, t) also enters through the non-canonical symplectic form ω. Equations ( 49) and ( 50) are more complicated than those of H ϕ .But the exact solution map of H x can also be written down explicitly [85-87, 90, 91], The exact solution maps φ y,∆t and φ z,∆t for the sub-systems H y and H z can be explicitly written down in a similar manner.
It is crucial to recognize that the exact solution maps φ x,∆t φ y,∆t , and φ z,∆t are written down in the Cartesian coordinates.It has been shown that for an orthogonal, uniform curvilinear coordinate system (q 1 , q 2 , q 3 ) in R 3 , the exact solution for sub-system H q i can also be written down explicitly [91].An orthogonal coordinate system is called uniform if where h i (i = 1, 2, 3) are the Lamé coefficients.If a curvilinear coordinate system is not orthogonal or non-uniformly orthogonal, the solution maps for the subsystem defined by H q i can't be written down explicitly in general [91].
Arbitrary high-order symplectic algorithms for the original Hamiltonian system ( 23)-( 25) can be constructed from different compositions of φ ϕ,∆t , φ x,∆t φ y,∆t , and φ z,∆t .Two of the 1st-order algorithms for H are ν * (1) disadvantages are also prominent.The 14MeV fusion neutron demands new technologies for energy conversion and hardware protection, and the fuel self-sustainability of tritium imposes a stringent constraint on the recycling and recovery rate.Advanced fuel fusion using D-He3 or P-B11 might be a viable alternative for commercial fusion energy production.However, to be practical, advanced fuel fusion needs to be in non-thermal regimes to avoid excessive radiation loss of energy because the fusion cross-sections for advanced fuel fusion peak at much higher energies.This leads to large recirculating power in the system.
Advanced fuel fusion trades the requirement of a large of amount recirculating tritium for that of large recirculating power.To meet the challenge of maintaining non-thermal particle distributions at a reasonable cost, phase space engineering technologies utilizing externally injected electromagnetic fields can be applied.For charged particles in ionized plasmas, injecting high-frequency electromagnetic fields seems to be the only way to manipulate phase space distributions faster than the collisional thermalization, and such methods have been extensively investigated for achieving hot ion modes via α-channeling in D-T fusion.
In the present study, we have studied the physical process of phase space engineering, such as the Maxwell demon and electromagnetic energy extraction, from a theoretical and algorithmic perspective.We demonstrated that the operational space of phase space engineering is limited by the underpinning symplectic dynamics of charged particles.Volume conservation, or incompressibility, according to the well-known Liouville theorem is just one of many phase space constraints.Gromov's non-squeezing theorem determines the minimum footprint of the charged particles on every conjugate phase space plane.Other constraints, such as symplectic capacities, need to be included in the phase space engineering designs as well.In this sense and level of sophistication, phase space engineering is abstractly symplectic topology.To calculate the minimum footprint of charged particles for phase space engineering and to accurately simulate the processes of phase space engineering, recently developed structure-preserving geometric algorithms can be used.The family of algorithms conserves exactly, on discretized spacetime, symplecticity and thus incompressibility, non-squeezability, and symplectic capacities.The algorithms are implicit but explicitly solvable and apply to finite dimensional non-canonical symplectic dynamics of charged particles under the influence of external electromagnetic fields, as well as to the infinite dimensional non-canonical symplectic charged particle-electromagnetic field system described by the Vlasov-Maxwell equations in the geometric form.

Figure 2 .
Figure 2. A Maxwell demon guarding a one-way door.The blue balls inside (representing backwardmoving particles) are bounced back at the door, and the red balls outside (representing forwardmoving particles) are allowed in.

Figure 3 .
Figure 3.Not all phase space manipulations by a Maxwell demon driven by electromagnetic fields are possible.The phase space volume of the ions needs to be conserved according to the Liouville theorem.If the operation is confined in the (x, v x ) plane, incompressibility rules out the possibility of forming a cold, focused, energetic beam as shown in (b), which would be ideal for fusion reactions.

Figure 4 .
Figure 4. Gardner ground state is reached by "doing a Lebesgue integral".

Figure 5 .
Figure 5. Gromov's non-squeezing theorem.No symplectic map φ can squeeze the ball B 2n (r) into the cylinder Z 2n j (R) when r > R.But for any r ≤ R, the ball B 2n (r) is already inside the cylinder Z 2n j (R).The ball B 2n (r) in phase space is "rigid".
), we select the orthonormal bases (a, b) for the (x 1 , p 1 ) coordinates.The images of a and b under S T are ā = S T a and b = S T b.Consider 1