Physics-Informed Neural Network Method for Space Charge Effect in Particle Accelerators

The electromagnetic coupling of a charged particle beam with vacuum chambers is of great interest for beam dynamics studies in the design of a particle accelerator. A deep learning-based method is proposed as a mesh-free numerical approach for solving the field of space charges of a particle beam in a vacuum chamber. Deep neural networks based on the physical model of a relativistic particle beam with transversally nonuniform charge density moving in a vacuum chamber are constructed using this method. A partial differential equation with the Lorentz factor, transverse charge density, and boundary condition is embedded in its loss function. The proposed physics-informed neural network method is applied to round, rectangular, and elliptical vacuum chambers. This is verified in comparison with analytical solutions for coupling impedances of a round Gaussian beam and an elliptical bi-Gaussian beam. The effects of chamber geometries, charge density, beam offset, and energy on the beam coupling impedance are demonstrated.


I. INTRODUCTION
Electromagnetic coupling of a charged particle beam with vacuum chambers is of great interest for beam dynamics studies in the design of a particle accelerator [1,2]. Such coupling effects are quantified using the concept of the beam coupling impedance [2] in the frequency domain. The following two approaches have been widely used to calculate the impedances of various accelerator vacuum chambers: analytical approaches such as the mode-matching method [3,4], the image charge method [5], and numerical methods such as the finite integration technique (FIT) [6], and finite element method (FEM) [7]. In the FIT and FEM, a computational domain is decomposed into a set of volume meshes. However, when using volume meshes for the impedance computations, the standard FIT with structured grids suffers from the so-called staircasing error of curved boundaries, and the FEM with unstructured grids allows accurate boundary modeling but may require the generation of dense meshes owing to a large variation in the field in the vicinity of the space charge (SC) of a relativistic beam traversing in a vacuum chamber [7].
To address the above problems, we propose a novel meshfree approach for solving the field of SCs of a beam moving in accelerator vacuum chambers. The proposed method is based on deep learning, for example, in [8][9][10][11]. In our approach, deep neural networks based on a physical model of the beam coupling effect due to the SC field are constructed. A partial differential equation (PDE) with a charge density distribution and boundary condition (BC) is embedded in the constructed networks. By virtue of the use of the NN, our approach includes no mesh in the computational domain. This offers the advantage of modeling the SC effect compared to the FIT [6] and FEM [7]. The main purpose of this study is to demonstrate the feasibility of the proposed method for modeling the SC effect.
It should be noted that the basic idea of this method is inspired by the physics-informed neural network (PINN) [8,9], which was recently presented in the research field of machine learning. For electromagnetics, there are PINNs for timedomain Maxwell's equations [10] in microwave engineering and the PINN for inverse problems in nano-optics and metamaterials [11]. To the best of the authors' knowledge, this work is the first application of the PINN to the SC effect in particle accelerators. In this study, we develop a PINN-based numerical method for computing the SC field and the beam coupling impedance of accelerator vacuum chambers. A relativistic charged particle beam with transversally nonuniform charge density moving at constant velocity on the axis of an infinitely long vacuum chamber with elliptical cross section.

