4.1 Force equilibrium

The electric current density of the warm plasma $\boldsymbol {J}_{{\rm w}}$, which appears in (3.4), can be found from the force balance equation for the warm plasmaFootnote ⁴

(4.1)\begin{equation} \boldsymbol{\nabla} p=c^{-1}\boldsymbol{J}_{{\rm w}}\times\boldsymbol{B},\end{equation}

where $p$ is the warm plasma pressure.Footnote ⁵ It can be seen that, due to the stationary axisymmetric equilibrium being considered, the sum of the azimuthal friction forces acting on the entire warm plasma is equal to zero.

It immediately follows from (4.1) that the warm plasma pressure $p$ is the constant on the flux surfaces $\psi ={\rm const}.$ and, therefore, is a function of the magnetic flux $\psi$ only: $p=p(\psi )$. In addition, due to the high longitudinal electron thermal conductivity, the temperature of the warm plasma $T$ appears to be constant along the field lines. Together with axial symmetry, this results in the temperature also being a flux function: $T=T(\psi )$. Eventually, assuming the pressure, temperature and densityFootnote ⁶ of the warm plasma to be related by the equation of state

(4.2)\begin{equation} p=\left(1+Z\right)n_{{\rm i}}T,\end{equation}

we arrive at the same for the warm plasma density: $n_{{\rm i}}=n_{{\rm i}}(\psi )$.

As a result, (4.1) allows the warm plasma current density to be expressed in terms of the magnetic field and pressure gradient

(4.3)\begin{equation} J_{{\rm w}\theta}=2{\rm \pi} rc\frac{{\rm d}p}{{\rm d}\psi}.\end{equation}

The magnetic flux distribution $\psi =\psi (r,z)$ is determined by Maxwell's equations (3.4), and the pressure profile $p=p(\psi )$ should be found from the solution of the mass and energy conservation equations for the warm plasma.

4.2 Mass conservation

To obtain the warm plasma equilibrium equation, we consider the material balance in a flux tube $\psi =\mathrm {const}$.

(4.4)\begin{equation} 2\varPhi_{{\rm i}\parallel{\rm m}}+\varPhi_{{\rm i}\perp}-W_{{\rm i}}=0,\end{equation}

where $\varPhi _{{\rm i}\parallel {\rm m}}$ is the axial loss of the warm ions from the flux tube, determined by the longitudinal flow through a mirror throat, $\varPhi _{{\rm i}\perp }$ is the warm ion flow transverse to the flux surface $\psi =\mathrm {const}.$ and $W_{{\rm i}}=W_{{\rm i}}(\psi )$ is the total warm ion source inside the flux tube. In other words, $W_{{\rm i}}$ represents the number of warm ion–electron pairs supplied by an external source to the interior of the surface $\psi =\mathrm {const}.$ per unit time. Longitudinal and transverse ion fluxes are determined as follows:

(4.5a,b)\begin{equation} \varPhi_{{\rm i}\parallel{\rm m}}=\int_{S_{{\rm m}}\left(\psi\right)}n_{{\rm i}}\boldsymbol{u}_{{\rm i}}\boldsymbol{\cdot} {\rm d}\boldsymbol{S},\quad\varPhi_{{\rm i}\perp}=\int_{S_{{\perp}}\left(\psi\right)}n_{{\rm i}}\boldsymbol{u}_{{\rm i}}\boldsymbol{\cdot} {\rm d}\boldsymbol{S}, \end{equation}

where $\boldsymbol {u}_{{\rm i}}$ is the warm ion flow velocity. Integration is carried out over the cross-section of the magnetic surface in the mirror $S_{{\rm m}}(\psi )$ and the flux surface $S_{\perp }(\psi )$ between the mirrors, respectively.

Gas-dynamic confinement implies the axial mirror loss being described by an outflow through the nozzle

(4.6)\begin{equation} \varPhi_{{\rm i}\parallel{\rm m}}=\int_{S_{{\rm m}}\left(\psi\right)}n_{{\rm i}}u_{{\rm m}}\,{\rm d}S,\end{equation}

where $u_{{\rm m}}\sim \sqrt {T/m_{{\rm i}}}$ is the warm plasma flow velocity in the mirror throats. The ion flow velocity transverse to the flux surface $\psi =\mathrm {const}.$ can be found from the generalized Ohm's law, which can be written in the form

(4.7)\begin{equation} \frac{1}{c}\left[\boldsymbol{u}_{{\rm i}}\times\boldsymbol{B}\right]_{\theta}=\frac{J_{{\rm w}\theta}}{\sigma _{{\rm w}}},\end{equation}

where $\sigma _{{\rm w}}$ is the electrical conductivity of the warm plasma.Footnote ⁷ The force balance (4.1) together with the azimuthal projection of Ohm's law (4.7) yields the warm plasma transverse flow velocity

(4.8)\begin{equation} \boldsymbol{u}_{{\rm i}\perp}=-\frac{c^{2}}{\sigma _{{\rm w}}B^{2}}\boldsymbol{\nabla} p.\end{equation}

Then the transverse ion flux is

(4.9a,b)\begin{equation} \varPhi_{{\rm i}\perp}=-\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi},\quad \varLambda_{{\rm i}\perp}=4{\rm \pi}^{2}c^{2}n_{{\rm i}}\int_{\gamma_{{\perp}}\left(\psi\right)}\frac{r^{2}\,{\rm d}l}{ \sigma _{{\rm w}}B},\end{equation}

where integration is carried out in $(r,z)$ space along the curve $\gamma _{\perp }(\psi )$ corresponding to the flux surface $\psi =\mathrm {const}$.

Finally, substituting (4.6) and (4.9a,b) into the material balance (4.4) and differentiating it with respect to $\psi$, we obtain the mass conservation equation for the warm plasma in the following form:

(4.10)\begin{equation} {-}\frac{{\rm d}}{{\rm d}\psi}\left(\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}\right)+\frac{2u_{{\rm m}}}{B_{{\rm m}}}n_{{\rm i}}= \frac{{\rm d}W_{{\rm i}}}{{\rm d}\psi},\end{equation}

where $B_{{\rm m}}$ is the magnetic field in the mirror throats; here, we also take into account that ${\rm d}\psi =B_{{\rm m}}\,{\rm d}S_{{\rm m}}$.

4.3 Energy conservation

For a given distribution of the magnetic flux, the equilibrium state of the warm plasma is fully described by three parameters: density $n_{{\rm i}}$, temperature $T$ and pressure $p$. A closed system of equations can be obtained by supplementing (4.2) and (4.10) with another equation relating these three parameters. In the previous MHD models (see Beklemishev Reference Beklemishev2016; Khristo & Beklemishev Reference Khristo and Beklemishev2019, Reference Khristo and Beklemishev2022), the plasma temperature $T$ is simply assumed to be constant. However, introducing an equation describing the energy balance of the warm plasma would be more proper.

The law of thermodynamics for a moving warm plasma element reads

(4.11)\begin{equation} \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\frac{3}{2}p\boldsymbol{u}_{{\rm i}}\right)= \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\varkappa_{{\rm w}}\boldsymbol{\nabla} T\right)+\frac{J_{{\rm w}}^{2}}{\sigma _{{\rm w}}}+q_{{\rm h}}-p\boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{u}_{{\rm i}},\end{equation}

where $-\varkappa _{{\rm w}}\boldsymbol {\nabla } T$ is the heat flux density due to thermal conductivity with a coefficient $\varkappa _{{\rm w}}$, $q_{{\rm h}}$ is the density of the heating power from the hot ions (see expression (5.18) in § 5). The left-hand side of the equation is the change in the internal energy of the moving element with time. The first three terms on the right-hand side describe the total heat release, and the last term is related to the mechanical work. Here, we consider the heating to be due to the neutral beam injection. If necessary, heating from other sources can also be taken into account by simply adding the corresponding terms to the right-hand side of the equation.

Equations (4.1) and (4.8) together yield the resistive heating power density

(4.12)\begin{equation} \frac{J_{{\rm w}}^{2}}{\sigma _{{\rm w}}}=-\boldsymbol{u}_{{\rm i}}\boldsymbol{\cdot}\boldsymbol{\nabla} p. \end{equation}

Substituting it into (4.11) we arrive at

(4.13)\begin{equation} \tfrac{5}{2}\boldsymbol{\nabla}\boldsymbol{\cdot}\left(\,p\boldsymbol{u}_{{\rm i}}\right)- \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\varkappa_{{\rm w}}\boldsymbol{\nabla} T\right)=q_{{\rm h}}. \end{equation}

Integrating this equation over the volume of the flux tube $\psi =\mathrm {const}.$ between the mirrors results in

(4.14)\begin{equation} 2\varPhi_{{\rm E}\parallel}+\varPhi_{{\rm E}\perp}=Q_{{\rm h}},\end{equation}

where $\varPhi _{{\rm E}\parallel }$ is the axial energy loss from the flux tube

(4.15)\begin{equation} \varPhi_{{\rm E}\perp}=\int_{S_{{\perp}}\left(\psi\right)}\frac{5}{2}p\boldsymbol{u}_{{\rm i}}\boldsymbol{\cdot} {\rm d}\boldsymbol{S}-\int_{S_{{\perp}}\left(\psi\right)}\varkappa_{{\rm w}}\boldsymbol{\nabla} T\boldsymbol{\cdot} {\rm d}\boldsymbol{S} \end{equation}

is the energy flow transverse to the flux surface $\psi =\mathrm {const}.$, and $Q_{{\rm h}}$ is the total heating power from the hot ions inside the flux tube.

The axial energy loss in a GDT can be described as follows:

(4.16)\begin{equation} \varPhi_{{\rm E}\parallel}=\alpha_{{\rm E}}\int_{S_{{\rm m}}\left(\psi\right)}Tn_{{\rm i}}u_{{\rm m}}\,{\rm d}S. \end{equation}

In other words, an ion–electron pair escaping the trap carries away energy $\alpha _{{\rm E}}T$ through the mirror (see Ryutov Reference Ryutov2005; Skovorodin Reference Skovorodin2019). The factor $\alpha _{{\rm E}}$ depends on the nature of the plasma outflow through the mirror; for the GDT facility in Budker INP, this parameter lies in the range of $6\div 8$ (see Soldatkina et al. Reference Soldatkina, Maximov, Prikhodko, Savkin, Skovorodin, Yakovlev and Bagryansky2020). Further, taking into account (4.8), the transverse energy flow is written in the form

(4.17)\begin{equation} \varPhi_{{\rm E}\perp}=\frac{5}{2}\frac{p}{n_{{\rm i}}}\varPhi_{{\rm i}\perp}-\varLambda_{{\rm E}\perp}\frac{{\rm d}T}{{\rm d}\psi}, \end{equation}

