Gyrotron is a microwave oscillators. It changes the DC energy to high power millimetre microwave. It has been successfully used in the magnetic confinement thermonuclear fusion research as the microwave source for the electron cyclotron resonance heating and current drive.
A typical gyrotron cavity consists of a central section and one taper section at both ends, as shown in Fig. 1. The central region is usually a regular hollow metallic cylinder with a circular cross section. The electron beam gives part of its kinetic energy to the electromagnetic field here. The input taper is a cutoff section which prevents the back propagation of the electromagnetic wave to the gun region. The output taper is a horn like. It allows part of the radiation to be coupled to the wave guide. Reflections at tapers lead to resonant. The Quality factor describes how damped the resonator is. The resonance frequency and the Q factor are determined by the radiation boundary conditions.
Vlasov Method
A traditional approach of the gyrotron cavity is developed by Vlasov et al[1]. It assumes that the transverse mode doesn’t changed and only its amplitude changes along the axis of the cavity. So the transverse field \(E_t\) could be expand as the superposition of all the possible TE modes and TM modes.
\(E_t=\sum V_{mp}(z) \vec e_{mp}(r,\theta)+ \sum U_{mp}(z) \vec h_{mp}(r,\theta)\)
Here, m and p are the azimuthal and radial mode number respectively. \( \vec e_{mp}(r,z,\theta)\) and \(\vec h_{mp}(r,z,\theta)\) are the TE modes and TM modes basis functions. \(V_mp(z)\) and \(U_mp(z)\) are the axial distribution of the complex field amplitude along z for each mode. The perfect electrical boundary condition on both taper results in the modes coupling. The mode coupling coefficients can be obtained and they are very complicated. For small taper angles, the modes decouple. Only a single mode exists in the cavity. This is called the Vlasov approximation. Then the auxiliary function \(h(r,\theta\,z)\) for a single mode mp is
\(h(r,\theta,z)=V(z)\vec e(r,\theta)\, z)
It is governed by the homogeneous Helmholtz equation
\(\nabla^2 V(z)e(r,\theta)+k^2 V(z) e(r,\theta)=0\)
evaluate the differential operator and divide \(V(z)e(r,\theta)\) on both sides, it turns
\(\frac{\nabla_t^2 e(r,\theta)}{e(r,\theta)}+\frac{d^2 V(z)}{dz^2}+k^2=0\)
where \(\nabla_t^2\) is the Laplacian operator in the transverse direction. \(k=\omega/c\) is the wave number, c is the speed of light.
The first term is only dependent on \(r, \theta\) and the second term is only dependent on z. So the equation is separable. Let the first term equal \(-k_{mp}^2\) and the second term equal \(-k_z^2\), so \(k_z^2=k^2-k_{mp}^2\), and
\(\nabla_t^2 e(r,\theta)+k_{mp}^2 e(r,\theta)=0\,\,\,\,\, \frac{d^2V(z)}{dz^2}+(k^2-k_{mp}^2)V(z)=0\)
The 1st equation is the traditional Helmholtz equation in the transverse direction. It is well known that the transverse wave number \(k_{mp} \) is \( \frac{\chi_{mp}}{R(z)}\), \(\chi_{mp}\) is the p-th root of \(J’_m(x)\) for TE mode. \( R(z)\) is the radius of the cavity.
Now we focus on the the second equation. It is the Helmholtz equation in the axial direction and analogous to the Schrodinger’s equation for a particle in a potential well with a potential \(k_z^2\). For a stable resonate mode, there is no reflection on both ends. The wave that reaches the output will continue propagate along z-axis without any reflection and the wave that reaches the input will continue propagate along -z axis without any reflection. It means
\(\frac{dV}{dz}-i k_z V=0\) at \(z=z_{in}\) and \(\frac{dV}{dz}+i k_z V=0\) at \(z=z_{out}\)
{k_z<0} at the input, the square root has two values. As it is an evanescent wave at \(z_{in}\), it is
\(\frac{dV}{dz}= |k_z| \) at \(z=z_{in}\)
There are some analytical solutions. However they are all too complex to be useful for me. For numerical calculations to find the frequency and the Q factor, the boundary condition at \(z_{in}\) is used to determine the starting value of \(dV/dz\) and the \(v(z_{in})\) is set to a small value. It is an inverse problem. It implies an iterative procedure in which the Helmholtz equation is integrated from the input taper to the output taper for a sequence of guess values and check the boundary condition at the output. There are two important procedures: the integration method and the searching method.
A fast eigen frequency search procedure
Most frequently, the search procedures are based on the MINUIT software package for a nonlinear optimization with continuous parameters or complex root finding code based on the Muller’s method. The search procedures used to calculate the complex frequency are in fact minimization algorithms and the iterative process continues until the reflection at the output port is small enough. This numerical approach has been implemented in many computer codes developed by different research groups worldwide and used successfully for analysis and design of gyrotron cavities.
Both the local optimization and the global optimization in the mathematica were tried. The local optimization methods with a start guess search the eigen frequency successfully in about 8s sometimes. The best algorithm is Principal Axis. It is impossible to calculate the symbolic derivatives of the reflection. The finite differences will be uesd in the Gauss-Newton and conjugate gradient methods. Computing derivatives with finite differences can impose a significant cost in this case and certainly affects the reliability of derivatives. The principal Axis method is a derivative-free algorithm. For an 2-variable problem, take a set of search direction and a start point. Move to the best direction and search again. However, the search result dependent strongly on the start point. The eigen frequency is not got sometimes. The best global optimization is Nelder-Mead. It usually takes about 30s to find the eigen frequency successfully. These methods are suite for the general optimization problem. They don’t know anything about this problem. We develop a fast search procedure based on the property of the Helmholtz equation.
The property of the reflection
The reflection at the output taper is a function of the complex frequency. The complex frequency is always description as the frequency and the quality factor. Here is a contour plot of the reflection. There are two eigen frequency in the plot region because there are two hole in the contour plot. The equal reflection line becomes curve around the eigen frequency. It seems all the curves have a peek at the eigen frequency.
In order to verify this guess, five reflection curves are plotted at difference quality factor. All the curves have the same local minimum value at the eigen frequency. The eigen frequency is decoupled from the quality factor. We search the eigen frequency at a very high quality factor 5000. The first local minimum frequency is the eigen frequency if we search from the cutoff frequency.
The important property is that the frequency and the quality factor are decoupled. The search procedure is
- Move from the cutoff frequency with a step=0.1GHz until the reflection increases.
- Use three-section method to find the local minimum very quickly. Cut the searching domain into 3 equal part, compare the reflection at the inner 2 points. Drop the part in which the reflection is greater and so on.
step=0.1; //frequency step
a=geometry[f]*10^-9; //Get the cutoff frequency
b=a+step;
While[reflect[a,3000]>reflect[b,3000],a=a+step;b=b+step;];//move until the curve goes up.
a=a-step;
eps=0.0002; //the maximum value of the frequency difference.
While[b-a>eps,a1=(2 a+b)/3;b1=(2 b+a)/3;If[reflect[a1,2344]<reflect[b1,2344],b=b1,a=a1];]; //three section method
f0=(a+b)/2The procedure to search the quality factor after the eigen frequency is obtained is very similar. However we start at a very high value and search the curve from right to left.
step=200;b=5000.0;a=b-step;While[reflect[f0,b]>reflect[f0,a],a=a-step;b=b-step;]; //search the local minimum roughly
b=b+step;
eps=0.1;
While[b-a>eps,a1=(2 a+b)/3.0;b1=(2 b+a)/3.0;If[reflect[f0,a1]<reflect[f0,b1],b=b1,a=a1];];// search the quality factor exactly.
q0=(a+b)/2;The new search procedure runs much more quickly than the traditional procedure because the algorithm knows the property of this problem. It usually find the eigen frequency and the quality factor in 2s.
The magnitude of the field in the cavity
\(V(z)\) shows the filed profile. However the output power depends on the product of \(V(z)\) and \(\vec e(r,\theta)\). It is convenient to work with a normalized field profile. This will be done by defining
\(V_{mn}(z)=V_{max}\hat f_{mn}(z)\)
Here, \(\hat f_{mn}(z)\) is the field profile normalized to a maximum absolute value 1. The quality factor Q is defined in the resonate cavity as \(2\pi\) times the ratio of the time-averaged energy W stored in the cavity to the energy loss per cycle:
\(Q=2\pi \frac{W}{T P_{out}}=2 \pi f \frac{W}{P_{out}}=\omega\frac{W}{P_{out}}\)
where \(P_{out}\) is the output power of the cavity. So,
\(QP_{out}=\omega W=\omega (W_e+W_m)=2\omega W_e=2\omega\frac{\epsilon}{4}V_{max}^2\int_{z_{in}}^{z_{out}} |\hat f(z)|^2 \int_0^{2\pi}\int_0^{R(z)}r \vec e_{mn}\cdot \vec e*_{mn}drd\theta dz\)
If \(I=\int_0^{2\pi}\int_0^{R(z)}r \vec e_{mn}\cdot \vec e*_{mn}drd\theta=1\), then the relationship between the output power, quality factor and maximum amplitude is
\(QP_{out}=\frac{\omega\epsilon}{2}V_{max}^2\int_{z_{in}}^{z_{out}} |\hat f(z)|^2 dz\)
Let the scalar function \(\phi=C J_m(k_c x)e^{-j m\theta}\), where \(Hz=\phi e^{-j k_z z}\) for the TE mode in the waveguide.
\(\vec e_{mn}=\nabla \phi=C k_c J_m'(k_c r)e^{-j m \theta}\hat e_r-\frac{j m}{r}C J_m(k_c r) e^{-j m \theta}\hat e_{\theta}\)
\(\vec e_{mn}^*=\nabla \phi^*=C k_c J_m'(k_c r)e^{j m \theta}\hat e_r+\frac{j m}{r}C J_m(k_c r) e^{j m \theta}\hat e_{\theta}\)
\(I=\iint_A\nabla \phi \cdot \nabla \phi^* da=\int_0^{2\pi}d\theta\int_0^{R(z)} r C^2 k_c^2 J_m’^2(k_c r)+\frac{m^2}{r}C^2 J_m^2(k_c r) dr\)
\(I=2\pi C^2 \int_0^{R(z)} r k_c^2 J_m’^2(k_c r)+\frac{m^2}{r} J_m^2(k_c r) dr\)
Let \(x=k_c r\) then \(dx=k_c dr\)
\(I=2\pi C^2 \int_0^{P’_{mn}} x J_m’^2 (x)+\frac{m^2}{x} J_m^2 (x) dx= \pi C^2 (P’^2_{mn}-m^2)J_m^2(P’_{mn})=1\)
So, \(C_{mn}=\frac{1}{\sqrt{\pi (P’^2_{mn}-m^2)}|J_m(P’_{mn})|}\)
New method from LU Bo
Although the Vlasov method got great success in the last 50 years, the duration to search the mode is not always the same because the inverse-problem depend strongly on the initial conditions. There is also the risk that the search stops on a local optimize value. The Vlasov equation is not linear because the boundary condition. If the nonlinear boundary condition would have been avoided, then the question becomes linear.
The domain is expanded by adding materials next to the two boundaries to absorb the wave that radiate from the cavity. Then metals is placed behind the material. The boundary condition becomes perfect electric wall. It is linear. The absorbing materials are the perfect match layers. It is an well developed imaginary material in the calculated electromagnetic. Its relative permittivity and relative permeability are the same so that the wave impedance is the same as the impedance in the vacuum.
\(Z=\sqrt{\frac{\mu_0 \mu_r}{\epsilon_0\epsilon_r}}=\sqrt{\frac{\mu_0}{\epsilon_0}}=Z_0\)
Its electrical conductivity is not zero. The wave damps in PMLs. However, the electrical conductivity caused the impedance changed. The anisotropy is used to cancel the effect to the impendance from the electrical conductivity and the the anisotropy is realized in two dimension. The r-z plane is solving domain. We separate the variables as
\(\vec E(r,\theta,z)=\vec E(r,z) M(\theta)\)
apply it into the Helmholtz equation and separate the variable,
\(\frac{\partial^2 \vec E}{\partial r^2}+\frac{\partial^2 \vec E}{\partial z^2}+\frac{1}{r} \frac{\partial \vec E}{\partial r}+(k^2-\frac{m^2}{r(z)}) \vec E(r,z)=0\)
It is a linear partial differential equation in r-z plane. Its boundary condition is PEC and it is linear too.
\(\vec n\times \vec E=0\)
This problem is solved by the finite element method.
Example
The calculus is always seems complex. However an example is always simple. We next try an example in reference 2. The geometry parameters is
| \(\theta_1\) | \(\theta_2\) | L1 | L2 | L3 | \(R_0\) | \(f_0\) |
| 1 degree | 3 degree | 10mm | 16mm | 15mm | 3.24mm | 150GHz |
The geometry is
And the axial wave number square \(k_z^2\) is
The axial wave number at both end is \(|k_z(0)|=1069 \) and \(kz(0.041)=1866.7\). Without a source, the eigen-frequencies are complex with positive imaginary part to account the damping due to the radiation out of the cavity. The amplitude is arbitrary. When a source is applied to the cavity, the amplitude is determined by the source. Here, we assume the \(V(0)=1\times 10^{-6}\). The wave length \(\lambda\) is 2mm, the step size \(\Delta z\) is set to one 20th wave length, so
\(\Delta z=\lambda/20=0.0001m\)
and use the 1st order Taylor series and the boundary condition at the input
\(V(\Delta z)=V(0)+V'(0)*\Delta z=V(0)+|k_z(0)|*V(0)*\Delta z=\)1.11e-6
Next we apply Numerov algorithm[3] to integrate
\(V_{n+1}={\frac {\left(2-{\frac {5\Delta z^{2}}{6}}kz_{n}\right)V_{n}-\left(1+{\frac {\Delta z^{2}}{12}}kz_{n-1}\right)V_{n-1}}{1+{\frac {\Delta z^{2}}{12}}kz_{n+1}}}\)
The iteration is carried out by Microsoft Excel. Here is the file.
It is clear that there is a evanescent wave at the input and a propagation wave at the output. The field is very strong at the uniform section. V(z) is real. It is impossible to meet the boundary at the output. Then a search algorithm will be used to search the eigen frequency and Quality factor.
Numerical Tools
I developed two numerical tools to obtain the eigenmodes. One is based on the finite element method on r-z plane. The other one is the traditional vlasov method. The new method is very quick and provide the the vlasov method the initial condition.
The traditional method is realized in Mathmatica. Here is the main procedure to obtain the eigenvalue.
Step 1 Define the geometry
geometry = <|\[Theta]1 -> 5.0 Degree, \[Theta]2 -> 3.0 Degree,
L1 -> 0.014, L2 -> 0.016, L3 -> 0.012, R0 -> 0.01557, m -> 22,
n -> 6, f -> 140*10^9|>Step 2 Load the package and some constant is initialized.
<<Gyrotron'Step 3 Search the eigen-mode.
search[reflect[f, q], {f, 139.6, 140.1, 0.05}, {q, 600, 1600, 200}]An example file is in D:\document\old\电磁场\gyrotron\vlasov.nb
Summary
The eigen-mode in a gyrotron cavity is introduced. The Helmholtz equations of both the transverse and axial are derived in a single mode approximation. An integration example is carried out by using excel and the numerov algorithm.