Electromagnetic Particle Algorithm for Beam–Wave Interaction in Traveling Wave Tube of Symmetry

Abstract

In many fields, such as space astrophysics, plasma and vacuum electronics, there are many nonlinear strong coupling physical problems. In order to solve the problem of electron beam–wave interaction in cylindrical Traveling wave tube (TWT) with symmetrical structure, a multi particle simulation algorithm for beam circuit is designed. The algorithm allows aperiodic time input, nonuniform linearity and large space diagnosis for different particles. In this algorithm, the field of coupled slow-wave transmission line is simulated by finite difference method. Assuming that there is strong coupling between the beam and the circuit, the space center equation of transmission along the line is obtained, and the space charge field is modeled considering the space charge effect, which can easily be ignored. The Particle In Cell (PIC) method of frog leaping step scheme is adopted to evaluate the electric field of each particle center, determine the circuit and space charge field, and design the termination part to compensate for the loss in order to avoid self-excited agitation. Finally, a simple numerical simulation is carried out for the electromagnetic problem and the accuracy of the code is checked against the electromagnetic simulator CHPIC. Therefore, the algorithm can be used to solve the problem of beam–wave interactions in a fixed structure (cylindrical) in multiple fields and can accurately record the data in the process.

1. Introduction

As the best combination of bandwidth and gain in microwave tube, the traveling wave tube (TWT) is still ahead of solid-state amplifier in bandwidth and efficiency. TWT has the advantages of high power, high efficiency, high gain, wide bandwidth, large dynamic range and low noise [1]. It is also widely used in radar, electronic countermeasures and satellite communication. In the era of the rapid development of computer technology, computer numerical simulation technology is also developing rapidly [2,3]. In terms of the problem of scientific research and practical application being difficult to link, the computer numerical simulation method can conduct numerical modeling of the studied system and then simulate from a given initial state and obtain the whole simulation process and meaningful results, which is a good solution [4,5,6]. The numerical simulation results can also be compared with the theoretical calculation results and experimental results to study the essence of the problem.

Particle simulation method is an excellent computer numerical simulation method, which can be used to study the interaction process of charged particles with electric field and magnetic field [7,8,9]. It plays an irreplaceable role in the application of electronic signal amplifier, space astrophysics, plasma and vacuum electronics and has important military application prospects [10,11]. The particle simulation method can fully consider the various factors of particle motion, truly reflect the particle motion, and obtain better results than the experiment [12,13,14,15]. At the same time, in the process of simulation calculation, all motion data of each simulated particle will be recorded in the computer, and then detailed data information will be obtained for observation and analysis.

The interaction between the magnetic field generated by the rapid movement of the electron beam in the aircraft and the electron beam body can accumulate enough energy to release, thus causing damage to the aircraft [16,17,18]. For example, the electrons moving to the space around the aircraft gather on the surface of the aircraft, which leads to surface electrification [19,20]. This charge and related magnetic field will lead to surface arc discharge, electronic blowout, enhanced pollution and the movement of the electrical ground of the aircraft. The continuous flow of electric charges entering and exiting the aircraft surface constitutes a current, which may lead to electric leakage inside or outside the aircraft and become a potential safety hazard [21,22,23,24,25]. Therefore, it is of great significance to study the electromagnetic radiation generated during the launch and movement of the aircraft, analyze the electromagnetic characteristics of the aircraft surface potential change in time and frequency domain, and then effectively protect and shield its impact [26,27,28]. Because of the background of high current electron beam and electromagnetic static field, many physical problems are nonlinear, and there are strong coupling problems [29,30]. It is difficult to obtain accurate results only by traditional analytical methods and linear analysis. The physical process must be simulated by the electromagnetic particle algorithm. Therefore, it is very necessary to design the electromagnetic particle algorithm for beam–wave interactions.