where

(4.18)\begin{equation} \varLambda_{{\rm E}\perp}=4{\rm \pi}^{2}\int_{\gamma_{{\perp}}\left(\psi\right)}\varkappa_{{\rm w}}Br^{2}\,{\rm d}l. \end{equation}

As a result, after taking the derivative of the energy balance (4.14) with respect to $\psi$, we arrive at the energy conservation equation for the warm plasma in the following form:

(4.19)\begin{equation} {-}\frac{{\rm d}}{{\rm d}\psi}\left(\frac{5}{2}\frac{p}{n_{{\rm i}}}\varLambda_{{\rm i}\perp} \frac{{\rm d}p}{{\rm d}\psi}\right)-\frac{{\rm d}}{{\rm d}\psi} \left(\varLambda_{{\rm E}\perp}\frac{{\rm d}T}{{\rm d}\psi}\right)+ \frac{2\alpha_{{\rm E}}u_{{\rm m}}}{B_{{\rm m}}}n_{{\rm i}}T=\frac{{\rm d}Q_{{\rm h}}}{{\rm d}\psi}.\end{equation}

4.4 Bubble core equilibrium

Strictly speaking, (4.3), (4.10) and (4.19) cannot be used to describe the equilibrium of the warm plasma in the bubble core, where $\psi \rightarrow 0$. In addition, the magnetized plasma approximation obviously ceases to be applicable there. For this reason, the equilibrium of the warm plasma inside the bubble core should be considered separately.

Since the magnetic field inside the bubble core is close to zero, the warm plasma transport in the core should increase significantly compared with the external transverse diffusion. Therefore, we assume the pressure, the temperature and the density of the warm plasma to be uniform inside the entire bubble core. In addition, we further suppose the warm plasma in the bubble core to be perfectly conducting, which results in the magnetic field there being constant and identically equal to zero

(4.20)\begin{equation} B\equiv0,\quad r< r_{0}\left(z\right). \end{equation}

This is possible only if there is no plasma current inside the core, i.e. the electric current of the hot ions is completely compensated by the inductive current of the warm plasma

(4.21)\begin{equation} J_{{\rm w}\theta}=-J_{{\rm h}\theta},\quad r< r_{0}\left(z\right). \end{equation}

It is shown by Chernoshtanov (Reference Chernoshtanov2020, Reference Chernoshtanov2022) that the so-called non-adiabatic loss should arise in the diamagnetic trap. In the case of a sufficiently collisional warm plasma with an isotropic distribution function, there are particles outside the absolute confinement region (2.2). These particles may reach the mirrors and escape the trap directly from the bubble core. The thermal equilibrium distribution of the warm ions in a diamagnetic trap can be approximately represented as follows (see Chernoshtanov Reference Chernoshtanov2020, Reference Chernoshtanov2022):

(4.22)\begin{equation} f_{{\rm i}}\left(\mathcal{E},\mathcal{P}\right)=n_{{\rm i}0}\left(\frac{m_{{\rm i}}}{2{\rm \pi} T_{0}}\right)^{{3}/{2}}{\rm e}^{-{\mathcal{E}}/{T_{0}}}\varTheta\left(a-\frac{\left|\mathcal{P}\right|}{\sqrt{2m_{{\rm i}} \mathcal{E}}}\right),\end{equation}

where $n_{{\rm i}0}$ and $T_{0}$ are the warm plasma density and temperature in the bubble core, $\varTheta (x)$ is the Heaviside step function, $a=\max [r_{0}(z)]=r_{0}(0)$ is the maximum radius of the core and $\mathcal {P}$ and $\mathcal {E}$ are the canonical angular momentum and the total energy, respectively, defined as in (2.1a,b). Such a distribution function corresponds to the collisional plasma when the characteristic Maxwellization time is much less than the confinement time. It can be obtained by direct averaging of the ‘ordinary’ Maxwellian distribution over the hypersurface of the constant $\mathcal {P}$ and $\mathcal {E}$ (see the definition of averaging (5.13) in § 5). The total axial loss from the core of the bubble is determined by the flow through the mirrors

(4.23)\begin{equation} \varPhi_{{\rm i}\parallel0}\simeq\frac{3}{4}n_{{\rm i}0}v_{T{\rm i}0}\frac{2{\rm \pi} a\rho_{{\rm i}0}}{\mathcal{R}},\end{equation}

where $v_{T{\rm i}0}=\sqrt {T_{0}/m_{{\rm i}}}$ is the warm ion thermal velocity, $\rho _{{\rm i}0}=v_{T{\rm i}0}m_{{\rm i}}c/eZB_{{\rm v}}$ is the characteristic warm ion Larmor radius in the vacuum field $B_{{\rm v}}$, $\mathcal {R}=B_{{\rm m}}/B_{{\rm v}}$ is the vacuum mirror ratio and $B_{{\rm m}}\simeq \mathrm {const}.$ is the mirror magnetic field. The expression (4.23) is obtained in the limit $a\gg \rho _{{\rm i}0}/\mathcal {R}$; for details, refer to Appendix A.

The non-adiabatic loss (4.23) is to be taken into account when setting the internal boundary condition for the mass conservation equation (4.10). Namely, the total transverse ion flow through the boundary of the bubble core is the total warm ion source inside the core $W_{{\rm i}0}$ minus the non-adiabatic loss through the mirrors $2\varPhi _{{\rm i}\parallel 0}$

(4.24)\begin{equation} \left.\left(-\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}\right)\right|_{\psi\rightarrow0}=W_{{\rm i}0}- 2\varPhi_{{\rm i}\parallel0}.\end{equation}

The corresponding axial energy loss is defined as follows:

(4.25)\begin{equation} \varPhi_{{\rm E}\parallel0}=\alpha_{{\rm E}0}T_{0}\varPhi_{{\rm i}\parallel0}, \end{equation}

where $\alpha _{{\rm E}0}$ is a coefficient similar to the previously introduced $\alpha _{{\rm E}}$ for gas-dynamic energy loss. Then, the internal boundary condition for the energy conservation equation (4.19) is

(4.26)\begin{equation} \left.\left(-\frac{5}{2}\frac{p}{n_{{\rm i}}}\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}- \varLambda_{{\rm E}\perp}\frac{{\rm d}T}{{\rm d}\psi}\right)\right|_{\psi\rightarrow0}=Q_{{\rm h}0}-2\alpha_{{\rm E}0}T_{0} \varPhi_{{\rm i}\parallel0},\end{equation}

where $Q_{{\rm h}0}$ is the total heat release from the hot ions inside the bubble core.

To complete the formulation of the boundary value problem, the external boundary conditions are to be set. As such, for example, the conditions

(4.27a,b)\begin{equation} \left.n_{{\rm i}}\right|_{r=a_{\mathrm{lim}}}=0,\quad\left.T\right|_{r=a_{\mathrm{lim}}}=0\end{equation}

can be chosen, which correspond to an external limiter located at the radius $a_{\mathrm {lim}}$. If the particle source is localized inside the bubble core, as we assume in the present paper, due to the large axial loss, the actual boundary of the warm plasma typically does not reach the limiter. This means that, in this case, the equilibrium appears to be weakly dependent on the external boundary conditions (4.27a,b).

5.1 Hot ion distribution function

As noted above, we consider the hot ion density to be small enough to neglect their collisions with each other. In addition, we neglect the angular scattering of the hot ions due to Coulomb collisionsFootnote ⁸ and take into account only the warm electron drag force

(5.7)\begin{equation} \boldsymbol{F}_{{\rm s}}\simeq-\nu_{\mathrm{ie}}m_{{\rm i}}\boldsymbol{v},\end{equation}

where $\boldsymbol {v}$ is the hot ion velocity,

(5.8)\begin{equation} \nu_{\mathrm{ie}}=1.6\times10^{-9}\mu_{{\rm i}}^{-1}Z^{3}\varLambda_{\mathrm{ie}}\frac{n_{{\rm i}}}{{\rm cm}^{-3}}\left(\frac{T}{{\rm eV}}\right)^{-{3}/{2}}{\rm s}^{-1},\end{equation}

is the inverse ion slowing down time (Trubnikov Reference Trubnikov1965), $\mu _{{\rm i}}=m_{{\rm i}}/m_{{\rm p}}$ is the ion mass $m_{{\rm i}}$ normalized to proton mass $m_{{\rm p}}$, $\varLambda _{\mathrm {ie}}$ is the ion–electron Coulomb logarithm. In what follows, we assume that the neutral beam is injected into the absolute confinement region (2.2); in this case, the loss of the hot ions can be ignored. Therefore, in the stationary case under consideration, the kinetic equation for the distribution function of the hot ions $f_{{\rm h}}$ has the following invariant form:

(5.9)\begin{equation} \boldsymbol{\nabla}_{\boldsymbol{X}}\boldsymbol{\cdot}\left(\dot{\boldsymbol{X}}f_{{\rm h}}\right)=g_{{\rm h}},\end{equation}

where $\boldsymbol {X}$ is the six-dimensional vector of generalized phase variables, $\dot {\boldsymbol {X}}$ is the corresponding phase velocity vector, $\boldsymbol {\nabla }_{\boldsymbol {X}}$ is the nabla differential operator in the $\boldsymbol {X}$-space and $g_{{\rm h}}=g_{{\rm h}}(\boldsymbol {X})$ is the phase density of the hot ion source intensity, the explicit form of which is determined by the injection. Specifying the variables as $\boldsymbol {X}=(r,\theta,z,P_{r},P_{\theta },P_{z})$, we obtain the equations of motion in the form

(5.10)\begin{equation} \left.\begin{gathered} \dot{r} =\frac{P_{r}}{m_{{\rm i}}}, \quad \dot{P}_{r} =\frac{P_{\theta}-\varPsi}{m_{{\rm i}}r^{2}} \left(\frac{P_{\theta}-\varPsi}{r}+\partial_{r}\varPsi\right)-\nu_{\mathrm{ie}}P_{r},\\ \dot{\theta} =\frac{P_{\theta}-\varPsi}{m_{{\rm i}}r^{2}}, \quad \dot{P}_{\theta} =- \nu_{\mathrm{ie}}\left(P_{\theta}-\varPsi\right),\\ \dot{z} =\frac{P_{z}}{m_{{\rm i}}}, \quad \dot{P}_{z} =\frac{P_{\theta}-\varPsi}{m_{{\rm i}}r^{2}}\partial_{z}\varPsi-\nu_{\mathrm{ie}}P_{z}. \end{gathered}\right\} \end{equation}

Since we assume axial symmetry, the distribution function $f_{{\rm h}}$ and the magnetic flux $\varPsi$ do not depend on the cyclic variable $\theta$.