A. PROBLEM FORMULATION
Consider a relativistic charged particle beam with total charge Q moving in an infinitely long vacuum chamber, as shown in Fig.1. We assume that the beam has a rigid charge density distribution ρ and a constant longitudinal velocity v=vez=βcez, where β=v/c, c is the speed of light in vacuum, and ez is the unit vector in the direction of the beam motion (z-direction). The beam current density J has only the longitudinal component, Jz=ρv. The charge and current densities obey the following continuity equation [12]: The influence of the field on the particle distribution is neglected in the field calculations. This is not self-consistent, but an excellent approximation for relativistic beams [1,2]. Here, we consider the transverse magnetic (TM) mode generated by a relativistic beam in a perfectly conducting vacuum chamber with a constant but arbitrarily shaped cross section, because the longitudinal component of the electric field is nonzero. We study a particular harmonic component with an angular frequency ω=2πf (or wave number k=ω/v). The charge and current densities ( , ) and longitudinal component Ez can be written as [1,2,4] ( , , , ) = ⊥ ( , ) ( − ) , = where ⊥ is the transverse charge density distribution function normalized by Substituting Eqs.
where γ=(1−β 2 ) −1/2 is the Lorentz factor, also known as the relativistic factor, ε0 and μ0 denote the permittivity and the permeability of vacuum, respectively. In particular, we are interested in the case of an elliptical bi-Gaussian charge density as [13][14][15] where (σx, σy) are the half values of the Gaussian distribution in the x-and y-directions, and (xc,yc) is the center position in the transverse plane. If we set σx=σy=σr, Eq. (7) can be reduced to a round Gaussian beam [16,17]: To find Ez inside the chamber cross-section, we solve Eq. (6) with a perfectly electric conductor (PEC) BC that is, the tangential component of the electric field vanishes on the surface of the chamber cross section. Here, we call Ez the space charge field because it is the field induced by the space charge ρ in a vacuum chamber. Note that, owing to 1/γ 2 in Eq. (6), Ez vanishes in the ultra-relativistic limit γ→∞.

B. DEEP LEARING APPROACH
We introduce a novel mesh-free approach based on deep learning [8,9] for solving electromagnetic boundary-value problems in the presence of a relativistic beam current. This method utilizes the universal function approximation capabilities of deep neural networks. The key idea of the method is to include PDE (6) and PEC-BC (9) into the loss function of a neural network (NN) using automatic differentiation. It is expected that this strategy works well for a transversally smooth charge density, as shown in Eqs. (7) and (8), respectively: A schematic of the proposed algorithm is illustrated in Fig.2. Note that γ and ρ⊥ are explicitly embedded in PINN. This feature originates from the physical model described in the previous section. Our method is summarized in the following list.
1) Set up a computational domain in a chamber cross section, wall surfaces imposing the PEC-BC, source domains related to ρ⊥, physical constants and parameters such as ε0, k, γ. Note that the beam traverses inside the chamber and the SC field is zero outside the chamber. 2) Generate randomly sampled points (or grid points on a lattice) within the vacuum domain surrounded by the inner PEC walls. Note that no sampling point is generated outside the chamber. The generated sampling points are used to train a NN as input. 3) Construct a neural network with output ̂( , ; ) as a surrogate of the PDE solution ( , , ), where θ is a vector containing all weights w and bias b in the neural network to be trained, and σ denotes an activation function. 4) Define the loss function L including Eqs. (6) and (9) 5) Train the constructed neural network to find the best parameters θ by minimizing the loss L via the L-BFGS algorithm [18] as a gradient-based optimizer, until the loss is smaller than a threshold .
In the above method, we use the loss function L defined by where p denotes the sampling point. NPDE and NBC are the numbers of sampling points in the computational domain and boundary condition, respectively. LPDE is the loss function related to the PDE (6), and its minimization (LPDE→0) enforces (6) at a set of finite sampling points in the computational domain. LBC is the loss function related to the PEC-BC (9), and its minimization (LBC→0) enforces (9) at a set of finite sampling points on the boundary surface. Throughout this study, we adapted a fully connected neural network and the tanh activation function. We used three hidden layers and 20 neurons per layer. For various chamber geometries, we chose NPDE and NBC, as shown in Table I. NPDE random points are generated inside a chamber and NBC grid points are generated on the chamber wall.
It is important that a dataset is properly scaled for training, as mentioned in the literatures, e.g. [19], [20]. In this study, the scaling of the input and output to a NN is performed as preprocessing and post-processing for the proposed method. See Appendix A for more details.
It is worth mentioning that the prediction accuracy depends on the NN architecture and hyperparameters in deep learning. As stated in [8], the general trend of the PINN on accuracy is that a good prediction accuracy can be achieved as a sufficiently expressive NN architecture and sufficient numbers of sampling points are given. This trend is confirmed in Appendix B. . Physics-informed neural network for space charge effect in particle accelerator. A partial differential equation with the Lorentz factor, elliptical bi-Gaussian charge density and boundary condition is embedded in its loss function.