In this algorithm, the finite difference method is used to simulate the field of the coupled slow-wave transmission line. The coupling between the beam and the transmission line is based on the finite difference equation proposed by Brillouin (1949) [31,32]. Based on the difference principle, it directly starts from Maxwell’s curl equation, which summarizes the general law of electromagnetic field, converts it into difference equations and samples the data of the continuous electromagnetic field at a certain volume and a certain period of time. Therefore, it is the most basic, essential and complete numerical simulation of electromagnetic field problems and has the widest applicability. It is also based on the characteristics of the FDTD algorithm, which it is mostly used as the time domain algorithm of electromagnetic field in electromagnetic PIC particle simulation. The electron beam emitted from the cathode of TWT is accelerated by the beam–wave interaction in the high-frequency system and emits high power electromagnetic waves. The study of beam–wave interaction in TWT is basically to solve Maxwell equation with a source (current source formed by electron beam) and the dynamic equation of electron motion.

2. Coupling between Beam and Circuit

It is assumed that the electron beam and the circuit will have a strong coupling effect. All beam charges in the traveling wave tube in Figure 1 generate charges on the equivalent circuit.

$I_{N-1/2}$ is the current between nodes $N-1$ and $N$; $V_N$ is the voltage at node $N$; $Q_{\mathrm{N}}$ is the charge on the capacitor at node $N$; $q_N$ is the line-deposited; $q_N^{\prime }$ is the beam-deposited; and $Q_{\mathrm{N}}=q_N+q_N^{\prime }$.

Using Brillouin notation, the following spatial center equation transmitted along the route can be obtained:

$$\begin{array}{cc}{I}_{N-1/2}-I_{N+1/2}=\frac{\partial q_N}{\partial t},& L_0(\frac{\partial I_{N-1/2}}{\partial t})=V_{N-1}-V_N\\ V_N=\frac{Q_N}{C_0},& Q_N=q_N+q_N^{\prime }=q_N+{\rho }_N{d}_0\end{array}$$

(1)

However, in our proposed algorithm, the beam charge $q_N^{\prime }$ is only a small fraction of the beam space $G$ charge in the guide section and depends on the wavenumber of the space charge wave. The capacitor charge is shown in Equation (2).

$$Q_N=q_N+\kappa (k)q_N^{\prime }$$

(2)

Beam–circuit coupling coefficient $\kappa (k)\le 1$ reflects the nature of the space charge coupling constant of a finite-diameter traveling wave tube, where beam motion generates displacement currents, which in turn generate conduction currents in slow waveguides. The functional dependence of $\kappa (k)$ on $k$ is determined by the relative diameter of the electron beam and slow waveguide, as well as the transverse field distribution for the different modes.

3. Space Charge Effect

Some early analysis models ignored the effect of space charge, which may be acceptable for wavelet amplitude, but it cannot be ignored for large amplitude signals close to saturation. The early space charge model follows the wavelength of a disk or a ring and usually only has an approximate space charge field. The total force acting on the electron is the sum of the forces generated by the circuit field and the space charge field. However, in our proposed algorithm, PIC technology (i.e., space grid) is used to model the space charge field, as shown in Figure 2.

In order to realize the difference calculation of spatial coordinates and considering the orthogonal and cross-linked relationship of electromagnetic fields in space, in the spatial grid, the field components in each coordinate axis direction are separated by half a grid space step, so the same field components are separated by one space step: the electric field and the magnetic field are separated by half a time step in time sequence, so that Maxwell’s equation can form an explicit difference equation after being dispersed, and then the electromagnetic problem can be solved iteratively in time after the initial value is given.

Assume that the radial distribution of the perturbed space charge field variable (such as $Ex$) is of the form $J_0\left(k_{vt}\mathrm{r}\right)$, where $k_{vt}=2.4/r_b$, $r_b$ is the beam radius. This mode indicates that there is a conductor at the edge of the beam. Then, we use Poisson’s equation to calculate the disturbing space charge potential $\varphi$, density $\rho$ and the cylinder is shown in Formula (3):

$$\frac{{\partial }^2\varphi }{\partial x^2}-k_{vt}^2\varphi =-\frac{\rho }{e_0}$$