As mentioned above, we consider that the dynamics of ions is globally chaotic, and the ergodization time is of the order of several periods of longitudinal oscillations of the ion $\tau _{{\rm b}}\sim L/\sqrt {\mathcal {E}/m_{{\rm i}}}$. Since the characteristic slowing down time $\nu _{\mathrm {ie}}^{-1}$ significantly exceeds the characteristic bounce time $\tau _{{\rm b}}$ (typically $\nu _{\mathrm {ie}}\tau _{{\rm b}}\sim 10^{-3}\div 10^{-4}$), the drag force (5.7) can be considered weak on time scales of the period of longitudinal oscillations $\tau _{{\rm b}}$. By means of the conventionally used Krylov–Bogoliubov–Mitropolsky averaging approach (Bogoliubov & Mitropolskii Reference Bogoliubov and Mitropolskii1961), the slow phase space diffusion associated with the drag force (5.7) can be separated from the fast longitudinal oscillations. Based on this, we replace the exact kinetic equation (5.9) by an averaged one, assuming that the averaging should be performed over the $\mathcal {E},\mathcal {P}=\textrm {const}.$ surface in the phase space.

Instead of $(P_{r},P_{\theta },P_{z})$, it is convenient to use the new variables $(\mathcal {E},\mathcal {P},\alpha )$, where $\alpha \in [0,2{\rm \pi} )$ is defined by

(5.11a,b)\begin{equation} P_{r}=\sqrt{2m_{{\rm i}}\left(\mathcal{E}-\varphi_{\mathrm{eff}}\right)}\sin\alpha,\quad P_{z}=\sqrt{2m_{{\rm i}}\left(\mathcal{E}-\varphi_{\mathrm{eff}}\right)}\cos\alpha. \end{equation}

Hence, we have

(5.12)\begin{equation} \frac{\partial\left(P_{r},P_{\theta},P_{z}\right)}{\partial\left(\mathcal{E},\mathcal{P},\alpha\right)}=m_{{\rm i}}, \end{equation}

and we should also keep in mind the additional condition $\mathcal {E}>\varphi _{\mathrm {eff}}$, which defines the region of permissible values of the new variables. Then averaging of some quantity $\mathcal {Q}$ over the surface $\mathcal {E},\mathcal {P}=\textrm {const}.$ is performed as follows:

(5.13)\begin{equation} \left\langle \mathcal{Q}\right\rangle _{\varGamma}\overset{\mathrm{def}}=\varGamma^{-1}\idotsint_{\mathcal{E}>\varphi_{\mathrm{eff}}} \mathcal{Q}m_{\mathrm{i}}\,{\rm d}r\,{\rm d}\theta \,{\rm d}z\,{\rm d}\alpha,\end{equation}

where

(5.14)\begin{equation} \varGamma=\idotsint _{\mathcal{E}>\varphi_{\mathrm{eff}}}m_{{\rm i}}\,{\rm d}r\,{\rm d}\theta \,{\rm d}z\,{\rm d}\alpha=4{\rm \pi}^{2}m_{{\rm i}}\iint _{\mathcal{E}>\varphi_{\mathrm{eff}}}\,{\rm d}r\,{\rm d}z, \end{equation}

is the phase volume of the hypersurface $\mathcal {E},\mathcal {P}=\mathrm {const}$. Finally, at the leading order with respect to $\nu _{\mathrm {ie}}\tau _{{\rm b}}\ll 1$, when the distribution function $f_{{\rm h}}$ approximately depends only on $\mathcal {E}$ and $\mathcal {P}$, averaging the perturbed system yields

(5.15)\begin{equation} \left.\begin{gathered} \varGamma^{-1}\partial_{\mathcal{E}}\left(\varGamma\left\langle \dot{\mathcal{E}}\right\rangle _{\varGamma}f_{{\rm h}} \right)+\varGamma^{-1}\partial_{\mathcal{P}}\left(\varGamma\left\langle \dot{\mathcal{P}}\right\rangle _{\varGamma}f_{{\rm h}} \right)=\left\langle g_{{\rm h}}\right\rangle _{\varGamma},\\ \left\langle \dot{\mathcal{E}}\right\rangle _{\varGamma}=-2\left\langle \nu_{\mathrm{ie}}\right\rangle _{\varGamma} \mathcal{E},\quad\left\langle \dot{\mathcal{P}}\right\rangle _{\varGamma}=-\left\langle \nu_{\mathrm{ie}}\right\rangle _{\varGamma}\mathcal{P}+\left\langle \nu_{\mathrm{ie}}\varPsi\right\rangle _{\varGamma}.\\ \left\langle \dot{r}\right\rangle _{\varGamma}=0,\quad\left\langle \dot{z}\right\rangle _{\varGamma}=0,\quad\left\langle \dot{\alpha}\right\rangle _{\varGamma}=0. \end{gathered}\right\}\end{equation}

Solving the kinetic equation (5.15) one can find the averaged distribution function of the hot ions $f_{{\rm h}}\simeq f_{{\rm h}}(\mathcal {E},\mathcal {P})$. This in turn enables the density of some quantity $\mathcal {Q}$ to be obtained

(5.16)\begin{equation} \left\langle \mathcal{Q}\right\rangle _{{\rm h}}\overset{\mathrm{def}}=\frac{1}{r}\iiint\mathcal{Q}f_{{\rm h}}\,{\rm d}P_{r}\,{\rm d}P_{\theta}\, {\rm d}P_{z}=\frac{m_{{\rm i}}}{r}\int_{0}^{2{\rm \pi}}{\rm d}\alpha\int _{-\infty}^{+\infty}{\rm d}\mathcal{P}\int_{\varphi_{\mathrm{eff}}}^{+\infty}{\rm d}\mathcal{E} \mathcal{Q}f_{{\rm h}}. \end{equation}

In particular, the azimuthal electric current density of the hot ions is given by

(5.17)\begin{equation} J_{{\rm h}\theta}=\left\langle Zev_{\theta}\right\rangle _{{\rm h}}=\frac{2{\rm \pi} Ze}{r^{2}}\int _{-\infty}^{+\infty}{\rm d}\mathcal{P}\int_{\varphi_{\mathrm{eff}}}^{+\infty}{\rm d}\mathcal{E}\left(\mathcal{P}-\varPsi\right)f_{{\rm h}}, \end{equation}

the power density of heating the warm plasma by the hot ions is defined as

(5.18)\begin{equation} q_{{\rm h}}=\left\langle 2\left\langle \nu_{\mathrm{ie}}\right\rangle _{\varGamma}\mathcal{E}\right\rangle _{{\rm h}}=\frac{4{\rm \pi} m_{{\rm i}}}{r}\int _{-\infty}^{+\infty}{\rm d}\mathcal{P}\int _{\varphi_{\mathrm{eff}}}^{+\infty}{\rm d}\mathcal{E}\left\langle \nu_{\mathrm{ie}}\right\rangle _{\varGamma}\mathcal{E}f_{{\rm h}}.\end{equation}

6.1 Magnetic field distribution

The equilibrium equation for the magnetic field (3.4) can only be considered in the central cylindrical region, and the magnetic flux in the mirrors is considered to be given: $\psi _{{\rm m}}=B_{{\rm m}}{\rm \pi} r^{2}$. Then, after substituting the current density of the warm plasma (4.3), (3.4) predictably reduces to the paraxial equilibrium with the Lorentz force from the hot ions

(6.1)\begin{gather} \frac{{\rm d}}{{\rm d}\psi}\left(\frac{B^{2}}{8{\rm \pi}}+p\right)=-\frac{J_{{\rm h}\theta}}{2{\rm \pi} rc}, \end{gather}

(6.2)\begin{gather}\frac{{\rm d}r}{{\rm d}\psi}=\frac{1}{2{\rm \pi} r}\frac{1}{B}. \end{gather}

In this equations we use the flux coordinate $\psi$, so $r=r(\psi )$ is the inverse function of $\psi =\psi (r)$ and is essentially the radius of the magnetic surface that corresponds to the magnetic flux $\psi$. In addition, the magnetic field $B$ and the current density of the hot ions $J_{{\rm h}\theta }$ should also be considered here as functions of $\psi$. The corresponding boundary conditions (3.7), (3.8) and (3.9) can be reduced to

(6.3a,b)\begin{gather} \left.B\right|_{\psi=0}=0,\quad\left.r\right|_{\psi=0}=a, \end{gather}

(6.4)\begin{gather}\left.B\right|_{\psi\rightarrow+\infty}=B_{{\rm v}}, \end{gather}

where $B_{{\rm v}}=\mathrm {const}.$ is the vacuum magnetic field at the central section of the trap.

6.2 Warm plasma equilibrium

Warm plasma equilibrium equations (4.10) and (4.19) mainly remain the same except for the transport coefficients $\varLambda _{{\rm i}\perp }$ and $\varLambda _{{\rm E}\perp }$, which are reduced to a simpler form in the cylindrical bubble approximation

(6.5)\begin{gather} -\frac{{\rm d}}{{\rm d}\psi}\left(\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}\right)+ \frac{2u_{{\rm m}}}{B_{{\rm m}}}n_{{\rm i}}=0, \end{gather}

(6.6)\begin{gather}-\frac{{\rm d}}{{\rm d}\psi}\left(\frac{5}{2}\frac{p}{n_{{\rm i}}}\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}\right)- \frac{{\rm d}}{{\rm d}\psi}\left(\varLambda_{{\rm E}\perp}\frac{{\rm d}T}{{\rm d}\psi}\right)+ \frac{2\alpha_{{\rm E}}u_{{\rm m}}}{B_{{\rm m}}}n_{{\rm i}}T=0, \end{gather}

(6.7a,b)\begin{gather}\varLambda_{{\rm i}\perp}\simeq4{\rm \pi}^{2}c^{2}n_{{\rm i}}\frac{r^{2}L}{\sigma _{{\rm w}}B}, \quad\varLambda_{{\rm E}\perp}\simeq4{\rm \pi}^{2}\varkappa_{{\rm w}}Br^{2}L. \end{gather}

Here, we also take into account that, due to the hot ions mainly slowing down inside the bubble core, the external hot ion heating power release is negligible, i.e. $Q_{{\rm h}}\equiv 0$ and $W_{{\rm i}}\equiv 0$ for $\psi >0$. The boundary conditions (4.24), (4.26) and (4.27a,b), in turn, remain exactly the same, except that $Q_{{\rm h}0}$ and $W_{{\rm i}0}$ now have the meaning of the total absorbed injection power and the total warm plasma source, respectively.

6.3 Hot ion equilibrium