C. BEAM COUPLING IMPEDANCE COMPUTATION USING TRAINED NEURAL NETWORK
After training the neural network, we can simulate the SC field ̂( , ; ) for any position (x,y) in a vacuum chamber crosssection. Note that no mesh exists in this method. From the simulated field ̂( , ; ) , we can compute the beam coupling impedance (per unit length) Z|| defined as [5] where I=Qv is the total beam current. Note that the impedance at any position inside the beam domain of (7) or (8) can be obtained from (14). If one obtains the impedance averaged over the beam domain, we use the following definition [7]: where the asterisk "*" denotes taking its complex conjugate. Although Eq. (15) can be made by numerical integrations, in this paper we focus on the local impedance defined in (14), because for comparison analytical solutions for (7) and (8) are available in the local one. The evaluation of the averaged impedance (15) is out of the present paper, and it is left for future work.
In the following, we show the applications of the PINN method to accelerator vacuum chambers.

III. NUMERICAL RESULTS
To verify the proposed PINN method, we show some benchmark examples, including transversally nonuniform charge density. Our simulated results are compared with the available analytical solutions. The beam parameters and vacuum chamber cross sections used are summarized in Tables II and III, respectively. The purpose of this section is to demonstrate its feasibility for SC problems in accelerators.
It should be noted that closed-form analytical impedances of transversally nonuniform beams for any frequency are published only for very limited cases [17]. Very recently, in [15], the analytical impedance formula of an elliptical bi-Gaussian beam in free space was derived as an extension of the impedance theory of a round Gaussian beam [17]. In this study, we combine the extended impedance theory [15] with the image charge method [5] for our purpose, and calculate the analytical impedances of a round Gaussian and an elliptical bi-Gaussian beam in parallel plate and rectangular chambers. In [15], although the use of the image charge method has already been mentioned, only the free space solution is derived and discussed when we submit the present paper to this journal. Throughout this study, the number of image charges in each direction is chosen as 40. This gives a convergent result, even for the off-axis beams tested here.
The application of the proposed method to an irregularlyshaped chamber is demonstrated in Appendix C.

A. CIRCULAR CHAMBER
To assess the proposed PINN method, we apply it to a round vacuum chamber with radius R=1 cm and a round Gaussian beam with Q=1pC, γ=100 and σr=1 mm, (xc,yc)=(0 mm,0 mm). PEC-BC (9) is imposed on the chamber wall. To the best of the authors' knowledge, the closed-form exact solution for any frequency in this case has not yet been published. Therefore, we verify the numerical results with the PINN in comparison with the analytical coupling impedances known for the following two cases: a beam with uniform charge density in the PEC round chamber [5] and a beam with a transversally nonuniform (Gaussian) density (8) in free space [15][16][17]. Here, we choose a uniform beam radius rb=1.747σr. This choice allows us to approximate the space charge impedance of the round Gaussian beam by that of the round uniform beam, up to kσr/γ ≃ 0.5, as suggested in [17]. This condition reads f<2.3THz in our case. Fig. 3 shows a comparison of the SC impedances obtained for the circular chamber with the PINN simulation and the above analytical solutions [5,17]. The PINN-simulated SC impedance is in good agreement with the analytical impedance (red line) of the uniform beam in the round chamber at f<2.3 THz. At higher frequencies f>2.3THz, the PINN-simulated SC impedance is in good agreement with the analytical one (black dashed line) of the Gaussian beam in free space. Some SC field distributions (Ez) are also displayed in the insets, where the field is strong near the center of the chamber owing to the existence of the space charge and approaches zero at the wall. This illustrates that in the trained deep neural network, the charge density is taken into account and the PEC BC (9) is properly imposed.  We discuss the relationship between the SC field and impedances in Fig.3. In the left inset (0.2THz), the nonzero field spread widely inside the chamber. This indicates that the shielding effect is relatively large and can contribute to a reduction in the impedance. In fact, the circular chamber impedance becomes smaller than the free space in the lowfrequency region f<2.3 THz. In the right inset (3.1 THz), the field becomes strong only in the vicinity of the beam. This indicates that the shielding effect is small and hardly contributes to impedance. Therefore, the circular chamber impedance of the Gaussian beam is very similar to the free space impedance of the same beam in the high-frequency region f>2.3THz.