(3)

In the formula $\varphi (\mathrm{x},\mathrm{r},\theta )=\varphi (\mathrm{x})J_0$. The next section discusses the finite difference decomposition of this equation using the PIC technique.

In Equation (3), choosing to propagate as $e^{j(kx-\omega t)}$, one can obtain

$$\varphi (k)=\frac{1}{{\epsilon }_0}\frac{1}{k^2}\left(\frac{k^2}{k^2+k_{vt}^2}\right)\rho (k)$$

(4)

This equation shows the usual space charge reduction factor, where the effective plasma frequency is $W_q=R{W}_q$,

$$R\left(k_{pl},k_{vt}\right)\equiv {\left(\frac{k_{pl}^2}{k_{pl}^2+k_{vt}^2}\right)}^{1/2}\le 1$$

(5)

This coefficient gives a simple explanation for the reduction of the axial field and the presence of the transverse field. The “best” choice of $k_{vt}$ depends on the type of slow-wave circuit chosen. As mentioned above, if the slow-wave circuit is a drift tube located at the edge of the beam, we choose to use $k_{vt}$, which is probably the best choice for a coupled-cavity circuit.

4. PIC Method

The wave-equation finite-difference time-domain method (WEFDTD) is one of the important methods to simulate electromagnetic wave propagation. Its calculation step is a time-stepping, that is, frog-leaping iterative calculation carried out along the time axis. Therefore, the calculation steps of WEFDTD can also be spatially stepped, i.e., a frog leaping iterative calculation along a certain-spatial axis can be performed. Obviously, the two algorithms have their own advantages, and each has its own best problem to solve, so the PIC method of the frog-leaping step scheme is adopted. Each electron advances at each time step through the discrete form of the equation of motion. The acceleration $a$, velocity $v$ and position equations used in the specification are similar to the following equations:

$$\begin{array}{cc}\frac{dv}{dt}=a\Rightarrow & v\left(t_{\mathrm{n}+1}\right)=v\left(t_{\mathrm{n}}\right)+a_{ar}\mathsf{\Delta }t\\ \frac{dx}{dt}=v\Rightarrow & x\left(t_{n+1 }\right)=x\left(t_n\right)+v_{ar}\mathsf{\Delta }t\end{array}$$

(6)

In it, $t_{n+1}=t_n+Ar$, $\mathsf{\Delta }t$ is the time step in seconds, the subscript “$Ar$” represents the average of the variable between time $t_n$ and time $t_{n+1}$. However, there are two differences between this equation and the equation used in the algorithm; first, the electron velocity and position are evaluated at different times, separated by $v_{ar}$, so if $\mathsf{\Delta }t$ is small enough, $v\left(t_n+\mathsf{\Delta }t/2\right)$ will be close to $v_{ar}$, and the same is true for $a_{ar}$. Second, the algorithm uses normalized units, so $x_{\mathrm{c}}=x/\mathsf{\Delta }x$ and $v_c=v(\mathsf{\Delta }t/\mathsf{\Delta }x)$, where $\mathsf{\Delta }x$ is the cell width in meters ($\mathsf{\Delta }x=$ Tube length, $L/$ number of cells, $N_C$) and $v_{\mathrm{c}}=v(\mathsf{\Delta }t/\mathsf{\Delta }x)$.

For electric potential and electric fields,

$$V_c=V\frac{e}{m}\frac{{(\mathsf{\Delta }t)}^2}{{(\mathsf{\Delta }x)}^2}\hspace{1em}{E}_c=E\frac{e}{m}\frac{{(\mathsf{\Delta }t)}^2}{\mathsf{\Delta }x}$$

(7)

where e is the magnitude of the electron charge and m is the electron mass. The equation of motion used in the algorithm is on the right ($E$ is derived from ${\mathrm{q}}_c=-e$):