Consider the injection of a monoenergetic neutral beam with the injection energy $\mathcal {E}_{\mathrm {NB}}$ and the total absorbed injection power $Q_{{\rm h}0}$. The injection is carried out at the angle $\xi _{\mathrm {NB}}\in (0,{\rm \pi} /2)$ to the axis of the trap, and the distance from the beam to the trap axis is finite and equal to $r_{\mathrm {NB}}< a$. In this case, the source in the kinetic equation (5.15) has the formFootnote ⁹

(6.8)\begin{equation} \left\langle g_{{\rm h}}\right\rangle _{\varGamma}=\frac{1}{\varGamma}\frac{Q_{{\rm h}0}}{\mathcal{E}_{\mathrm{NB}}} \delta\left(\mathcal{E}-\mathcal{E}_{\mathrm{NB}}\right)\delta\left(\mathcal{P}-\mathcal{P}_{\mathrm{NB}}\right),\end{equation}

where $\delta (x)$ is the Dirac delta function

(6.9)\begin{equation} \mathcal{P}_{\mathrm{NB}}=-\sqrt{2m_{{\rm i}}\mathcal{E}_{\mathrm{NB}}}r_{\mathrm{NB}}\sin\xi_{\mathrm{NB}}<0,\end{equation}

is the angular momentum of the injected ions. Due to the beam slowing down mainly in the core of the bubble, averaging in (5.15) yields

(6.10a,b)\begin{equation} \left\langle \nu_{\mathrm{ie}}\right\rangle _{\varGamma}\simeq\nu_{\mathrm{ie}0},\quad\left\langle \nu_{\mathrm{ie}}\varPsi\right\rangle_{\varGamma}\simeq0, \end{equation}

where $\nu _{\mathrm {ie}0}=\nu _{\mathrm {ie}}|_{T=T_{0},n_{{\rm i}}=n_{{\rm i}0}}$. The resulting kinetic equation

(6.11)\begin{equation} {-}2\nu_{\mathrm{ie}0}\partial_{\mathcal{E}}\left(\mathcal{E}\varGamma f_{{\rm h}}\right)-\nu_{\mathrm{ie}0}\partial_{\mathcal{P}}\left(\mathcal{P}\varGamma f_{{\rm h}}\right)\simeq\frac{Q_{{\rm h}0}}{\mathcal{E}_{\mathrm{NB}}}\delta \left(\mathcal{E}-\mathcal{E}_{\mathrm{NB}}\right)\delta\left(\mathcal{P}-\mathcal{P}_{\mathrm{NB}}\right), \end{equation}

is satisfied by

(6.12a,b)\begin{equation} f_{{\rm h}}\simeq\frac{1}{2\nu_{\mathrm{ie}0}}\frac{1}{\varGamma}\frac{Q_{{\rm h}0}}{\mathcal{E}_{\mathrm{NB}}} \frac{1}{\mathcal{E}^{{3}/{2}}}\delta\left(\frac{\mathcal{P}}{\sqrt{\mathcal{E}}}- \frac{\mathcal{P}_{\mathrm{NB}}}{\sqrt{\mathcal{E}_{\mathrm{NB}}}}\right),\quad \mathcal{E}<\mathcal{E}_{\mathrm{NB}},\end{equation}

where the phase volume $\varGamma$ is approximately equal to

(6.13)\begin{equation} \varGamma=4{\rm \pi}^{2}m_{{\rm i}}L\int _{\mathcal{E}>\varphi_{\mathrm{eff}}}\,{\rm d}r, \end{equation}

and $\varphi _{\mathrm {eff}}$ is the effective potential (5.4) in the central plane.

The remaining integral in $\varGamma$ can be represented in a simpler form. To do this, we further assume that, first, there is no reversed field and, second, the magnetic field either increases with the radius or decreases not faster than $r^{-1}$. In other words, the following conditions are met:

(6.14a,b)\begin{equation} B\geq0,\quad\frac{r}{B}\frac{{\rm d}B}{{\rm d}r}\geq -1. \end{equation}

It can be shown that, in this case, the effective potential $\varphi _{\mathrm {eff}}$ has only one minimum. Consequently, the equation $\mathcal {E}=\varphi _{\mathrm {eff}}$ has no more than two roots: $r_{\min }=r_{\min }(\mathcal {E},\mathcal {P})$ and $r_{\max }=r_{\max }(\mathcal {E},\mathcal {P})$, which are the radial boundaries of the integration domain $\mathcal {E}>\varphi _{\mathrm {eff}}$. These roots are essentially the minimum and maximum distance from the axis available to a hot ion with energy $\mathcal {E}$ and angular momentum $\mathcal {P}$. In what follows, $r_{\min }$ and $r_{\max }$ are shortly referred to as the minimum radius and the maximum radius, respectively, and we also define the corresponding minimum and maximum radii for the injected ions: $\bar {r}_{\min }=r_{\min }(\mathcal {E}_{\mathrm {NB}},\mathcal {P}_{\mathrm {NB}})$ and $\bar {r}_{\max }=r_{\max }(\mathcal {E}_{\mathrm {NB}},\mathcal {P}_{\mathrm {NB}})$. Eventually, the phase volume can be represented in the following form:

(6.15)\begin{equation} \varGamma=4{\rm \pi}^{2}m_{{\rm i}}L\left(r_{\max}-r_{\min}\right). \end{equation}

It is worth noting that the minimum radius for the ions passing through the interior of the bubble core can be found explicitly. Indeed, since $\varPsi \equiv 0$ in the bubble core, for $r_{\min }< a$ we have

(6.16)\begin{equation} \mathcal{E}=\varphi_{\mathrm{eff}}\left(r_{\min}\right)=\frac{\left|\mathcal{P}\right|^{2}}{2m_{{\rm i}} r_{\min}^{2}},\quad\Rightarrow\quad r_{\min}=\frac{\left|\mathcal{P}\right|}{\sqrt{2m_{{\rm i}}\mathcal{E}}}. \end{equation}

Note that the solution (6.12a,b) also implies that, in the approximation considered, the hot ion energy and angular momentum are related as follows: $\mathcal {P}\sqrt {\mathcal {E}_{\mathrm {NB}}}=\mathcal {P}_{\mathrm {NB}}\sqrt {\mathcal {E}}$. This means that the minimum radius for the ions injected in the bubble core ($r_{\mathrm {NB}}< a$) remains constant

(6.17)\begin{equation} r_{\min}=\bar{r}_{\min}=\frac{\left|\mathcal{P}_{\mathrm{NB}}\right|}{\sqrt{2m_{{\rm i}}\mathcal{E}_{\mathrm{NB}}}}. \end{equation}

Given (6.9), the minimum radius can also be expressed in terms of the injection radius $r_{\mathrm {NB}}$ and injection angle $\xi _{\mathrm {NB}}$: $\bar {r}_{\min }=r_{\mathrm {NB}}\sin \xi _{\mathrm {NB}}$. In addition, it follows that the sign of the angular momentum does not change during the ion slowing down, and the ions initially injected into the absolute confinement region (2.2) always stay in it

(6.18)\begin{equation} {-}\mathcal{R}\varOmega_{{\rm i}}\mathcal{P}_{\mathrm{NB}}>\mathcal{E}_{\mathrm{NB}},\quad\Rightarrow \quad-\mathcal{R}\varOmega_{{\rm i}}\mathcal{P}>\sqrt{\mathcal{E}\mathcal{E}_{\mathrm{NB}}}\geq \mathcal{E}. \end{equation}

Therefore, the axial loss of the hot ions can indeed be neglected.Footnote ¹⁰

Eventually, using the resulting distribution function (6.12a,b) and the definition (5.17), we obtain the azimuthal diamagnetic current density of the hot ions

(6.19)\begin{align} \frac{J_{{\rm h}\theta}}{2{\rm \pi} rc}& =-\frac{\bar{r}_{\min}}{a}\frac{\varPi_{{\rm h}}}{\psi_{{\rm h}}} \frac{a^{3}}{r^{3}}\varTheta\left(\frac{r-\bar{r}_{\min}}{a-\bar{r}_{\min}}- \frac{\psi}{\psi_{{\rm h}}}\right)\nonumber\\ & \quad \times\int_{(({a-\bar{r}_{\min}})/({r-\bar{r}_{\min}}))({\psi}/{\psi_{{\rm h}}})}^{1}\frac{a-\bar{r}_{\min}}{ r_{\max}^{*}\left(\eta\right)-\bar{r}_{\min}}\left(1+\frac{a-\bar{r}_{\min}}{\bar{r}_{\min}} \frac{\psi}{\psi_{{\rm h}}}\frac{1}{\eta}\right){\rm d}\eta, \end{align}

where $\psi _{{\rm h}}=2{\rm \pi} B_{{\rm v}}(a-\bar {r}_{\min })\rho _{\mathrm {NB}}$ is the characteristic magnetic flux induced by the hot ion current, $\rho _{\mathrm {NB}}=\sqrt {2m_{{\rm i}}\mathcal {E}_{\mathrm {NB}}}c/ZeB_{{\rm v}}$ is the characteristic Larmor radius of injected ions, the coefficient

(6.20)\begin{equation} \varPi_{{\rm h}}=\frac{Q_{{\rm h}0}}{\nu_{\mathrm{ie}0}{\rm \pi} a^{2}L}, \end{equation}

is the characteristic energy density of the hot ions and

(6.21)\begin{equation} r_{\max}^{*}\left(\eta\right)\overset{\mathrm{def}}=r_{\max}\left(\mathcal{E}_{\mathrm{NB}}\eta^{2},\mathcal{P}_{\mathrm{NB}}\eta\right), \end{equation}

is the maximum radius of the hot ions with the distribution (6.12a,b), which is to be found from the equation $\mathcal {E}=\varphi _{\mathrm {eff}}$ reduced to explicit form

(6.22)\begin{equation} r_{\max}^{*}\left(\eta\right)=\left.r\left(\psi\right)\right|_{\psi=\psi_{{\rm h}}\eta (({r_{\max}^{*}\left(\eta\right)-\bar{r}_{\min}})/({a-\bar{r}_{\min}}))}.\end{equation}

In particular, the outer boundary is $\bar {r}_{\max }=r_{\max }^{*}(1)$.

7.1 Equilibrium magnetic field distribution

The magnetic field equilibrium is determined by (6.1) and (6.2) with the boundary conditions (6.3a,b) and (6.4). When solving these equations, we consider the pressure profile of the warm plasma $p=p(\psi )$ given, implying that it can be found from the corresponding equilibrium equations (6.5) and (6.6). At the same time, the diamagnetic current of the hot ions is determined by the expression (6.19).

It is further convenient to use the following dimensionless quantities:

(7.2a-d)\begin{equation} \phi=\frac{\psi}{\psi_{{\rm h}}},\quad x=\frac{r}{a},\quad R=\frac{B}{B_{{\rm v}}},\quad\beta_{{\rm w}}=\frac{8{\rm \pi} p}{B_{{\rm v}}^{2}}, \end{equation}