B. SQUARE CHAMBER
In order to carefully check the accuracy of the PINN method, we simulate the space charge field of an on-axis round Gaussian beam in a square vacuum chamber with a width and height of 2 cm at 0.2THz, and compare it with the corresponding analytical field. The analytical calculation can be performed by combining the free-space solution of the round Gaussian beam [15][16][17] and the image charge method [5]. Fig. 4 demonstrates comparison of the space charge field of the on-axis round Gaussian beam in the square chambers simulated at 0.2THz with the PINN method with analytical results for some cases: free space, the parallel plate and the square chambers. It is obvious that the free space SC field (red dotted line) is not zero on the chamber walls at x=±0.01m and y=±0.01m, thus it does not satisfy the PEC-BCs at all the walls. The parallel plate SC field (blue dashed line) is not zero at x=±0.01m, while it approaches zero at y=±0.01m. This is considered to be reasonable because of the effect of the parallel plate chamber walls. The analytical square chamber field is zero for all walls; thus, it satisfies the PEC-BCs. As expected, the PINN-simulated square chamber field is zero for all walls; thus, it satisfies the PEC-BCs. Its field profile is consistent with the analytical profile inside the chamber.

C. RECTANGULAR AND ELLIPTICAL CHAMBERS
Here, we demonstrate the feasibility of the PINN method for other chamber geometries. The fields of the round Gaussian beam (Q=1pC, γ=100, σr=1 mm) with an offset (xc,yc)=(rs,0) to the axis in the x-direction are simulated for rectangular and elliptical vacuum chambers at 0.2THz. rs=6 mm is chosen in this study. A discrepancy in the field distributions for the two geometries is clearly observed in Fig. 5. The simulated impedances for the different geometries are listed in Table  IV. The elliptical chamber impedance is slightly smaller than the rectangular one, but it is larger than the circular impedance. This indicates that by approximating the elliptical and rectangular chambers as a similar circular chamber, an underestimated impedance can be obtained.
The results of Sections III-A, B, and C show that the PINN method can model various vacuum chamber geometries and evaluate its effect on the impedance.

D. BEAM OFFSET
We apply the PINN to the modeling of an off-axis round Gaussian beam. The impedance of a round Gaussian beam in a square vacuum chamber with width and height of 2 cm at 0.2THz is simulated for different beam offsets. As a reference, we calculate the corresponding analytical impedances using the image charge method [5]. A comparison of the simulated and analytical impedances is presented in Table V. Excellent agreement between the two results is observed. The relative error is less than 0.3%. Here, we return to Table IV, and discuss the impact of beam offset on impedance. The impedances of the round Gaussian beam in the circular, rectangular, and elliptical chambers are simulated for different offsets with the PINN. As the offset increases, the impedances decrease in all chambers. Fig. 6 shows the SC field of an off-axis round Gaussian beam in a circular chamber. Shifting of the field peak to the chamber axis is clearly shown because of the off-axis beam.