$$\begin{array}{l}{v}_m^{t+1/2}=v_m^{t-1/2}-\frac{e}{m}{E}_m^t\frac{{(\mathsf{\Delta }t)}^2}{\mathsf{\Delta }x}\Rightarrow \hspace{1em}{v}_m^{t+1/2}\equiv v_m^{t-1/2}-E_m^t\hfill \\ x_m^{t+1}=x_m^t+v_m^{t+1/2}\mathsf{\Delta }t\hspace{1em}\Rightarrow \hspace{1em}{x}_m^{t+1}\equiv x_m^t+v_m^{t+1/2}\hfill \end{array}$$

(8)

where “$m$“ represents the particle and “$t$“ represents the time step.

Variables are evaluated at the boundary or center of the $N_C$ unit in the simulation. Approximate the space charge density $\rho$ at the cell boundaries by assigning the charge on each electron to the two nearest cell boundaries by linear interpolation (linear weighting or intracellular particles). If particle $m$ is located between grid points $i$ and $i+1$, then the charge $Q_m$ on particle $m$ will be distributed so that the charge $Q_i\equiv Q_m\frac{x_{i+1}-x_m}{x_{i+1}-x_i}=Q_m\frac{x_{i+1}-x_m}{\mathsf{\Delta }x}$ assigned to the $i$ grid point, and $Q_{i+1}\equiv Q_m\frac{x_m-x_i}{\mathsf{\Delta }x}$ is assigned to $(i+1)$ grid point. Figure 3 shows the weights. After weighing all the particles in the tube, normalize the space charge to the standard or the number of electrons injected per $\mathsf{\Delta }t$ and per particle $N_E$. A uniform positive ionic space charge is then added to yield a net charge per unit length. Since $Q$ is the charge per cell, it means the integral of the radial dimension. Therefore, the difference in the radial density distributions of the RF and DC components of the space charge does not change $Q_i$ due to charge conservation.

A similar weighting method is used to approximate the beam current $I_b$; in the same way as the circuit current, $I_b$ is obtained at the center of each particle, at the $(i+1/2)-th$ grid point and half the time $(t-1/2)\Delta t$. At same as speed, each electron supplies a current to the two nearest cell centers via bilinear interpolation. If at time $(t-1/2)\Delta t$, particle $m$ is located between grid points $(i-1/2)$ and $(i+1/2)$, the current produced by that particle is:

$$I_{b,i\pm 1/2}^{t-1/2}\equiv Q_m{v}_m^{t-1/2}\left[\frac{1}{2}\pm \frac{x_m^{t-1/2}-x_i}{\mathsf{\Delta }x}\right]$$

(9)

$$x_m^{t-1/2}\equiv \frac{x_m^t+x_m^{t-1}}{2}=x_m^t-\frac{v_m^{t-1/2}\mathsf{\Delta }t}{2}$$

(10)

The weights are given in Figure 4. After the contributions of all electrons are accumulated, the current is normalized to the unit of the injected beam current.

The space charge potential is found by solving a diagonal matrix given by the finite difference form of Equation (3):

$${\varphi }_{i+1}-2{\varphi }_i+{\varphi }_{i-1}-\frac{{2.405}^2}{r_{bm}^2}{\varphi }_i=-\frac{{\rho }_{b,i}}{{\epsilon }_0}{(\mathsf{\Delta }x)}^2$$

(11)

$${\varphi }_{i+1}-2{\varphi }_i+{\varphi }_{i-1}-\frac{{2.405}^2}{r_6^2}{\varphi }_i=-{\left({\omega }_0\mathsf{\Delta }t\right)}^2{Q}_{b,i}$$

(12)

Due to the continuity of the beam charge, different drive terms can be used,

$$\frac{\partial I_b}{\partial x}=-\frac{\partial Q_b}{\partial t}$$

(13)

In our proposed algorithm, the telegraph operator equation has been transformed into a finite difference equation, which can also be solved by leap-frog time stepping. Use the following:

$$I_{\mathrm{c}}\equiv I\frac{e}{m{\epsilon }_0}{\left(\frac{\mathsf{\Delta }t}{\mathsf{\Delta }x}\right)}^3$$

(14)