where $\beta _{{\rm w}}$ is the warm plasma pressure normalized to the energy density of the vacuum magnetic field $B_{{\rm v}}^{2}/8{\rm \pi}$, also referred to simply as the relative pressure of the warm plasma. According to this, we also define the dimensionless minimum and maximum radii as follows:

(7.3a,b)\begin{equation} \bar{x}_{\min}=\frac{\bar{r}_{\min}}{a},\quad x_{\max}^{*}\left(\eta\right)=\frac{r_{\max}^{*}\left(\eta\right)}{a}, \end{equation}

where the latter is determined by the equation

(7.4)\begin{equation} x_{\max}^{*}\left(\eta\right)=\left.x\left(\phi\right)\right|_{\phi=\eta (({x_{\max}^{*}\left(\eta\right)-\bar{x}_{\min}})/{(1-\bar{x}_{\min}}))},\end{equation}

which results from the corresponding normalization of (6.22). Then, (6.1) and (6.2) with substituted hot ion current (6.19) take the following dimensionless form:

(7.5)\begin{gather} \frac{{\rm d}}{{\rm d}\phi}\left(R^{2}+\beta_{{\rm w}}\right)\nonumber\\ =\frac{\varLambda_{{\rm h}}}{x^{3}}\varTheta\left(\frac{x-\bar{x}_{\min}}{1-\bar{x}_{\min}}-\phi\right) \int_{(({1-\bar{x}_{\min}})/({x-\bar{x}_{\min}}))\phi}^{1} \frac{1-\bar{x}_{\min}}{x_{\max}^{*}\left(\eta\right)-\bar{x}_{\min}}\left(1+ \frac{1-\bar{x}_{\min}}{\bar{x}_{\min}}\frac{\phi}{\eta}\right){\rm d}\eta, \end{gather}

(7.6)\begin{gather} \frac{{\rm d}x}{{\rm d}\phi}=\frac{\epsilon}{x}\frac{1}{R},\end{gather}

where

(7.7a,b)\begin{equation} \varLambda_{{\rm h}}=\bar{x}_{\min}\frac{8{\rm \pi}\varPi_{{\rm h}}}{B_{{\rm v}}^{2}},\quad \epsilon=\left(1-\bar{x}_{\min}\right)\frac{\rho_{\mathrm{NB}}}{a}. \end{equation}

Normalizing the boundary conditions (6.3a,b) and (6.4) yields

(7.8a,b)\begin{gather} \left.R\right|_{\phi=0}=0,\quad\left.x\right|_{\phi=0}=1, \end{gather}

(7.9)\begin{gather}\left.R\right|_{\phi\rightarrow+\infty}=1. \end{gather}

As mentioned above, the hot ion transition layer thickness $\lambda _{{\rm h}}$ appears to be of the order of the Larmor radius of the injected ions $\rho _{\mathrm {NB}}$, which means that, in the thin transition layer limit before us $\lambda _{{\rm h}}\ll a$, we should also expect the ratio $\rho _{\mathrm {NB}}/a$ to be small. Given this, it seems appropriate to apply the method of successive approximations by considering the coefficient $\epsilon \propto \rho _{\mathrm {NB}}/a$ on the right-hand side of (7.6) as a small expansion parameter. Namely, we further assume that the solution of the system (7.5) and (7.6) can be represented in the form of an asymptotic expansion.

For a given magnetic field profile $R=R(\phi )$, by integrating (7.6) and taking into account the boundary condition (7.8a,b), the radius as a function of the magnetic flux can be explicitly expressed

(7.10)\begin{equation} x\left(\phi\right)=\sqrt{1+2\epsilon\int_{0}^{\phi}\frac{{\rm d}\phi^{\prime}}{R\left(\phi^{\prime}\right)}}.\end{equation}

At the leading order of the approximation $\epsilon \ll 1$, when small corrections are completely neglected, the expression (7.10) reduces to

(7.11)\begin{equation} x\left(\phi\right)\simeq1+o\left(\epsilon^{0}\right).\end{equation}

Substituting (7.11) into (7.4) yields the corresponding approximation for the maximum radius

(7.12)\begin{equation} x_{\max}^{*}\left(\eta\right)\simeq1+o\left(\epsilon^{0}\right).\end{equation}

The obtained formulas (7.11) and (7.4) correspond to $r\simeq a$. In other words, on the scale of the transition layer, the radius $r$ can be considered almost constant and equal to the bubble radius $a$, while the magnetic flux $\psi$ varies greatly. At a qualitative level, this is associated with the high density of the diamagnetic current inside the layer.

Taking into account (7.11) and (7.12), the integral in the right-hand side of (7.5) is evaluated explicitly

(7.13)\begin{align} \frac{{\rm d}}{{\rm d}\phi}\left(R^{2}+\beta_{{\rm w}}\right)& \simeq \varLambda_{{\rm h}}\varTheta\left(1-\phi\right)\int _{\phi}^{1}\left[1+\left(\bar{x}_{\min}^{-1}-1\right) \frac{\phi}{\eta}\right]{\rm d}\eta+o\left(\epsilon^{0}\right)\nonumber\\ & =\varLambda_{{\rm h}}\left[1-\phi-\left(\bar{x}_{\min}^{-1}-1\right)\phi\ln\phi\right] \varTheta\left(1-\phi\right)+o\left(\epsilon^{0}\right). \end{align}

As can be seen, the outer boundary of the hot ions at this order of approximation corresponds to $\phi =1$. Given the boundary conditions (7.8a,b) and (7.9), the resulting equation yields

(7.14)\begin{equation} R\left(\phi\right)\simeq R^{\left(0\right)}\left(\phi\right)+o\left(\epsilon^{0}\right),\end{equation}

where

(7.15)\begin{align} R^{\left(0\right)}\left(\phi\right)=\begin{cases} \sqrt{\beta_{{\rm w}0}-\beta_{{\rm w}}\left(\phi\right)+\varLambda_{{\rm h}}\phi\left[1- \dfrac{\phi}{2}-\left(\bar{x}_{\min}^{-1}-1\right)\dfrac{\phi}{2}\left(\ln\phi-\dfrac{1}{2}\right)\right]}, & \phi<1,\\[10pt] {}\sqrt{1-\beta_{{\rm w}}\left(\phi\right)}, & \phi\geq1, \end{cases} \end{align}

and $\beta _{{\rm w}0}=\beta _{{\rm w}}(0)=8{\rm \pi} p_{0}/B_{{\rm v}}^{2}$ is the relative pressure of the warm plasma inside the bubble core. In the absence of surface currents, the solution (7.14) should be continuous at the boundary $\phi =1$, from which we also obtain the following condition:

(7.16)\begin{equation} 1\simeq\beta_{{\rm w}0}+\beta_{{\rm h}0}+o\left(\epsilon^{0}\right).\end{equation}

Here, the quantity

(7.17)\begin{equation} \beta_{{\rm h}0}=\frac{1}{4}\left(\bar{x}_{\min}^{-1}+1\right)\varLambda_{{\rm h}}= \left(1+\bar{x}_{\min}\right)\frac{2{\rm \pi}\varPi_{{\rm h}}}{B_{{\rm v}}^{2}}, \end{equation}

can be interpreted as the characteristic relative energy density of the hot ions in the bubble. The formula (7.16) essentially expresses the balance between the plasma pressure inside the bubble and the pressure of the external magnetic field.

By applying the successive approximations to find higher orders of the expansion, the solution can be further refined. In particular, taking into account first-order corrections in the expression (7.10) yields

(7.18)\begin{equation} x\left(\phi\right)\simeq1+\epsilon\int_{0}^{\phi}\frac{{\rm d}\phi^{\prime}}{R^{\left(0\right)}\left(\phi^{\prime}\right)}+o\left(\epsilon^{1}\right). \end{equation}

Using the obtained dependence, the corresponding correction to the maximum radius is also found from (7.4)

(7.19)\begin{equation} x_{\max}^{*}\left(\eta\right)\simeq1+\epsilon\int_{0}^{\eta}\frac{{\rm d}\phi}{R^{\left(0\right)} \left(\phi\right)}+o\left(\epsilon^{1}\right). \end{equation}

For $\eta =1$ this formula determines the thickness of the hot ion transition layer

(7.20)\begin{equation} \frac{\lambda_{{\rm h}}}{a}=x_{\max}^{*}\left(1\right)-1\simeq\epsilon\int_{0}^{1}\frac{{\rm d}\phi}{R^{ (0)}\left(\phi\right)}+o\left(\epsilon^{1}\right).\end{equation}

Further, in this paper, we consider the equilibrium configuration of the magnetic field to be approximately defined by (7.11) and (7.14) at $\epsilon \rightarrow 0$. In other words, we take into account only the leading order of the asymptotic expansion. The applicability criterion for this approximation is determined by the limit

(7.21)\begin{equation} \frac{1-x_{\max}^{*}\left(1\right)}{1-\bar{x}_{\min}}=\frac{\lambda_{{\rm h}}}{a-\bar{r}_{\min}}\ll1, \end{equation}

which is actually assumed in the derivation of (7.14). After substituting (7.20), we arrive at the following condition:

(7.22)\begin{equation} \int _{0}^{1}\frac{{\rm d}\phi}{R^{\left(0\right)}\left(\phi\right)}\ll\frac{a}{\rho_{\mathrm{NB}}}. \end{equation}

An upper bound for the integral involved can be obtained as follows:

(7.23)\begin{equation} \int_{0}^{1}\frac{{\rm d}\phi}{R^{\left(0\right)}\left(\phi\right)}\leq \frac{1}{\sqrt{1-\beta_{{\rm w}0}}} \mathcal{W}\left(\frac{a}{\bar{r}_{\min}}\right), \end{equation}

where

(7.24)\begin{equation} \mathcal{W}\left(z\right)\overset{\mathrm{def}}=\frac{\sqrt{z+1}}{2}\int_{0}^{1} \dfrac{{\rm d}\phi}{\sqrt{\phi\left[1-\dfrac{\phi}{2}-\left(z-1\right) \dfrac{\phi}{2}\left(\ln\phi-\dfrac{1}{2}\right)\right]}}.\end{equation}

As can be seen from figure 2, for all reasonable $a/\bar {r}_{\min }$ the function $\mathcal {W}$ is of the order of $2\div 3$.

Figure 2.

Function $\mathcal {W}=\mathcal {W}(z)$ defined by the expression (7.24).