E. TRANSVERSE Bi-GAUSSIAN CHARGE DENSITY
We apply the PINN to the modeling of the bi-Gaussian charge density with α=σy/σx≠1. The impedance of an on-axis elliptical bi-Gaussian beam with σx=1 mm in a square vacuum chamber with width and height of 2 cm at 0.2THz is simulated for α=0.5 and 2. Again, we calculated the corresponding analytical impedances using the image charge method [5].
A comparison of the simulated and analytical impedances is shown in Fig. 7. For both cases (α=0.5, and 2), the simulated impedance is in good agreement with the analytical impedance of the same square chamber. We find that the square chamber impedances of α=2 are less than those of α=0.5. A discrepancy compared to the free-space impedances is clearly observed at low frequencies. This is due to the shielding effect of the square chamber wall, as discussed in Sections III-A, III-B, and III-C. The insets demonstrate the SC field distributions simulated at 3.1THz for α=0.5 and 2. A discrepancy in the field distributions due to the different density profiles in the transverse directions is observed. VOLUME XX, 2021 9

F. BEAM ENERGY (LORENTZ FACTOR)
Finally, we show that the PINN method can consider finite γ, that is, the beam energy. We simulate the impedances of an on-axis elliptical bi-Gaussian beam with σx=1 mm and α=0.5, in a square vacuum chamber with a width and height of 2 cm at 0.2THz for a wide range of beam energies. We also calculate the corresponding analytical impedances in free space and the same square chamber using the image charge method [5]. The dependencies of the simulated and analytical impedances on the finite γ are shown in Fig. 8. The simulated square chamber impedance is in good agreement with the analytical square chamber impedance. Both impedances decrease as γ increases. This is due to the feature of the SC field, as mentioned in Section II. In contrast, for large γ, a discrepancy of the square chamber impedance to the freespace impedance can be seen. This effect is related to k/γ on the left-hand side of the PDE (6). As γ increases, k/γ decreases. This is similar to the situation when k=2πf/v, i.e. f decreases, as discussed in Section III-B. Therefore, as γ increases, the shielding of the chamber wall becomes more effective, resulting in a reduction in the impedance. A similar discussion was made in the impedance theory of a round uniform beam [5].
The results of Sections III-D, E, and F show that the PINN method can model various charge densities, beam offsets, and energies, and can also evaluate their effects on the impedance.

IV. CONCLUSION
The PINN method for the calculation of beam coupling impedances of accelerator vacuum chambers is proposed as a mesh-free approach. The method has the capability of modeling the SC field of an elliptical bi-Gaussian charge density for any beam energy and for calculating the impedances of various vacuum chamber cross sections for any frequency. The proposed method was successfully applied to circular, rectangular, and elliptical chambers. The SC fields and impedances simulated for different charge densities, different beam offsets, and different beam energies were demonstrated. Their numerical results are verified in comparison with the analytical impedances of a round Gaussian and an elliptical bi-Gaussian beam in parallel plates and rectangular chambers based on the image charge method.
Although the PEC walls are assumed, the proposed method can be extended to resistive walls. This extension will be discussed in the near future. Generalization of the proposed PINN method for transverse coupling impedances will be considered for future publication.