$$I_{b0}\equiv I_{\mathrm{bm}}\frac{e}{m{\epsilon }_0}{\left(\frac{\mathsf{\Delta }t}{\mathsf{\Delta }x}\right)}^3$$

(15)

$$Q_{b0}\equiv \frac{I_{\mathrm{bm}}}{v_{\mathrm{bm}}}\frac{e}{m{\epsilon }_0}{\left(\frac{\mathsf{\Delta }t}{\mathsf{\Delta }x}\right)}^2$$

(16)

$$C_{\mathrm{c}}\equiv C\equiv \frac{C_0}{{\epsilon }_0}$$

(17)

$$L_{\mathrm{c}}\equiv L\equiv L_0{\epsilon }_0{\left(\frac{\mathsf{\Delta }t}{\mathsf{\Delta }x}\right)}^2$$

(18)

The equations used in the algorithm are:

$$V_i^t=V_i^{t-1}-\frac{1}{C}\left(I_{i+1/2}^{t-1/2}-I_{i-1/2}^{t-1/2}\right)+\kappa (k)\frac{Q_{b0}}{C}\left(Q_{b,i}^t-Q_{b,i}^{t-1}\right)$$

(19)

$$I_{i-1/2}^{t+1/2}=I_{i-1/2}^{t-1/2}-\frac{1}{L}\left(V_i^t-V_{i-1}^t\right)$$

(20)

Alternatively, other driver items can be used:

$$-\kappa \frac{I_{b0}}{C}\left(I_{b,i+1/2}^{t-1/2}-I_{b,i-1/2}^{t-1/2}\right).$$

(21)

Figure 5 shows the overall model in time and space.

The wave phase velocity $v_p$ is determined by $\Delta x/\Delta t$, so

$$L_{\mathrm{c}}=Z_{\mathrm{Oc}}={\epsilon }_0{v}_p{Z}_0$$

(22)

Figure 6 shows the model of transmitted beam in TWT. Due to the use of the finite differential, there may be chromatic dispersion in the transmission line. Remove the driving term from Equation (19) and assume that the propagation of $V_C$ and $J_C$ is $e^{j(kx-\omega t)}$.

In the latter, the subscript of $r_b=r_{\mathrm{bm}}/x$, $Q_{b,i}$ indicates that it is the beam space charge, not the charge in the transmission line, and ${\omega }_0$ is the electron plasma frequency of the beam:

$$I_{\mathrm{bm}}=e{n}_0{v}_{\mathrm{bm}}\pi r_{\mathrm{bm}}^2,\hspace{1em}\mathrm{or}\hspace{1em}{n}_0\cong 2\times {10}^{18}\frac{I_{\mathrm{bm}}}{v_{\mathrm{bm}}{r}_{\mathrm{bm}}^2}$$

(23)

Both the circuit and the space charge field are determined by the following finite difference method:

$$E=-\nabla V\hspace{1em}\Rightarrow \hspace{1em}{E}_{i+1/2}=\frac{V_i-V_{i+1}}{\mathsf{\Delta }x}\hspace{1em}\Rightarrow \hspace{1em}{E}_{i+1/2}=V_i-V_{i+1}$$

(24)

Among them, the last is used in the algorithm and the same is used for $\varphi$, each variable has a corresponding scale. The area sum is given by:

$$E_{(sum)i+1/2}=E_{(cc)\mathrm{i}+1/2}+SCcp\ne E_{(\mathrm{oc}-\mathrm{cg})i+1/2}$$

(25)

The field $E_m$ for each electron is then obtained from the two nearest particle centers $E_{i-1/2}$ and $E_{i+1/2}$ using linear weighting, so:

$$E_m=E_{i-1/2}\frac{x_{i+1/2}-x_m}{\mathsf{\Delta }x}+E_{i+1/2}\frac{x_m-x_{i-1/2}}{\mathsf{\Delta }x}$$

(26)

This is the field used to move electrons in the next time step.