Before moving on to considering the equilibrium of the warm plasma, it is also useful to examine the expression (7.16) in greater detail. For given parameters of the warm plasma inside the bubble core: density $n_{{\rm i}0}$ and temperature $T_{0}$, which are to be found from the solution of the equilibrium equations (6.5) and (6.6), the expression (7.16) defines the relation between the radius of the core $a$ and the fixed external parameters: absorbed heating power $Q_{{\rm h}0}$, vacuum magnetic field $B_{{\rm v}}$ and geometric parameters $L$ and $\bar {r}_{\min }$. After expressing the warm plasma density $n_{{\rm i}0}$ from the equation of state (4.2), the formula (7.16) is reduced to a cubic equation for the normalized core radius $k=a/\bar {r}_{\min }=\bar {x}_{\min }^{-1}$

(7.25)\begin{equation} k^{3}-\left(4\tau_{{\rm s}}\nu_{\mathrm{ie}0}\right)^{-1}\left(k+1\right)\simeq0,\end{equation}

where

(7.26)\begin{equation} \tau_{{\rm s}}=\beta_{{\rm h}0}\frac{B_{{\rm v}}^{2}}{8{\rm \pi}}\frac{{\rm \pi}\bar{r}_{\min}^{2}L}{Q_{{\rm h}0}}. \end{equation}

Provided that the core radius $a$ exceeds the minimum radius $\bar {r}_{\min }$, i.e. $k\geq 1$, this equation has exactly one real root for a fixed $\beta _{{\rm w}0}\in [0,1]$

(7.27)\begin{equation} k = \begin{cases} \dfrac{3}{\mathfrak{D}}\cos\left(\dfrac{1}{3}\arccos\mathfrak{D}\right), & 0<\mathfrak{D}\leq1,\\[10pt] \dfrac{3}{\mathfrak{D}}\cosh\left(\dfrac{1}{3}\mathrm{arccosh}\,\mathfrak{D}\right), & 1<\mathfrak{D}, \end{cases} \end{equation}

where $\mathfrak {D}=\sqrt {27\tau _{{\rm s}}\nu _{\mathrm {ie}0}}$ and $\mathrm {arccosh}\,x=\ln (x+\sqrt {x^{2}-1})$. This means that, for a given temperature in the core $T_{0}$, this root specifies the functional dependency of the core radius on the pressure in the core: $k=k(\beta _{{\rm w}0})$.

7.2 Warm plasma equilibrium

The equilibrium of the warm plasma is described by (6.5) and (6.6) with the boundary conditions (4.24), (4.26) and (4.27a,b). In what follows, the electrical and thermal conductivity of the warm plasma are considered classical (Braginskii Reference Braginskii1965)

(7.28)\begin{equation} \left.\begin{gathered} \sigma _{{\rm w}}=\frac{Z^{2}e^{2}n_{{\rm i}}}{m_{{\rm i}}\nu_{\mathrm{ie}}}=8.7\times10^{13}Z^{-1}\varLambda_{\mathrm{ie}}^{-1} \left(\frac{T}{\mathrm{eV}}\right)^{{3}/{2}}\ \mathrm{s}^{-1},\\ \varkappa_{{\rm w}}=10^{-3}\left(\mu_{{\rm i}}^{{1}/{2}}+0.077Z\right)Z^{2}\varLambda_{\mathrm{ie}} \left(\frac{B}{\mathrm{G}}\right)^{-2}\left(\frac{n_{{\rm i}}}{{\rm cm}^{-3}}\right)^{2}\left(\frac{T}{{\rm eV}}\right)^{-{1}/{2}}{\rm cm}^{-1}\ {\rm s}^{-1}. \end{gathered}\right\} \end{equation}

In addition, the warm plasma flow velocity in the mirror is assumed to be

(7.29)\begin{equation} u_{{\rm m}}=\sqrt{\frac{2}{\rm \pi}}v_{T{\rm i}}=\sqrt{\frac{2T}{{\rm \pi} m_{{\rm i}}}}, \end{equation}

which corresponds to the gas-dynamic loss in the case of short mirrors and a filled loss cone (Ivanov & Prikhodko Reference Ivanov and Prikhodko2017). In what follows, we also consider $r\simeq a$, assuming the leading order of the thin transition layer approximation $\lambda _{{\rm w}}\ll a$.

Let us show that in this case there exists an equilibrium of the warm plasma such that the temperature can be considered slowly varying on the scale of the warm plasma inhomogeneity

(7.30)\begin{equation} \frac{n_{{\rm i}}}{T}\frac{{\rm d}T}{{\rm d}n_{{\rm i}}}\simeq\frac{p}{T}\frac{{\rm d}T}{{\rm d}p}=\gamma_{T}\ll1.\end{equation}

Subtracting (6.5) multiplied by $\alpha _{{\rm E}}T$ from (6.6) yields

(7.31)\begin{equation} \alpha_{{\rm E}}T{\rm d}\left(\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}\right)- {\rm d}\left(\frac{5\left(1+Z\right)}{2}T\varLambda_{{\rm i}\perp}\frac{{\rm d}p}{{\rm d}\psi}\right)- {\rm d}\left(\varLambda_{{\rm E}\perp}\frac{{\rm d}T}{{\rm d}\psi}\right)=0. \end{equation}

Taking into account the condition (7.30), the left-hand side of the resulting expression is approximately reduced to an exact differential

(7.32)\begin{equation} {\rm d}\left(\left(\alpha_{{\rm E}}-\frac{5\left(1+Z\right)}{2}\right)T\varLambda_{{\rm i}\perp} \frac{{\rm d}p}{{\rm d}\psi}\right)-{\rm d}\left(\varLambda_{{\rm E}\perp}\frac{{\rm d}T}{{\rm d}\psi}\right)\simeq0, \end{equation}

which further allows the equation to be explicitly integrated. The constant of integration should be set equal to zero, since at the boundary of the warm plasma, where $n_{{\rm i}}\rightarrow 0$, the transverse fluxes of mass and energy should vanish. Finally, we find the slope factor

(7.33)\begin{align} \gamma_{T}=\frac{p}{T}\frac{{\rm d}T}{{\rm d}p}\simeq\left(\alpha_{{\rm E}}- \frac{5\left(1+Z\right)}{2}\right)\frac{\varLambda_{{\rm i}\perp}}{\varLambda_{{\rm E}\perp}} p\simeq\frac{1}{60}\left(\alpha_{{\rm E}}- \frac{5\left(1+Z\right)}{2}\right)\frac{1+Z}{\left(\mu_{{\rm i}}^{{1}/{2}}+0.077Z\right)Z}.\end{align}

For hydrogen plasma $Z=1$, $\mu _{{\rm i}}=1$ and $\alpha _{{\rm E}}=8$, which is typical for GDT (Ivanov & Prikhodko Reference Ivanov and Prikhodko2017; Skovorodin Reference Skovorodin2019; Soldatkina et al. Reference Soldatkina, Maximov, Prikhodko, Savkin, Skovorodin, Yakovlev and Bagryansky2020), the slope factor proves to be quite small: $\gamma _{T}\sim 10^{-1}\ll 1$.

When the condition (7.30) is met, the equilibrium equations (6.5) and (6.6) are equivalent and reduced to

(7.34a,b)\begin{equation} \frac{{\rm d}\mathcal{F}}{{\rm d}\chi}+2\beta_{{\rm w}}\simeq0,\quad \mathcal{F}=-\frac{\beta_{{\rm w}}}{R}\frac{{\rm d}\beta_{{\rm w}}}{{\rm d}\chi},\end{equation}

where we use the dimensionless magnetic flux $\chi =\psi /\psi _{\mathrm {GD}}$ normalized to $\psi _{\mathrm {GD}}=2{\rm \pi} a\lambda _{\mathrm {GD}}B_{{\rm v}}$, and $\lambda _{\mathrm {GD}}$ is defined by (7.1). We also introduce the quantity $\mathcal {F}$, which has the meaning of the warm plasma transverse flow. The magnetic field distribution is given by (7.14) taking into account the corresponding renormalization of the magnetic flux: $\phi =\chi /\chi _{{\rm h}}$, where $\chi _{{\rm h}}=\psi _{{\rm h}}/\psi _{\mathrm {GD}}=(1-k^{-1})\rho _{\mathrm {NB}}/\lambda _{\mathrm {GD}}$ corresponds to the boundary of the hot ions $\phi =1$. When solving (7.34a,b), according to the approximation (7.30), the temperature should be considered constant and equal to $T_{0}$. Further, a weak temperature dependence can be found from the equality (7.33)

(7.35)\begin{equation} \frac{T}{T_{0}}=\left(\frac{p}{p_{0}}\right)^{\gamma_{T}}.\end{equation}

Internal boundary conditions (4.24) and (4.26) also merge into one

(7.36)\begin{gather} \left.\mathcal{F}\right|_{\chi\rightarrow0}\simeq\mathcal{S}, \end{gather}

(7.37)\begin{gather}\mathcal{S}=\beta_{{\rm w}0}\frac{3}{2}\sqrt{\frac{\rm \pi}{2}}\frac{\rho_{{\rm i}0}}{ \lambda_{\mathrm{GD}}}\left(\frac{W_{{\rm i}0}}{2\varPhi_{{\rm i}\parallel0}}-1\right), \end{gather}

with an additional condition

(7.38)\begin{equation} \alpha_{{\rm E}}T_{0}\left(W_{{\rm i}0}-2\varPhi_{{\rm i}\parallel0}\right)=Q_{{\rm h}0}- 2\alpha_{{\rm E}0}T_{0}\varPhi_{{\rm i}\parallel0},\end{equation}

which defines the relationship between the temperature and the pressure of the warm plasma in the bubble core. For simplicity, we further assume $\alpha _{{\rm E}0}=\alpha _{{\rm E}}$, then the temperature of the warm plasma $T_{0}$ can be explicitly found from the equality (7.38). It proves to be independent of the pressure and is determined by the ratio of the sources

(7.39)\begin{equation} T_{0}\simeq\frac{Q_{{\rm h}0}}{\alpha_{{\rm E}}W_{{\rm i}0}}.\end{equation}

External boundary conditions (4.27a,b) are reduced to the following:

(7.40)\begin{equation} \left.\beta_{{\rm w}}\right|_{\chi=\chi_{{\rm w}}}=0,\quad\Leftrightarrow\quad\left. \mathcal{F}\right|_{\chi=\chi_{{\rm w}}}=0,\end{equation}

where $\chi =\chi _{{\rm w}}$ corresponds to the outer boundary of the warm plasma $r=a+\lambda _{{\rm w}}$. The second condition is essentially the vanishing of the transverse flow at the boundary, which is consistent with the pressure gradient and the current density being finite. It is also worth noting that the position of the boundary $\chi _{{\rm w}}$ is not specified and should be found self-consistently from the solution of the equilibrium equation (7.34a,b). This means that the two conditions (7.40) are related by $\chi _{{\rm w}}$ and are formally reducible to one. In other words, one of the expressions (7.40) can be considered as a boundary condition, and the other as a definition of the boundary position $\chi =\chi _{{\rm w}}$.