A. Data and Equation Scaling
In many cases, we will have a dataset including different magnitudes of input for calculating the coupling impedance. The magnitude of the SC field may vary largely within a vacuum chamber under consideration or in a frequency range of interest. In such situations, we have to deal with the dataset carefully. Directly using a raw dataset may result in a slow convergence of the gradient-based optimizer. To avoid this problem, we employ the following special scaling scheme for the input, output and PDE (6).
First, consider modeling a vacuum chamber geometry in Cartesian coordinates. The x-axis and y-axis are scaled with a typical chamber length s0 (e.g., radius, height and width) as where ⊥ is given by (7). Finally, by introducing a new parameter B and choosing E0 as we can derive a scaled PDE as follows: where = ⊥ / denotes the normalized charge density distribution. Throughout this paper, B=100 is empirically chosen. Note that γ and ρn are included even in the scaled PDE (23). From (9) and (19), the PEC-BC of the scaled SC field is given by = 0 (24) It should be mentioned that the input (x,y) is scaled to (X,Y) with s0 in (17), and the output Ez is scaled to ez with E0 in (19), and the absolute value of the right-hand side in the PDE (23) is scaled to B. This scaling is performed as pre-processing.
After the data and equation scaling described above, in the same way to step 1)−5) described in Section II-B, we construct a neural network with output ̂( , ; ) as a surrogate of the scaled PDE solution ( , , ) , define the loss function including (23) and (24), and train the constructed NN using the gradient-based optimizer. Finally, we obtain the original output (unscaled) from ̂( , ; ) = 0̂( / 0 , / 0 ; ) (25) This step is performed as post-processing.
The implementation of the above scaling scheme needs only minor modifications of (10)-(13) and the PINN (Fig.2) in the original scheme described in Section II-B. It is worth noting that the methodology of this scaling scheme can be applied to other PDEs in electromagnetics, as addressed in [10], [11]. This application is left for a future study.
Here, we shall compare Eq. (23) with Eq. (6). The righthand side (RHS) of the original PDE (6) will have a relatively large value owing to the existence of ε0≈8.854×10 −12 [F/m]. This is undesirable for using the gradient-based optimizer. In addition, the RHS depends on the wavenumber k (or the frequency f). On the other hand, when our scaling scheme (17)-(25) is used, the RHS of the scaled PDE (23) has only a constant value over a frequency range of interest. This feature is preferable for calculating the beam coupling impedance in the frequency domain and for maintaining the accuracy of trained NNs over a frequency range. In fact, our scaling scheme has been successfully applied to all benchmark examples shown in Section III. Note that the concept of our scaling scheme is not described in [10], [11].
In summary, the input, output and PDE are properly scaled in the proposed method with the above scaling scheme, which is performed as pre-processing and post-processing.

B. Impact of Hyperparameters on Prediction Accuracy
Here we confirm the general trend of the PINN on prediction accuracy in [8]. It is that a good prediction accuracy can be achieved as a sufficiently expressive NN architecture and sufficient numbers of sampling points are given. As an example, we use the same square vacuum chamber as that described in Section III-B. We define the relative error as where Zanaly denotes the analytical impedance of the on-axis round Gaussian beam in the square vacuum chamber. The analytical and numerical impedances are evaluated at the origin (x,y)=(0,0). First, we change the number of hidden layers and the number of neurons in each layer while keeping all other parameters fixed. The relative error of the impedance simulated at 0.2THz is listed in Table VI. As expected, we observe that the prediction accuracy increases, as the number of hidden layers and neurons increases.
Second, we change the number of sampling points NPDE within the chamber cross section while keeping all other parameters fixed. Fig.9 shows the influence of sampling points on the relative errors of the impedance simulated at 0.2THz. As the number of sampling points increases, the relative error decrease. The error fluctuates near 0.2%, when NPDE≥2000.
These results demonstrate the accuracy and convergence of the proposed PINN method on the impedance computation.

C. Application to L-shaped Vacuum Chamber
We demonstrate the applicability of the proposed PINN method to an irregularly shaped chamber. For this purpose, we model and simulate the SC field of an on-axis round Gaussian beam with (xc, yc) = (−5mm, −5mm) in a L-shaped vacuum chamber as shown in Fig.10. Note that it is difficult to analytically determine the SC field in this case. We used three hidden layers and 20 neurons per layer. We chose NPDE=3750 and NBC=400. NPDE random points are generated inside a chamber and NBC grids points are generated on the chamber wall. We additionally multiply the loss term LBC by a factor of 1000 to enhance the effect of the PEC-BC in training, as used in [11]. Fig.10 shows the field distribution in the L-shaped vacuum chamber simulated at 0.2THz. We find that the simulated field approaches zero on the edges and even at the corner near the origin.