Finally, the boundary conditions for the particles are the source and collector absorption. For example, if a particle passes the collector coordinates, it will be removed from the array that stores the particle’s position and velocity, and the number of electrons will be reduced by one. The beam particles are injected into the source uniformly in time. Furthermore, the space charge potentials of both the source and collector are fixed to zero.

5. Transmission Line Algorithms

The circuit in Figure 1 can be written with loop and nodal equations, resulting in the telegraph operator’s equation plus the driving term, under the constraints of small $\Delta x$ and $\Delta t$:

$$\begin{array}{l}\frac{\partial V}{\partial t}=-\frac{1}{C_0}\frac{\partial I}{\partial x}+\frac{1}{C_0}\frac{\partial Q_b}{\partial t}\hfill \\ \frac{\partial I}{\partial t}=-\frac{1}{L_0}\frac{\partial V}{\partial x}\hfill \end{array}$$

(27)

At the source location and collector location:

$$I_{N_C+1/2}^{t+1/2}=\left(V_{\mathrm{o}}-V_{N_C}-\left(R_{\mathrm{o}}-L\right)I_{N_C+1/2}^{t-1/2}\right){\left(R_{\mathrm{o}}+L\right)}^{-1}$$

(28)

where $R_{\mathrm{i}}={ϵ}_0{v}_{\mathrm{p}}{R}_{\mathrm{i}}$, $R_{\mathrm{i}}={ϵ}_0{v}_{\mathrm{p}}{R}_{\mathrm{i}}$. However, these terminations are matched (reflection free) only if there is no beam circuit coupling. The greater the coupling between the beam and the circuit, the greater the reflected signal, which pushes the loop gain beyond unity and the tube becomes an oscillator. To avoid self-oscillation, a segment of the transmission line is lossily compensated using the circuit segment shown in Figure 7.

There are two options included in the algorithm. The first scheme has a smooth section located near the center of the line, starting from grid point $N_0$, the spatial variation of the resistance is as follows:

$$R_i=R_0(1-\cos \theta )$$

(29)

In Formula (29), $\theta =2\pi \frac{i-N_0}{N_C-2N_0}$ and $N_0<i<N_C-N$. This loss is to attenuate the reflected signal before it reaches the input of the tube, while reducing the gain by a small amount. As the reflected wave approaches the lossy cross-section, the wave reflected from the collector end of the tube (which may be up to 10% of the incident signal) becomes larger relative to the forward (growing) wave and these two waves produce interesting effects in the lossy section the interferogram.

To overcome these limitations, a nearly completely absorbing or non-reflective termination section was designed. The termination is a lossy line at the pipe end of length $N_1$ units, where

$$R_i=R_t\left({(1+\gamma )}^{i-N_C+N_1}-1\right)$$

(30)

where $R_t\le {10}^{-5}$, $\gamma \le 0.2$ and $i>N_C-N_1{}_0$. In this section, the coupling between the beam and the circuit is also reduced in a tapered manner. It is not possible to abruptly terminate the coupling between the beam and the circuit, as this would also introduce reflections. The current decreasing loss form is as follows:

$$\cos \left(2\pi \left(i-N_C+N_1\right)/4N_1\right)$$

(31)

where $x=\Delta x$ and $t=\Delta t$, and Equations (18) and (19) become:

$$\mathsf{\Omega }{V}_c=\frac{K}{C}{I}_c$$

(32)

$$\mathsf{\Omega }{I}_c=\frac{K}{C}{V}_c$$

(33)

In the above formula, $\mathsf{\Omega }=\omega \sin c(\omega \Delta t/2)$ and $K=k\sin c(k\Delta x/2)(\sin cy\equiv \sin y/y)$. The elimination of $V_{\mathrm{C}}$ and $I_{\mathrm{C}}$ produces transmission line dispersion relations:

$$K^2={\mathsf{\Omega }}^{2}LC$$

(34)

This equation shows two results. First, it should be written as:

$${\sin }^2(\omega \mathsf{\Delta }t/2)=\frac{\mathsf{\Delta }{t}^2}{LC\mathsf{\Delta }{x}^2}{\sin }^2(k\mathsf{\Delta }x/2)$$