7.3 Numerical solution

The boundary value problem (7.34a,b), (7.36) and (7.38) described above turns out to be rather complex, and to find its exact solution we use numerical methods.

As mentioned above, provided that the temperature is given by the expression (7.39), (7.25) relates the radius of the bubble core and the pressure of warm plasma in the core: $k=k(\beta _{{\rm w}0})$. This means that $\beta _{{\rm w}0}$ remains the only free (i.e. unknown yet) parameter that determines the equilibrium. In that regard, the boundary value problem before us is convenient to reformulate as the following Cauchy problem. The equilibrium equations (7.34a,b) can be considered as the system of ordinary differential equations for the functions $\beta _{{\rm w}}=\beta _{{\rm w}}(\chi )$ and $\mathcal {F}=\mathcal {F}(\chi )$ with the initial conditions (7.36) and $\beta _{{\rm w}}|_{\chi =0}=\beta _{{\rm w}0}$. At the same time, the quantity $\beta _{{\rm w}0}$ should be treated as a variable parameter, which is found from the external boundary condition (7.40).

To solve the formulated problem, we apply an approach similar to the shooting method. For a given parameter $\beta _{{\rm w}0}$, the equations (7.34a,b) are integrated by means of the Runge–Kutta methods, starting from the corresponding initial conditions at $\chi =0$. Computing continues until either $\beta _{{\rm w}}=\beta _{{\rm w}}(\chi ;\beta _{{\rm w}0})$ or $\mathcal {F}=\mathcal {F}(\chi ;\beta _{{\rm w}0})$ hits zero at some point $\chi =\chi _{{\rm w}}(\beta _{{\rm w}0})$.Footnote ¹¹ Thus, by varying the parameter $\beta _{{\rm w}0}$, we can formally obtain the functional dependence $\chi _{{\rm w}}=\chi _{{\rm w}}(\beta _{{\rm w}0})$. On the other hand, for a given relation $\chi _{{\rm w}}=\chi _{{\rm w}}(\beta _{{\rm w}0})$, the conditions (7.40) are actually reduced to an equation for $\beta _{{\rm w}0}$ on the finite interval $0\leq \beta _{{\rm w}0}\leq 1$. The root of this equation, which can be found using standard numerical methods, is a true value of the parameter $\beta _{{\rm w}0}$ corresponding to the real equilibrium profile $\beta _{{\rm w}}=\beta _{{\rm w}}(\chi )$.

Figure 3 shows the numerical solution of the warm plasma equilibrium equation, found using the procedure described above. There are also plotted the corresponding analytical asymptotics in the vicinity of the bubble core and the outer boundary of the warm plasma, respectively,

(7.41)\begin{align} \beta_{{\rm w}}& \sim\beta_{{\rm w}0}-\frac{4\mathcal{S}\mathfrak{C}_{{\rm h}}}{3\beta_{{\rm w}0}}\chi^{{3}/{2}}-\frac{\mathcal{S}^{2}}{6\beta_{{\rm w}0}^{2}}\chi^{2}\nonumber\\ & \quad +\frac{\mathfrak{C}_{{\rm h}}}{5}\left\{ 8-\frac{\mathcal{S}^{3}}{36\mathfrak{C}_{{\rm h}}^{2}\beta_{{\rm w}0}^{3}}+\frac{\mathcal{S}}{\beta_{{\rm w}0}\chi_{{\rm h}}} \left[1+\left(k-1\right)\left(\ln\left(\frac{\chi}{\chi_{{\rm h}}}\right)-\frac{9}{10}\right)\right]\right\} \chi^{{5}/{2}}\nonumber\\ & \quad +o\left(\chi^{{5}/{2}}\right),\quad\chi\rightarrow0, \end{align}

(7.42a,b)\begin{align} \beta_{{\rm w}}& \sim\tfrac{1}{3}R\left(\chi_{{\rm w}}\right)\left(\chi_{{\rm w}}-\chi\right)^{2}+o \left[\left(\chi_{{\rm w}}-\chi\right)^{2}\right],\quad\chi\rightarrow\chi_{{\rm w}},\end{align}

where

(7.43)\begin{equation} \mathfrak{C}_{{\rm h}}=\sqrt{\frac{1}{\chi_{{\rm h}}}\frac{\beta_{{\rm h}0}}{k+1}}. \end{equation}

The second term in the expansion (7.41) includes the contribution from the hot ions – it corresponds to the pressure profile of the warm plasma in the magnetic field induced by the hot ion current only. The third term, in turn, takes into account the diamagnetic current of the warm plasma. The non-polynomial logarithmic contribution of the hot ion current appears in the fourth term of the expansion.

Figure 3.

Profile of the warm plasma pressure depending on the magnetic flux: the numerical solution (black solid curve), left (blue dashed curve) and right (red dash-dotted curve) asymptotics given by the expressions (7.41) and (7.42a,b), respectively. Simulation corresponds to the parameters of the GDMT10 regime (see table 1).

Table 1.

Parameters of the numerical simulations. The thicknesses of transition layers $\lambda _{{\rm h}}$ and $\lambda _{{\rm w}}$ are the widths of the current profiles shown in figure 5. The warm plasma thermodynamic parameters: relative pressure $\beta _{{\rm w}0}$, density $n_{{\rm i}0}$ and temperature $T_{0}$ are the maxima of the corresponding quantities shown in figure 4. The parameter $\eta _{\parallel 0}$ is found from numerical simulations.

$^{a}$Thickness of the warm plasma transition layer according to MHD equilibrium (Beklemishev Reference Beklemishev2016): $\lambda _{{\rm w}\,\mathrm {MHD}}=7\lambda _{\mathrm {GD}}$. See the beginning of § 7.

$^{b}$Proportion of the non-adiabatic loss (4.23) in the total axial loss of the warm plasma: $\eta _{\parallel 0}=2\varPhi _{{\rm i}\parallel 0}/W_{{\rm i}0}$.

Figure 4.

Radial profiles of the warm plasma (a) relative pressure, (b) density and (c) temperature outside the bubble core $r\ge a$ for the regimes GDMT05 (blue dashed curves), GDMT10 (black solid curves) and GDMT20 (red dash-dotted curves) specified in table 1.

8.1 Analysis of the solutions

1. (i) For all solutions the relative energy density of the hot ions is significantly greater than the warm plasma relative pressure (see table 1 and figure 4)

   (8.1)\begin{equation} \beta_{{\rm w}0}\ll\beta_{{\rm h}0}\sim1.\end{equation}
   Thus, it turns out that the energy content of the plasma as well as the diamagnetic current almost entirely correspond to the hot ions, while the contribution of the warm plasma is negligible. Since the warm plasma density in this case should be relatively small, the drag force (5.7) acting on the injected ions appears to be weak as well. Therefore, this results in the energy being accumulated by the hot component rather than being transferred to the warm plasma, which is actually consistent with (8.1).

2. (ii) In table 1, listed are the simulated values of the transition layer thickness for the warm plasma $\lambda _{{\rm w}}$ and the hot ions $\lambda _{{\rm h}}$, determined from the corresponding current profiles plotted in figure 5. As expected, the total transition layer thickness $\lambda =\max \{ \lambda _{{\rm w}},\lambda _{{\rm h}}\}$ is manly determined by the hot component, since the transition layer for the warm plasma $\lambda _{{\rm w}}$ proves to be much thinner than that for the hot ions $\lambda _{{\rm h}}$. The thickness of the transition layer for the hot ions naturally turns out to be of the order of the Larmor radius of the injected ions in the vacuum magnetic field: $\lambda _{{\rm h}}\sim \rho _{\mathrm {NB}}$. At the same time, the warm plasma transition layer thickness proves to be surprisingly close to the estimate made for MHD equilibrium (Beklemishev Reference Beklemishev2016): $\lambda _{{\rm w}\,\mathrm {MHD}}=7\lambda _{\mathrm {GD}}$. This effect is probably due to a combination of the following two factors. On the one hand, the warm plasma transition layer should apparently be wider in the presence of injected ions, since the characteristic scale of the warm plasma transverse diffusion proves to be greater in the magnetic field weakened by the diamagnetism of the hot component. On the other hand, for a fixed warm plasma source, the lower maximum relative pressure of the warm plasma (8.1) seems to correspond to a radial pressure drop in a more narrow region. A more clear explanation of this phenomenon and a proper quantitative estimate of the thickness $\lambda _{{\rm w}}$ should be made based on analysis the warm plasma equilibrium equations, which is planned to be addressed in future work.

3. (iii) The equilibria presented in figures 4 and 5 are obtained within the thin transition layer limit $\lambda \ll a-\bar {r}_{\min }$ (see § 7). For all the solutions found, the expansion parameter $\lambda /(a-\bar {r}_{\min })$ proves to be less than unity. However, in the regime GDMT05 with relatively weak vacuum magnetic field (figure 5a), the expansion parameter is not too small: $\lambda /(a-\bar {r}_{\min })\simeq 7/20$. A proper description of such equilibria should be obtained by accurate numerical modelling that takes into account the finite thickness of the transition layer.

4. (iv) To simplify the system of the warm plasma equilibrium equations presented in § 7.2, we assume the temperature of the warm plasma to be approximately constant: $T\simeq \mathrm {const}$. As a result, the temperature profiles determined by the formula (7.35) actually prove to be relatively flat, but in addition, the corresponding pressure and density profiles turn out to be identically shaped (see figure 4). It is worth noting, however, that the approximation of a slowly varying temperature $T\simeq \mathrm {const}.$ is valid in the particular case of the hydrogen plasma with the energy loss factor $\alpha _{{\rm E}}$ corresponding to the GDT (Ivanov & Prikhodko Reference Ivanov and Prikhodko2017; Skovorodin Reference Skovorodin2019; Soldatkina et al. Reference Soldatkina, Maximov, Prikhodko, Savkin, Skovorodin, Yakovlev and Bagryansky2020), i.e. $Z=\mu _{{\rm i}}=1$ and $\alpha _{{\rm E}}=8$. At the same time, the nature of energy loss in the diamagnetic confinement mode may differ significantly from the gas-dynamic regime. In particular, the presence of the non-adiabatic loss (Chernoshtanov Reference Chernoshtanov2020, Reference Chernoshtanov2022) could greatly affect both longitudinal and transverse equilibrium.