(35)

The right-hand side of the equation must be less than 1 for $\omega$ to be correct; therefore, for a stable solution, we must make $\mathsf{\Delta }t\le {(LC)}^{1/2}\mathsf{\Delta }x$ as a condition. Formula (5) is then:

$$v_p^2\equiv \frac{{\omega }^2}{k^2}=\frac{1}{LC}\frac{sin{c}^2(\omega \mathsf{\Delta }t/2)}{sin{c}^2(k\mathsf{\Delta }x/2)}$$

(36)

This is the case if and only if $v_p=\Delta x/\Delta t$, $v_p={(LC)}^{-1/2}$.

Termination and Transmission Line Loss

The transmission line is terminated with a resistor and a voltage source at each end. Since the input and output source voltages $V_i$ and $V_o$, respectively, are user-defined functions of time, only the updated circuit current needs to be calculated. The passive components connected to the transmission line are resistors $R_i$ and $R_o$, which are input and output, respectively, and an inductor, $L_0/2$, is connected in series with the input and output, respectively. The resulting equation for the second half is:

$$I_{-1/2}^{t+1/2}=\left(V_i-V_0-\left(R_{\mathrm{i}}-L\right)I_{-1/2}^{t-1/2}\right){\left(R_{\mathrm{i}}+L\right)}^{-1}$$

(37)

However, when $N_C=300$, a gain of about 60 dB is achievable before the reflected signal becomes comparable to the input signal. The typical parameter values used are $N_1=80$, $R_t={10}^{-6}$ and $\gamma =0.1$.

6. Case Simulation

The code is mainly developed for beam–wave interaction analysis in cylindrical TWT with symmetrical structure using MATLAB software. Because the applicable structure of this algorithm is fixed, it does not have the simulation function for irregular TWT, and the fixed structure makes its simulation speed faster.

Next, the electromagnetic solution algorithm in TWT is tested, and the electron beam with beam voltage of 20 kV and 50 kV and input power of 100 W is added to the cylindrical TWT with symmetrical structure and center frequency of 2 × 10¹¹ Hz for simulation. The same input conditions are set and the relatively mature CHPIC electromagnetic simulation software for simulation analysis and comparison are used.

As shown in Figure 8, when the beam voltage in the algorithm is 20 kV, the circuit electric field distribution diagram and CHPIC code results generated at different positions in the traveling wave tube are compared.

However, if the beam voltage is increased to 50 kV, the electric field distribution law has changed greatly. As shown in Figure 9, when the beam voltage is 50 kV, the circuit electric field distribution diagram and CHPIC code results generated at different positions in the traveling wave tube are compared.

When the beam voltage in the algorithm is 20 kV, Figure 10 show the electron tracking using the space charge field and circuit field. It can be seen intuitively that the track of 20 kV electrons with beam voltage changes with time compared with the CHPIC code results.

Figure 11 also shows the electron tracking trajectory of space charge field and circuit field in the calculation example when the beam voltage is 50 kV, and the comparison effect with CHPIC code.

Figure 12 shows the simulation results of output voltage changes with different positions in TWT when the beam voltage is 20 kV, and the results are compared with the simulation results of CHPIC code.

Table 1. Comparison of the output voltage results of the proposed algorithm and CHPIC code generated by a beam–wave interaction simulation under 20 kV RF voltage.

	Distance (mm)	Output Voltage of TWT Head End (kV)	Distance (mm)	Maximum Output Voltage in TWT (kV)	Distance (mm)	Minimum Output Voltage in TWT(kV)
Algorithm testing	0	7.2930	0.1143	13.1300	0.0726	0.3168
CHPIC code generation result	0	9.1650	0.1126	14.3400	0.0699	0.3641

compares the algorithm test results and CHPIC code results in detail, the output voltage value generated at the TWT inlet, and the maximum and minimum output voltage values corresponding to different positions in the TWT.