5. (v) In § 6.3, we assume that the injected hot ions release the energy mainly inside the bubble core. As already mentioned, the warm plasma temperature radial profile turns out to be almost flat, while the density of the warm plasma $n_{{\rm i}}$ sharply drops just beyond the bubble core on the scale of $\lambda _{{\rm w}}$ (figure 4). In the thin transition layer approximation $\lambda _{{\rm w}}\leq \lambda \ll a-\bar {r}_{\min }$, this results in the ion–electron drag $\nu _{\mathrm {ie}}\propto n_{{\rm i}}T^{-{3}/{2}}$ indeed acting mainly inside the core. However, this appears to take place only when the source of the warm plasma is entirely contained inside the bubble core. Otherwise, the maximum density of the warm plasma could be located outside the core, which, apparently, could significantly increase the energy loss.

6. (vi) The simulations also yield the parameter $\eta _{\parallel 0}=2\varPhi _{{\rm i}\parallel 0}/W_{{\rm i}0}$ (see table 1), which represents the total proportion of the non-adiabatic loss (4.23) in the total warm plasma loss. It can be seen that the warm plasma loss turns out to be almost entirely non-adiabatic for all the considered regimes

   (8.2)\begin{equation} W_{{\rm i}0}\sim2\varPhi_{{\rm i}\parallel0},\end{equation}
   while the loss at the periphery of the bubble core $2\varPhi _{{\rm i}\parallel \mathrm {out}}=W_{{\rm i}0}-2\varPhi _{{\rm i}\parallel 0}$ proves to be negligible: $2\varPhi _{{\rm i}\parallel \mathrm {out}}\ll W_{{\rm i}0}$. One can also observe that the parameter $\eta _{\parallel 0}$ grows with increasing vacuum magnetic field. The reason for this is probably related to the high energy density of hot ions (8.1). The strong diamagnetism induced by the hot ions extending beyond the bubble core leads to a considerable weakening of the magnetic field $B$ in the vicinity of the bubble core $r\simeq a$. Therefore, the effective mirror ratio $\mathfrak {R}=B_{{\rm m}}/B$ inside the thin transition layer $a< r< a+\lambda _{{\rm w}}$ is greatly increased compared with the vacuum value $\mathcal {R}=B_{{\rm m}}/B_{{\rm v}}$. This should eventually result in the axial loss at the periphery $\varPhi _{{\rm i}\parallel \mathrm {out}}\sim n_{{\rm i}}u_{{\rm m}}S_{{\rm m}}$ being suppressed due to a decrease in the cross-section of the warm plasma flow in the mirrors $S_{{\rm m}}\sim a\lambda _{{\rm w}}/\mathfrak {R}$. At the same time, the non-adiabatic loss $\varPhi _{{\rm i}\parallel 0}\sim n_{{\rm i}0}u_{{\rm m}0}a\rho _{{\rm i}0}/\mathcal {R}$, being determined by the vacuum mirror ratio, $\mathcal {R}$, is not affected by the diamagnetism of the hot ions.

7. (vii) The warm plasma density $n_{{\rm i}0}$ turns out to be considerably lower in the regimes with a higher temperature of the warm plasma $T_{0}$ (see figure 4 and table 1). This appears to be due to the improved confinement resulting from the suppression of the axial loss in the regimes with a lower warm plasma temperature. When the non-adiabatic loss dominates, the warm plasma density can be approximately estimated from (8.2)

   (8.3)\begin{equation} n_{{\rm i}0}\sim\frac{W_{{\rm i}0}\tau_{{\rm i}\parallel0}}{{\rm \pi} a^{2}L}\propto\frac{W_{{\rm i}0}B_{{\rm m}}}{a}\frac{1}{T_{0}},\end{equation}
   where
   (8.4)\begin{equation} \tau_{{\rm i}\parallel0}\sim\frac{a}{\rho_{{\rm i}0}}\tau_{\mathrm{GD}}\propto aLB_{{\rm m}}\frac{1}{T_{0}}, \end{equation}
   is the characteristic time of the non-adiabatic loss. It turns out that the efficiency of the axial confinement does not depend explicitly on the mirror ratio $\mathcal {R}$, as for a ‘classical’ GDT $\tau _{\mathrm {GD}}\propto \mathcal {R}$, but it is enhanced with increasing the absolute value of the mirror magnetic field $B_{{\rm m}}$. In addition, the warm plasma density $n_{{\rm i}0}$ and the confinement time $\tau _{{\rm i}\parallel 0}$ indeed prove to be greater at a lower warm plasma temperature $T_{0}$. The inverse dependence of $n_{{\rm i}0}$ on $T_{0}$ is the result of the non-adiabatic loss (4.23) being proportional to $\rho _{{\rm i}0}v_{T{\rm i}0}\propto T_{0}$.

8. (viii) It can be observed from the simulations (see figure 4 and table 1) that the warm plasma temperature $T_{0}$ proves to be higher in the regimes with a stronger vacuum magnetic field $B_{{\rm v}}$. Indeed, sustaining the diamagnetic confinement equilibrium with stronger external field $B_{{\rm v}}$ (at fixed mirror field $B_{{\rm m}}$, length $L$ and radius $a$ of the bubble core) requires greater absorbed injection power $Q_{{\rm h}0}$, while the increased heating power $Q_{{\rm h}0}$, in turn, should result in a corresponding rise in the equilibrium temperature of the warm plasma $T_{0}$. This relation can be clarified by considering the force balance in the case of beam plasma (8.1)

   (8.5)\begin{equation} \varPi_{{\rm h}}\sim\varPi_{{\rm M}}, \end{equation}
   where
   (8.6a,b)\begin{equation} \varPi_{{\rm h}}\sim-\frac{Q_{{\rm h}0}}{\left\langle \dot{\mathcal{E}}\right\rangle _{\varGamma}V}\sim\frac{Q_{{\rm h}0}}{\nu_{\mathrm{ie}0}{\rm \pi} a^{2}L},\quad\varPi_{{\rm M}}\sim\frac{B_{{\rm v}}^{2}}{8{\rm \pi}}, \end{equation}
   are the characteristic energy densities of the hot ions and the vacuum magnetic field, respectively. This expression roughly corresponds to (7.25) in the limit: $\beta _{{\rm h}0}\rightarrow 1$, $\beta _{{\rm w}0}\rightarrow 0$ and $\bar {r}_{\min }/a\rightarrow 0$. Taking into account the estimate (8.3), we obtain the following relation:
   (8.7)\begin{equation} B_{{\rm v}}^{2}\sim\frac{8Q_{{\rm h}0}}{\nu_{\mathrm{ie}0}a^{2}L}\propto\frac{1}{aLB_{{\rm m}}}T_{0}^{{7}/{2}}.\end{equation}
   Therefore, the temperature of warm plasma indeed proves to increase with the vacuum magnetic field as $T_{0}\propto B_{{\rm v}}^{{4}/{7}}$.

9. (ix) Having fixed the mirror field $B_{{\rm m}}$, the warm plasma source $W_{{\rm i}0}$, the radius $a$ and the length $L$ of the bubble core, from the considered above qualitative estimates (8.3) and (8.7), we get

   (8.8a,b)\begin{equation} n_{{\rm i}0}\propto B_{{\rm v}}^{-{4}/{7}},\quad T_{0}\propto B_{{\rm v}}^{{4}/{7}}. \end{equation}
   Given the error due to the corresponding inaccuracy of the estimate (8.2), these relations agree well with the results of the simulations shown in figure 4 and table 1.

10. (x) In the present paper, we assume the classical model of transverse transport for the warm plasma, which corresponds to the Spitzer conductivity (7.28). At the same time, large gradients of equilibrium parameters, such as the magnetic field and the warm plasma density, inside the transition layer of the diamagnetic bubble could probably lead to non-classical transport, which corresponds to a much greater diffusion coefficient. In turn, a significant modification of the transverse transport could apparently have a qualitative impact on the equilibrium of the warm plasma. However, the issue concerning non-classical diffusion in the diamagnetic confinement mode remains poorly studied so far.

11. (xi) The model of the warm plasma equilibrium, constructed in § 4 within the framework of MHD, assumes isotropic pressure and gas-dynamic axial loss. This is known to correspond to the collisional regime with a filled loss cone. In the case of GDT (Ivanov & Prikhodko Reference Ivanov and Prikhodko2017), such a regime is realized when the gas-dynamic time $\tau _{\mathrm {GD}}\sim \mathcal {R}L/v_{T{\rm i}}$ significantly exceeds the kinetic time $\tau _{\mathrm {kin}}\sim \tau _{\mathrm {ii}}\ln \mathcal {R}$, where $\tau _{\mathrm {ii}}$ is the mean free time of ion–ion Coulomb collisions (Trubnikov Reference Trubnikov1965). However, when considering a diamagnetic trap, one should also take into account the effective increase in the mirror ratio $\mathfrak {R}\geq \mathcal {R}$, which results from the magnetic field weakening induced by the strong diamagnetic current of the high-pressure plasma. In addition, the exotic conditions of the diamagnetic confinement mode are likely to lead to anomalous collisionality with an effective mean free time $\tau _{\mathrm {eff}}$, which typically proves to be considerably shorter than the ‘classical’ time $\tau _{\mathrm {ii}}$. Thus, the warm plasma loss cone in the transition layer of a diamagnetic bubble can be considered filled when

    (8.9)\begin{equation} \frac{L\mathfrak{R}}{v_{T{\rm i}}}\gg\tau_{\mathrm{eff}}\ln\mathfrak{R}. \end{equation}
    This condition seems to be always satisfied at least in some neighbourhood of the bubble core $r\gtrsim a$, where the magnetic field approaches zero $B\rightarrow 0$ and, accordingly, the effective mirror ratio is greatly increased $\mathfrak {R}\rightarrow +\infty$. At the periphery of the warm plasma $r\lesssim a+\lambda _{{\rm w}}$, the loss cone may turn out to be only partially filled ($L\mathfrak {R}/v_{T{\rm i}}\sim \tau _{\mathrm {eff}}\ln \mathfrak {R}$) or even empty ($L\mathfrak {R}/v_{T{\rm i}}\ll \tau _{\mathrm {eff}}\ln \mathfrak{R}$), depending on the specific value of the anomalous mean free time $\tau _{\mathrm {eff}}$.

9.1 Future work

This work should be considered as a basis for further expansion of the theoretical model constructed here. In the near future, the equilibrium model presented in this paper will be used for a detailed investigation of the diamagnetic confinement regime in GDMT. In particular, it is planned to study the dependence of the bubble equilibrium on external conditions, such as heating power, vacuum magnetic configuration, warm plasma source, etc. In addition, the finite absorption efficiency of the injected neutral beam can be taken into account. The constructed model can also be refined by eliminating the thin transition layer approximation; this, however, would require more complex numerical simulations. In perspective, it is also planned to include in the model the effects of the warm plasma rotation and the inhomogeneity of the electrostatic potential. In the long term, collisional angular scattering of the hot ions can also be taken into account.