Figure 13 shows the simulation results of output voltage changes with different positions in TWT when the beam voltage is 50 kV, and the results are compared with the simulation results of CHPIC code.

Table 2. Comparison of the output voltage results of the proposed algorithm and CHPIC code generated by a beam–wave interaction simulation under 50 kV RF voltage.

Distance (mm)	Output Voltage of TWT Head End (kV)	Distance (mm)	Maximum Output Voltage in TWT (kV)	Distance (mm)	Minimum Output Voltage in TWT(kV)	Distance (mm)
Algorithm testing	0	15.5800	0.0739	22.0200	0.1473	0.5740
CHPIC code generation result	0	17.2200	0.0734	22.1300	0.1484	0.8620

compares the algorithm test results and CHPIC code results in detail, the output voltage value generated at the TWT inlet, and the maximum and minimum output voltage values corresponding to different positions in the TWT.

Figure 14 shows the simulation results that the output power caused by beam-wave action changes with different positions in the TWT when the beam voltage is 20 kV, and the results are compared with the CHPIC code simulation results.

Table 3. Comparison of the output power results of the proposed algorithm and CHPIC code generated by a beam–wave interaction simulation under 20 kV RF voltage.

	Distance (mm)	Output Power of TWT Head End (W)	Distance (mm)	Maximum Output Power in TWT (W)	Distance (mm)	Minimum Output Power in TWT (W)
Algorithm testing	0	169.5	0.0301	495.3	0.1152	549.3
CHPIC code generation result	0	267.6	0.0272	592.1	0.1132	655.0

compares the algorithm test results and CHPIC code results in detail, the output power value generated at the TWT inlet, and the maximum and minimum output power values corresponding to different positions in the TWT.

Figure 15 shows the simulation results that the output power caused by beam-wave action changes with different positions in the TWT when the beam voltage is 50 kV, and the results are compared with the CHPIC code simulation results.

Table 4. Comparison of the output power results of the proposed algorithm and CHPIC code generated by a beam–wave interaction simulation under 50 kV RF voltage.

		Distance (mm)	Output Power of TWT Head End (W)	Distance (mm)	Maximum Output Power in TWT (W)	Distance (mm)
Algorithm testing	0	794.8	0.0745	1590.0	0.1473	0.574
CHPIC code generation result	0	964.4	0.0736	1597.0	0.1484	2.491

compares the algorithm test results and CHPIC code results in detail, the output power value generated at the TWT inlet, and the maximum and minimum output power values corresponding to different positions in the TWT.

It can be seen from Figure 8, Figure 9, Figure 10 and Figure 11 that the electric field distribution produced by the proposed algorithm and the commercial CHPIC code under the beam–wave action of 20 kV and 50 kV RF voltage coincides well, and the consistency of electron trajectory changes is good. Table 1 and Table 2 show some significant data of the output voltage shown in Figure 12 and Figure 13. No matter the voltage at the head of the TWT, or the highest and lowest voltage values generated at different locations, they are all in good agreement. The same is true for the output power data shown in Figure 14 and Figure 15 and given in Table 3 and Table 4.

7. Conclusions

In order to solve many nonlinear and strong coupling problems in physics based on the background of high current electron beam and electromagnetic static field, a nonlinear beam–wave analysis method in the frequency domain is proposed for analyzing the electron beam–wave problem of cylindrical traveling wave tube. The basic working equation of the coupling between the beam and the circuit is given and considering the space charge effect caused by the large-amplitude signal close to saturation, the electric field at the center of each particle is given by the PIC method using the leap-frog step scheme. The specific steps for evaluating and determining the circuit and space charge field are presented, and the design principle of the termination part of the algorithm for lossy compensation in order to avoid self-excitation is introduced. This algorithm code generates the distribution diagram of electric field, electron trajectory, output voltage and output power variables in TWT, providing data for solving problems in physics. The design of the algorithm only uses a single cylindrical structure. Although it speeds up the simulation, the simulation of irregular structures is still lacking and needs further improvement.