Plane thermonuclear detonation waves initiated by proton beams and quasi-one-dimensional model of fast ignition

Abstract

The one-dimensional problem on bilatiral irradiation by proton beams of the plane layer of condensed DT mixture with length 2H and density ρ₀ ≤ 100ρ_s, where ρ_s is the fuel solid-state density at atmospheric pressure and temperature of 4 K, is considered. The proton kinetic energy is 1 MeV, the beam intensity is 10¹⁹ W/cm² and duration is 50 ps. A mathematical model is based on the one-fluid two-temperature hydrodynamics with a wide-range equation of state of the fuel, electron and ion heat conduction, DT fusion reaction kinetics, self-radiation of plasma and plasma heating by α-particles. If the ignition occurs, a plane detonation wave, which is adjacent to the front of the rarefaction wave, appears. Upon reflection of this detonation wave from the symmetry plane, the flow with the linear velocity profile along the spatial variable x and with a weak dependence of the thermodynamic functions of x occurs. An appropriate solution of the equations of hydrodynamics is found analytically up to an arbitrary constant, which can be chosen so that the analytical solution describes with good accuracy the numerical one. The gain with respect to the energy of neutrons G ≈ 200 at Hρ₀ ≈ 1 g/cm², and G > 2000 at Hρ₀ ≈ 5 g/cm². To evaluate the ignition energy E_ig of cylindrical targets, the quasi-1D model, limiting trajectories of α-particles by a cylinder of a given radius, is suggested. The model reproduces the known theoretical dependence E_ig ~ ρ₀⁻² and gives E_ig = 160 kJ for ρ₀ = 100ρ_s ≈ 22 g/cm³.

1. INTRODUCTION

High-gain targets for the inertial confinement fusion (ICF) are meant for the ignition of a small part of the fuel with the following propagation of the thermonuclear burn wave on its main part. As the deuterium–tritium (DT) reaction rate at the temperatures of about 10 keV is much higher than the rates of other thermonuclear reactions, the ignition of the equimolar DT mixture is considered commonly as the primary task.

There are two approaches to ICF. One of the two is based on compression of a spherical layer of the fuel by a single driver providing sufficiently high values of density and temperature (Lindl et al., 2004). At that, the thermonuclear burn wave arises in the target volume. Another approach known as the fast ignition is based on using two drivers (Basov et al., 1992; Tabak et al., 1994; Gus'kov, 2013). The first driver compresses the target up to necessary value of density, while the second driver provides for fast rise of temperature. As a variant of such approach, one can consider ignition at the target center by a converging shock wave (Shcherbakov, 1983; Betti et al., 2007). A possibility of using multi-shock wave compression for fast ignition is considered by Eliezer and Martinez Val (2011).

The present paper is devoted to the variant of the fast ignition for which the burn wave arises near the target surface and propagates inside the target. Laser-generated beams of electrons (Tabak et al., 1994), and ions (Roth et al., 2001; Gus'kov, 2001; Caruso & Strangio, 2003; Honrubia et al., 2010) as well as heavy ion beams (Churazov et al., 2001; Medin et al., 2002; Gus'kov et al., 2010), macroparticles (Caruso & Strangio, 2001), and high velocity flows of matter (Gus'kov, 2013; Gus'kov & Zmitrenko, 2012) are considered as a driver for the fast heat of the highly-compressed fuel. Results of many theoretical works show initiation of burn waves propagating inside the target at certain values of beam parameters.

The contemporary concept of fast ignition assumes that the fuel is compressed up to the initial density ρ₀ ~ 10³ρ_s, where ρ_s is the density of the solid-state DT-mixture at atmospheric pressure and temperature of 4 K. Such preliminary compression of the fuel is a serious technical problem due to, particularly, Rayleigh–Taylor instability (Anisimov et al., 2013) arising at the stage of deceleration of a heavy shell while compressing the fuel. In the present paper, we study a possibility to use targets with the fuel density ρ₀ ≤ 100ρ_s. In this connection, the work by Avrorin et al. (1984) should be mentioned, where propagation of thermonuclear detonation on unevenly compressed and heated non-spherical target is studied numerically, in particular, from a small part of the hot (T = 10 keV) and dense (ρ₀ ≈ 900ρ_s) DT mixture to a cold and much less compressed part (ρ₀ ≈ 45ρ_s).

The first theoretical estimate of the ignition energy threshold for DT fuel of the density ρ₀ was obtained by Tabak et al. (1994) in the form E _ig ~ ρ₀⁻². The similar estimate with an improved value of the constant of proportionality was obtained by Atzeni (1999). In many papers (see, for example, Caruso & Pais, 1996; Atzeni, 1999; Churazov et al., 2001; Caruso & Strangio, 2003), the ignition energy threshold is determined from numerical solutions of the two-dimensional (2D) axially symmetric problem with given radius and time dependencies of the particle beam intensity. The improved dependency E _ig ~ ρ₀^−1.85 based on results of many computations for ρ₀ ≥ 50 g/cm³ ≈ 230ρ_s is proposed by Atzeni (1999). Note also a number of simple burn models based on ordinary differential equations of evolution, which allow estimating of ignition parameters for uniformly heated volumes of fuel (Nayak & Menon, 2012).

In the following three works, shell cylindrical targets with pre-compressed DT fuel irradiated from the target end are studied using 2D numerical codes. The initiation of a thermonuclear burn wave in DT fuel of the density about 500ρ_s by a heavy ion beam was considered by Churazov et al. (2001). A possibility to decrease the ignition energy by heating a small part of the fuel near the symmetry axis before the fuel density reaches its maximal value was considered by Caruso and Strangio (2003). Propagation of thermonuclear burn waves within a gold shell for DT fuel density of 50, 100, and 200 g/cm³, as well as the ignition of such a target by heavy ion beams were considered by Ramis and Meyer-ter-Vehn (2014).

Symmetrically converging plane thermonuclear burn waves initiated by laser pulses for the initial density ρ₀ = ρ_s and 5ρ_s have been studied by Khishchenko and Charakhch'yan (2013) using the model of total absorption of laser radiation in the point with the critical density. For such kind of model, epithermal particles generated by laser-plasma interaction and heating the domain of supercritical density are ignored. A slow combustion wave arises only after at least one interaction of the shock wave reflected from the symmetry plane with the ablation front and moves in the fuel preliminary compressed and heated by several shock waves. The slow wave generates before itself a compression velocity profile that increases rapidly the fuel density. As a result, the slow combustion wave can transform into two detonation waves moving in opposite directions. After reflecting the slow or detonation wave from the symmetry plane, the intensive burning continues and increases considerably the burn-up factor of the target.

Further (Charakhch'yan & Khishchenko, 2013), the study of thermonuclear burn waves at ρ₀ = ρ_s and 5ρ_s was supplemented by the waves generated by proton beams. It was shown that in spite of different ways of ignition, various models of α-particle heat, whether the burn wave remains slow or transforms into the detonation wave, and regardless of way of such a transformation, the final value of the burn-up factor depends essentially on the only parameter Hρ₀, as in the known approximate formulas for the fuel expansion in spherical geometry (Basko, 2009).

In the present paper, we extend our study to the more wide range of the initial density 5ρ_s ≤ ρ₀ ≤ 100ρ_s. The one-dimensional (1D) problem to be considered can be treated as a rough model of burning the target shown in Figure 1. The fuel of the initial density ρ₀ is in the cylindrical channel with the length of 2H surrounded by a shell of a heavy material. Following, for example, Caruso and Strangio (2003), the monoenergetic beam of protons of kinetic energy 1 MeV is considered as the ignition driver. The beam intensity is 10¹⁹ W/cm² and the duration is 50 ps. The beam of the same energy with the intensity 10¹⁸ W/cm² and the duration 500 ps is also considered. The burn-up factor and the gain obtained from the 1D calculation should be considered as the maximum possible values for the target of the same initial density and length.

Fig. 1. The target scheme.

For the first time such type of target was suggested by Pashinin and Prokhorov (1972), for compression of a gaseous deuterium by laser pulses. Note that the target was experimentally studied by Stöckl and Tsakiris (1991).

The heavy shell can be compressed by a magnetic field (Prut et al., 1979) or by heavy ions depositing energy to a certain shell of a multilayer target (Churazov et al., 2001; Dolgoleva & Zabrodin, 2006). So we suppose that it is possible to create the configuration close to that shown in Figure 1 with a cylinder of the fuel compressed to the necessary density and with two holes for injecting the proton beams.

Let us assume that the key effect determining the fuel ignition is the heat by α-particles while the lateral expansion of the fuel, the heat transfer from the fuel to the shell and the self-radiation of plasma are insignificant. Such assumption has certain grounds. For the initial density ρ₀ = 100ρ_s ≈ 22 g/cm³, the shell density in its compression can exceed the above value more than ten times. For example, Ramis and Meyer-ter-Vehn (2014) considered the gold shell with a density of 400 g/cm³. For the time 150 ps, during which a detonation wave arises in our calculations (see Section 3), the fuel radius in the target can be changed insignificantly. In the case of compression by a strong magnetic field, the heat flux between the fuel and the shell can be essentially diminished (Pashinin & Prokhorov, 1972). Finally, we refer to the work by Churazov et al. (2001), where the ignition of a similar target by a heavy ion beam is simulated numerically. The authors conclude that the radiative transfer is not determinative when the ignition.

If the above assumption is true, there is an interesting possibility to evaluate the ignition energy of the cylindrical target E _ig using 1D calculations. To ignite the target, it is necessary to leave a significant part of the α-particles energy within the fuel. To determine approximately the ignition energy of the cylindrical target in Figure 1, we introduce the parameter R _α and solve the 1D problem taking into account the escape of α-particles from a cylinder of the radius R _α identifying the latter with the beam and the fuel radiuses. If the solution of the 1D problem contains the thermonuclear burn wave, we set E _ig = πR _α²I(∞), where I = I(t) is the proton beam energy inserted at the time t per unit of cross-sectional area

$$ I\lpar t\rpar = \int\limits_0^t \, {J_b}\lpar t^{\prime}\rpar {\rm d}t^{\prime}\comma \;$$

where J _b(t) is the beam intensity given in the 1D problem.

The trajectories of α-particles are relying straight lines. A possibility to use strongly magnetized cylindrical shells confining the trajectories of α-particles and therefore reducing the beam area necessary for the ignition (Kemp et al., 2003) is not considered here. For simulation of α-particle transport, we use the track method (Brueckner & Jorna, 1973; Duderstadt & Moses, 1982) and its modification that is a numerical method for the known Cauchy problem for the kinetic steady-state homogeneous equation in Fokker–Plank approximation (Gus'kov & Rozanov, 1982).

The argument in favor of the quasi-1D model described above is its qualitative agreement with theoretical estimates and numerical simulations of the fuel ignition by the particle beam of the radius much less than the fuel size. Apart from numerical results presented in Section 5, one can see such agreement from the following simple reasoning. The mean free-path length of α-particles during the isochoric heating is about l _α ~ ρ₀⁻¹. The energy fraction of α-particles that remains within the channel of the radius R _α can be estimated as β ≈ R _α/l _α ~ R _αρ₀. Assuming that β is independent of ρ₀, obtain R _α ~ ρ₀⁻¹, and, setting that I(t) is also independent of ρ₀, obtain the known estimate E _ig ~ ρ₀⁻² (Tabak et al., 1994).

A mathematical model is presented in Section 2. The target ignition by proton beams of different intensity and the same energy, as well as generation of a detonation wave are considered in Section 3. The detonation wave reflection from the symmetry plane is studied in Section 4. Integrated characteristics for the initial density 5ρ_s, 25ρ_s, and 100ρ_s are discussed in Section 5. Conclusions are summarized in Section 6. Additionally, the track method for simulation of α-particle transport and its modification using the known Cauchy problem for the kinetic steady-state homogeneous equation in the framework of Fokker–Plank approximation are described in Appendix.

2. PROBLEM STATEMENT AND NUMERICAL METHOD

We use an equation of state (EOS) of hydrogen p = p _H(ρ, T), ε = ε_H(ρ, T) based on the wide-range semiempirical EOS model (Khishchenko, 2008). Here p is the pressure, ε is the specific internal energy, ρ is the density, and T is the temperature. The EOS provides for a good agreement with the Thomas–Fermi model with quantum and exchange corrections (Kalitkin, 1960; Kalitkin & Kuzmina, 1975) at high densities (ρ≳2 g/cm³ for H) and low temperatures, as well as gets a form of model for ideal-gases mixture of electrons and nuclei at moderate densities and high temperatures.

At the initial point of time t = 0, a motionless plane layer of the equimolar DT mixture occupies the domain 0 ≤ x ≤ H. The one-temperature EOS of the medium is described by formulas

(2.1)$$ p = p\lpar {\rm \rho} \comma \; T\rpar = {\,p_{\rm H}}\lpar {A^{ - 1}}{\rm \rho} \comma \; T\rpar \comma \; \quad {\rm \varepsilon} = {\rm \varepsilon} \lpar {\rm \rho} \comma \; T\rpar = {A^{ - 1}}{{\rm \varepsilon} _{\rm H}}\lpar {A^{ - 1}}{\rm \rho} \comma \; T\rpar \comma \;$$

where A = 2.5 is the atomic weight of the mixture.

The initial density of the fuel is specified as ρ₀ = 5ρ_s, 25ρ_s, and 100ρ_s, where ρ_s ≈ 0.22 g/cm³. The initial temperature whether equals 300 K while the initial pressure p ₀ is determined from the EOS, or is chosen (together with p ₀) on the isentrope passing through the point (ρ_s, p _a), where p _a = 0.1 MPa. Numerical results in both cases are very close to each other.

A free boundary with the pressure p _a is initially at the point x = H. The monoenergetic beam of protons of kinetic energy 1 MeV acts upon this boundary during τ_pb = 50 ps. With the exception of the short time interval δτ_pb = 0.02τ_pb, the beam has the constant intensity J ₀ = 10¹⁹ W/cm²:

$${J_b}\lpar t\rpar = \left\{{\matrix{ {{J_0}t/{\rm \delta} {{\rm \tau} _{\,pb}}\comma \; }&{t \le {\rm \delta} {{\rm \tau} _{\,pb}}\semicolon \; } \cr {{J_0}\comma \; } &{{\rm \delta} {{\rm \tau} _{\,pb}} \lt t \le {{\rm \tau} _{\,pb}}\semicolon \; } \cr {0\comma \; } &{t \gt {{\rm \tau} _{\,pb}}.} \cr } } \right.$$

The less intensive beam of the same energy with J ₀ = 10¹⁸ W/cm² and τ_pb = 500 ps is also considered. At the point x = 0, we set the symmetry condition that is equivalent to action of the identical proton beam on a symmetrical layer of the fuel.

Only the primary fusion reaction between deuteron and triton with α-particle and neutron as the reaction products

(2.2)$$ {\rm D} + {\rm T} \to {\rm \alpha} \lpar 3.5\,{\rm MeV}\rpar + {\rm N}\lpar 14\,{\rm MeV}\rpar \comma \;$$

is taken into account. Neutrons are supposed to be escaped from the fuel without interaction.

In computing of all of the coefficients entering into the model, the plasma is supposed to be completely ionized.

A mathematical model is based on the equations of one-fluid two-temperature hydrodynamics. Electron and ion heat conduction, self-radiation of plasma and plasma heating by both the proton beam and α-particles are taken into account (Afanas'ev et al., 1982):

(2.3)$$\displaystyle{{{\rm d} \rho} \over {{\rm d}t}} = - {\rm \rho} \displaystyle{{\partial u} \over {\partial x}} \comma \;$$

(2.4)$${\rm \rho} \displaystyle{{{\rm d}u} \over {{\rm d}t}} = - \displaystyle{{\partial p} \over {\partial x}}\comma \;$$

(2.5)$${\rm \rho} \displaystyle{{{\rm d}{{\rm \varepsilon} _e}} \over {{\rm d}t}} = - {\,p_e}\displaystyle{{\partial u} \over {\partial x}} + \displaystyle{\partial \over {\partial x}}{{{\varkappa} }_e}\displaystyle{{\partial {T_e}} \over {\partial x}} + \displaystyle{3 \over 2}{n_i}{k_{\rm B}}\displaystyle{{{T_i} - {T_e}} \over {{{\rm \tau} _T}}} \!+\! {D_e} + {W_e} + R\comma \,$$

(2.6)$${\rm \rho} \displaystyle{{{\rm d}{{\rm \varepsilon} _i}} \over {{\rm d}t}} = - {\,p_i}\displaystyle{{\partial u} \over {\partial x}} + \displaystyle{\partial \over {\partial x}}{{{\varkappa} }_i}\displaystyle{{\partial {T_i}} \over {\partial x}} + \displaystyle{3 \over 2}{n_i}{k_{\rm B}}\displaystyle{{{T_e} - {T_i}} \over {{{\rm \tau} _T}}} + {D_i} + {W_i}\comma \;$$

where u is the mass velocity, d/dt = ∂/∂t + u∂/∂x is the Lagrangian derivative with respect to time, p _e and p _i are the electron and ion pressure, p = p _e + p _i is the total pressure, ε_e and ε_i are the electron and ion specific internal energy, T _e and T _i are the electron and ion temperature, ϰ_e and ϰ_i are the electron (Kalitkin & Kostomarov, 2006) and ion (Silin, 1971) heat conductivity. The third term in the right parts of Eqs. (2.5) and (2.6) defines the energy exchange between electrons and ions, n _i = ρ/(Am _u) is the ion number density, A is the atomic weight, m _u is the atomic mass unit, k _B is the Boltzmann constant, τ_T is the temperature relaxation time (Kalitkin & Kostomarov, 2006). Apart from Kalitkin & Kostomarov (2006) and Silin (1971), for ϰ_e, ϰ_i and τ_T, we used formulas from Charakhch'yan et al. (2011). The rest terms in Eqs. (2.5) and (2.6) define the heating of electrons and ions by the proton beam (D _e and D _i) and by α-particles (W _e and W _i) as well as the energy exchange between electrons and self-radiation of plasma (R). The radiation pressure and the momentum transfer under deceleration of α-particles are neglected in the equation of motion (2.4). We also leave out of account terms describing the change of plasma composition induced by reaction (2.2) in the equation of discontinuity (2.3).

EOS for electrons in Eqs. (2.3)–(2.6) is taken in a form, which corresponds to thermal contribution of ideal electronic gas in completely ionized plasma of hydrogen isotopes (with taking into account degeneracy):

(2.7)$${\,p_e}\lpar {\rm \rho} \comma \; T\rpar = \displaystyle{2 \over 3}{\rm \rho} {{\rm \varepsilon} _e}\lpar {\rm \rho} \comma \; T\rpar \comma \;$$

(2.8)$${{\rm \varepsilon} _e}\lpar {\rm \rho} \comma \; T\rpar = \displaystyle{3 \over 2}{R_A}T\displaystyle{{{{\rm \beta} _e}T} \over {3{{\rm \rho} ^{2/3}} + {{\rm \beta} _e}T}}\comma \;$$

where R _A = k _BN _A(Am _u)⁻¹, N _A is the Avogadro constant,

$$ {{\rm \beta} _e} = {\lpar {\rm \pi} /3\rpar ^{2/3}}\displaystyle{{{m_e}{k_{\rm B}}} \over {{{\rm \hbar }^2}}}{\lpar A{m_u}\rpar ^{2/3}}\comma \;$$

m _e is the electron mass, ħ is the Planck constant.

EOS for ions in Eqs. (2.3)–(2.6) is taken as follows,

(2.9)$${\,p_i}\lpar {\rm \rho} \comma \; T\rpar = {\,p_{\rm H}}\lpar {A^{ - 1}}{\rm \rho} \comma \; T\rpar - {\,p_e}\lpar {\rm \rho} \comma \; T\rpar \comma \;$$

(2.10)$${{\rm \varepsilon} _i}\lpar {\rm \rho} \comma \; T\rpar = {A^{ - 1}}{{\rm \varepsilon} _{\rm H}}\lpar {A^{ - 1}}{\rm \rho} \comma \; T\rpar - {{\rm \varepsilon} _e}\lpar {\rm \rho} \comma \; T\rpar.$$

So, in the case of temperature equality T = T _e = T _i, the EOS (2.1) is satisfied, p _e + p _i = p(ρ, T), ε_e + ε_i = ε(ρ, T).

It should be stressed that, at high temperatures, T≫3β_e⁻¹ρ^2/3, relations (2.7) and (2.8) have forms of pressure and internal energy of mono-particle ideal gas of Boltzmann,

$$ {\,p_e}\lpar {\rm \rho} \comma \; T\rpar = {\rm \rho} {R_A}T\comma \; \quad {{\rm \varepsilon} _e}\lpar {\rm \rho} \comma \; T\rpar = \displaystyle{3 \over 2}{R_A}T.$$

The same relations take place at high temperatures for ionic components (2.9) and (2.10) in the EOS model from Khishchenko (2008):

$$ {\,p_i}\lpar {\rm \rho} \comma \; T\rpar = {\rm \rho} {R_A}T\comma \; \quad {{\rm \varepsilon} _i}\lpar {\rm \rho} \comma \; T\rpar = \displaystyle{3 \over 2}{R_A}T.$$

This implies in particular that for the problems in which the temperature difference of electrons and ions occurs at a high degree of plasma heating when the Boltzmann ideal gas approximation is applicable for all kinds of particles, and ionization is complete, EOS of hydrogen isotopes in Eqs. (2.3)–(2.6) can be taken as follows (Charakhch'yan et al., 2011),

(2.11)$$ {\,p_e} = \displaystyle{1 \over 2}p\lpar {\rm \rho} \comma \; {T_e}\rpar \comma \; \quad {{\rm \varepsilon} _e} = \displaystyle{1 \over 2}{\rm \varepsilon} \lpar {\rm \rho} \comma \; {T_e}\rpar \comma \; \quad {\,p_i} = \displaystyle{1 \over 2}p\lpar {\rm \rho} \comma \; {T_i}\rpar \comma \; \quad {{\rm \varepsilon} _i} = \displaystyle{1 \over 2}{\rm \varepsilon} \lpar {\rm \rho} \comma \; {T_i}\rpar \comma \;$$

where functions of p and ε correspond to the one-temperature case (2.1). Such an approach to the separation of pressure and internal energy on the electronic and ionic parts was used previously (Khishchenko & Charakhch'yan, 2013; Charakhch'yan et al., 2013; Charakhch'yan & Khishchenko, 2013) with the EOS model from Khishchenko (2008). Calculations presented further with taking into account the degeneracy of electrons by Eqs. (2.7) and (2.8) did not show significant differences between the obtained results and the case of using Eqs. (2.11).

The number of events of reaction (2.2) per unit time and unit volume is as follows from Afanas'ev et al. (1982):

$$ F = {n_{\rm D}}{n_{\rm T}}{\left\langle {{\rm \sigma} v} \right\rangle _{{\rm DT}}}\comma \;$$

where n _D and n _T are the deuteron and triton number density respectively, 〈σv〉_DT is the ion-temperature dependence of the reaction rate averaged over Maxwellian distribution of ions (Brueckner & Jorna, 1973). The burnout of fuel nuclei is described by the equation

(2.12)$$ \displaystyle{{{\rm d}{n_j}} \over {{\rm d}t}} = - {n_j}\displaystyle{{\partial u} \over {\partial x}} - F\comma \;$$

where subscripts j = D and T correspond to the cases of deuterium and tritium.

Our model of α-particle heat considers the fuel burnout, but ignores the change of the density and EOS due to the change of the plasma composition. The computations for the model of local heat by α-particles performed previously (Khishchenko & Charakhch'yan, 2013) show, that the difference between numerical results of simulations with and without taking into account the change of the plasma composition turns out insignificant.

Following Lindl (1995); Basko (2009), we introduce the local burn-up factor

$${B_{{\rm loc}}} = \displaystyle{{{n_{\rm R}}} \over {{n_{\rm R}} + {n_{\rm D}}}}\comma \;$$

where the number density of the fusion reaction events n _R is defined by the equation

(2.13)$$ \displaystyle{{{\rm d}{n_{\rm R}}} \over {{\rm d}t}} = - {n_{\rm R}}\displaystyle{{\partial u} \over {\partial x}} + F\comma \;$$

and by the initial condition n _R = 0 at t = 0. Excluding the derivative ∂u/∂x from (2.12) and (2.13), one can obtain the following ordinary differential equation along trajectories of Lagrangian particles

(2.14)$$ \displaystyle{{{\rm d}{B_{{\rm loc}}}} \over {{\rm d}t}} = {{\rm \chi} _B}\lpar 1 - {B_{{\rm loc}}}\rpar \comma \;$$

where the function χ_B = F/n _D = n _T〈σv〉_DT is the burn-up rate. The solution of (2.14) satisfying the initial condition B _loc = 0 at t = 0 has the form

(2.15)$$ {B_{{\rm loc}}}\lpar s\comma \; t\rpar = 1 - \exp\left({ - \int\limits_0^t \, {{\rm \chi} _B}\lpar s\comma \; t^{\prime}\rpar {\rm d}t^{\prime}} \right)\comma \;$$

where

$$ s\lpar x\comma \; t\rpar = \int\limits_0^x \, {\rm \rho} \lpar x^{\prime}\comma \; t\rpar {\rm d}x^{\prime}$$

is the Lagrangian coordinate.

We use two models of α-particle heat. The first one is the track method (Brueckner & Jorna, 1973; Duderstadt & Moses, 1982) based on simple physical reasons applying to a discrete media. The domain of x is divided into cells by a numerical grid. The number of α-particles per unit cross-section area which are produced within a grid cell is divided into several angular groups basing on the uniform directional distribution of α-particles. The α-particles of one angular group move along the ray intersecting one or more grid cells and are decelerated in compliance with the equation

(2.16)$$ v\displaystyle{{{\rm d}v} \over {{\rm d}{\rm \xi} }} = a\lpar x\comma \; v\rpar \comma \; \quad v\lpar {x_0}\rpar = {v_0}\comma \; $$

where v is the α-particle velocity, v ₀ ≈ 13 Mm/s is the initial α-particle velocity, a(x, v) is the deceleration (negative acceleration) of α-particles in plasma, x = x ₀ + ξμ, x ₀ is the center position of the grid cell where the group of α-particles is produced, ξ is the coordinate along the ray, μ is the cosine of the angle between the x-axis and the ray direction.

The deceleration a(x, v) = a _e(T _e(x), ρ(x), v) + a _i(T _i(x), ρ(x), v), where the first term in the right-hand part is related to electrons (Vygovskii et al., 1990) as well as the second term—to ions (Sivukhin, 1964). The functions a _e(T _e, ρ, v) and a _i(T _i, ρ, v) from Gus'kov and Rozanov (1982) were also used in our simulations and gave close results.

Solution of Eq. (2.16) together with condition v ≥ v _th(x), where v _th = (3k _BT _i/m _α)^1/2 is the velocity of α-particle thermalization, m _α is the α-particle mass, enables one to determine the contributions of the group to the right-hand parts W _e and W _i of Eqs. (2.5) and (2.6) for the respective grid cells (see Appendix).

If plasma is located within a cylinder, the α-particles escaping from the cylinder do not contribute to the plasma heating. The 1D track method enables us to take into account approximately the above three-dimensional effect. We introduce a length parameter R _α and consider the solution of Eq. (2.16) along the bounded interval 0 ≤ ξ ≤ ξ_max, where

(2.17)$$ {{\rm \xi} _{\rm max}} = {R_{\rm \alpha} }{\lpar 1 - {{\rm \mu} ^2}\rpar ^{ - 1/2}} $$

corresponds to the intersection point of the ray and the lateral boundary of the cylinder of the radius R _α, which symmetry axis coincides with the x-axis.

Let R _α → 0. Then ξ_max → 0 for all the interval −1 ≤ μ ≤ 1 with the exception of its extreme points μ = ±1. The right parts of Eqs. (2.5) and (2.6), which are integrals along μ, W _e → 0, W _i → 0, and the ignition is impossible. By increasing R _α, one can determine its value, starting from which the burn wave arises, and interpret it as the beam radius necessary for the ignition.

The second model of α-particle heat is based on the kinetic steady-state equation for the distribution function f(x, v, μ) in Fokker–Plank approximation. If all of the produced α-particles have the same initial velocity v ₀, and the diffusion of the function f(x, v, μ) in the velocity space can be ignored, the known Cauchy problem for the kinetic homogeneous equation arises (Gus'kov & Rozanov, 1982):

(2.18)$$ {\rm \mu} v\displaystyle{{\partial f} \over {\partial x}} + \displaystyle{{\partial af} \over {\partial v}} = 0\comma \; \quad {v_{{\rm th}}}\lpar x\rpar \le v \le {v_0}\comma \; \quad f\lpar x\comma \; {v_0}\comma \; {\rm \mu} \rpar = - \displaystyle{{\tilde F\lpar x\comma \; {\rm \mu} \rpar } \over {a\lpar x\comma \; {v_0}\rpar }}\comma \; $$

where $\tilde F\lpar x\comma \; {\rm \mu}\rpar $ is the distribution function upon μ of the production rate of α-particles in a unit volume near the point x. Since

$$ \int\limits_{ - 1}^1 \, \tilde F\lpar x\comma \; {\rm \mu} \rpar {\rm d}{\rm \mu} = F\lpar x\rpar \comma \;$$

for isotropic distribution of produced α-particles, $\tilde F\lpar x\comma \; {\rm \mu}\rpar = F\lpar x\rpar /2$.

The right-hand parts of Eqs. (2.5) and (2.6) have the form

$${W_{e\comma i}}\lpar x\rpar = - {m_{\rm \alpha} } \int\limits_{ - 1}^1 \, \int\limits_{{v_{{\rm th}}}\lpar x\rpar }^{{v_0}} \, f\lpar x\comma \; v\comma \; {\rm \mu} \rpar {a_{e\comma i}}\lpar x\comma \; v\rpar v{\rm d}v{\rm d}{\rm \mu} .$$

As shown in Appendix, a minor modification of the track method is a numerical method for the Cauchy problem (2.18). The quasi-1D model is similar to that described above for the track method.

The intensity of a monoenergetic proton beam is determined by the proton velocity v as J = n _pvm _pv ²/2, where n _p is the proton number density, m _p is the proton mass. Protons are supposed to be decelerated in the plasma according to Eq. (2.16) for μ = −1 and x ₀ = x _b. The value of n _p is determined by given values of the initial proton velocity and the boundary beam intensity J _b(t). Assuming that n _p is independent of x, the function J(x) is determined by the solution of Eq. (2.16) v(x). The right-hand parts of Eqs. (2.5) and (2.6) are

$$ {D_e} = \displaystyle{{{a_e}} \over {{a_e} + {a_i}}}\displaystyle{{\partial J} \over {\partial x}}\comma \; \quad {D_i} = \displaystyle{{{a_i}} \over {{a_e} + {a_i}}}\displaystyle{{\partial J} \over {\partial x}}.$$

Self-radiation of plasma is described by the steady-state transfer equation in the diffusion approximation by solid angle (Zel'dovich & Raizer, 1967). Similarly to Marchuk et al. (2009), we take into account the cooling of electrons by the inverse Compton effect using the known approximate formula (Zel'dovich, 1975; Basko, 2009). Resulting equations are as follows,

(2.19)$$ \displaystyle{{\partial {q_{\rm \nu} }} \over {\partial x}} = {\rm \kappa} \lpar {B_{\rm \nu} }\lpar {T_e}\rpar - {u_{\rm \nu} }\rpar \comma \; \quad \displaystyle{{\partial {u_{\rm \nu} }} \over {\partial x}} = - 3{\rm \kappa} {q_{\rm \nu} }\comma \;$$

where ν is the frequency, B _ν(T _e) is the Planck function, κ = κ(ρ, T _e, ν) is the absorption coefficient with accounting for the induced emission. The term in the right-hand part of Eq. (2.5) has the form

(2.20)$$ R = - \displaystyle{{\partial Q} \over {\partial x}} - \displaystyle{{4{{\rm \sigma} _{\rm T}}{n_e}U} \over {{m_e}{c^2}}}{k_{\rm B}}\lpar {T_e} - {T_r}\rpar \comma \;$$

where

$$ Q = \int\limits_0^\infty \, {q_{\rm \nu} }{\rm d}{\rm \nu} \comma \; \quad U = \int\limits_0^\infty \, {u_{\rm \nu} }{\rm d}{\rm \nu} \comma \;$$

n _e = zn _i is the electron number density, m _e is the electron mass, c is the speed of light, σ_T is the Thomson scattering cross-section of photons by free electrons, T _r is the photon temperature determined by the equality

$$ \int\limits_0^\infty \, {B_{\rm \nu} }\lpar {T_r}\rpar {\rm d}{\rm \nu} = U.$$

Numerical method for the simulations is based on splitting into physical processes. The Godunov first order accurate method in Lagrangian variables (Godunov, 1959) is used for the hydrodynamics equations. For the heat conduction and the energy exchange between electrons and ions, the implicit over time method is used. The uniform Lagrangian grid contains from 350 to 700 nodes. Number of angular groups in the track method varies from 12 to 24. The results of calculations vary insignificantly when changing to double the number of grid points and angular groups.

The grid on the frequency ν for solving of Eqs. (2.19) occupies the range from 4 to 8 decimal exponents. Number of grid nodes per a decimal exponent varies from 5 to 20. The integration with respect to the frequency in Eq. (2.20) is performed by the trapezium method.

3. IGNITION

Let us consider the initial stage of the target ignition, which is limited to the proton beam duration, t ≤ τ_pb. Suppose the heating by the proton beam is so fast that the motion of the fuel and its density change can be ignored. Such kind of heating will be referred to as isochoric. If the thermonuclear reaction and the plasma self-radiation are also ignored, the sum of Eqs. (2.5) and (2.6) takes the form

$$ {{\rm \rho} _0}\displaystyle{{\partial {\rm \varepsilon} \lpar x\comma \; t\rpar } \over {\partial t}} = \displaystyle{{\partial \lpar J + {q_e} + {q_i}\rpar } \over {\partial x}}\comma \; \quad {q_e} = {{{ \varkappa} }_e}\displaystyle{{\partial {T_e}} \over {\partial x}}\comma \; \quad {q_i} = {{{\varkappa} }_i}\displaystyle{{\partial {T_i}} \over {\partial x}}\comma \;$$

ε = ε_e + ε_i. Integrating this equation with respect to time t from 0 to τ_pb and with respect to the spatial coordinate x from H − l _p to H, where l _p is the mean free-path length of protons, setting the heat fluxes q _e = q _i = 0 at the boundaries of the integration with respect to x and ignoring the initial energy of the fuel, obtain

(3.1)$$ \bar {\rm \varepsilon} = \displaystyle{1 \over {{l_p}}} \int\limits_{H - {l_p}}^H \, {\rm \varepsilon} \lpar x\comma \; {{\rm \tau} _{\,pb}}\rpar {\rm d}x = \displaystyle{I \over {{{\rm \rho} _0}{l_p}}}\comma \; \quad I = \int\limits_0^{{{\rm \tau} _{\,pb}}} \, J\lpar t^{\prime}\rpar {\rm d}t^{\prime}\comma \;$$

where $\bar {\rm \varepsilon}$ is the average specific internal energy of plasma in the heating region at t = τ_pb, I is the beam energy per unit cross-sectional area. Since the dependence of the proton deceleration in plasma a on the density is close to the linear one, the free-path length of protons l _p ~ ρ₀⁻¹. Therefore the averaged internal energy $\bar {\rm \varepsilon}$ is independent of ρ₀ and is determined by only the beam energy I and by the initial proton energy (1 MeV) which determines the free-path length.

Epithermal protons give up the greater part of their energy to electrons, while ions are heated by both epithermal protons and more hot electrons. For the given beam energy I, the ion temperature at t = τ_pb increases both with increasing τ_pb and with increasing ρ₀ because the relaxation time of the electron and ion temperatures τ_T ~ ρ⁻¹. If the electron temperature at t = τ_pb is much greater than the ion temperature, the latter can be sufficient for the ignition after stopping the action of the proton beam at t > τ_pb.

The isochoric heating violation is connected with the rarefaction wave moving at the speed of sound from the free boundary into the target, and with the shock wave, arising due to the rapid growth of pressure in the heating region. As the condition of the isochoric heating violation, we take the inequality c_sτ_pb > l _p, where c _s is the speed of sound. Since l _p ~ ρ₀⁻¹, with increasing ρ₀ one should decrease the beam duration τ_pb to keep close to the isochoric heating, which, as it follows from Eq. (3.1), needs the appropriate increase in the beam intensity J ₀ to conserve the average internal energy in the heating region.

The ion temperature and the mass velocity profiles by the Lagrangian variable with respect to the free boundary s − s _b at the time t = τ_pb for the beams of the same energy and different durations τ_pb = 50 and 100 ps, as well as for the two values of the initial density ρ₀ = 25ρ_s and 100ρ_s are presented in Figure 2. The parameter R _α limiting the trajectory of α-particles is selected near the ignition boundary of the target by at least one of the beams (R _α = 0.4 and 0.1 mm for ρ₀ = 25ρ_s and 100ρ_s respectively).

Fig. 2. The ion temperature (a, b) and the mass velocity (c, d) as functions of the Lagrangian variable with respect to the free boundary s − s _b at t = τ_pb: ρ₀ = 25ρ_s, R _α = 0.4 mm (a, c) and ρ₀ = 100ρ_s, R _α = 0.1 mm (b, d) for different beams of the same energy, J ₀ = 10¹⁹ W/cm², τ_pb = 50 ps (solid lines) and J ₀ = 10¹⁸ W/cm², τ_pb = 500 ps (dashed lines).

First, we examine the case of ρ₀ = 25ρ_s (see Figs. 2a and 2c).

At τ_pb = 50 ps, as follows from the corresponding velocity profile (the solid line in Fig. 2c), the heating is almost isochoric: the formation of the shock wave has not started yet, and the rarefaction wave occupies only a small part of the heating region. Although the ion temperature (the solid line in Fig. 2a) is relatively small, its subsequent growth due to heating by more hot electrons leads to the target ignition.

With increasing the duration τ_pb up to 500 ps, the proton beam heating is no longer isochoric. As seen in the velocity profile (the dashed line in Fig. 2c), a continuous compression wave is formed inside the heating region. Then this wave should transform into a shock wave. A rarefaction wave adjoins the compression one, which in the absence of heat by α-particles should lead to a rapid drop in the amplitude of the shock wave. However, the ion temperature (dashed line in Fig. 2a) is sufficient to ignite the target. Analysis of the calculation results shows that a significant heating of the plasma by α-particles takes place in the field of the compression wave in Figure 2c. As a result, the compression wave is rapidly converted to the detonation wave of the well-known type (Landau & Lifshitz, 1987) with a rarefaction wave, which is adjacent to the detonation wave front.

We now turn to the case ρ₀ = 100ρ_s (see Figs. 2b and 2d). When τ_pb = 50 ps, the ion temperature (the solid line in Fig. 2b) is much greater than in the case of the same beam and ρ₀ = 25ρ_s due to the decrease of the temperature relaxation time τ_T ~ ρ⁻¹. As can be seen from the velocity profile (the solid line in Fig. 2d), the shock wave begins to form, and the rarefaction wave takes until a small portion of the heating region. The ignition mechanism for these parameters will be discussed below.

For the beam with τ_pb = 500 ps, the target does not ignite. As can be seen from the velocity profile (dashed line in Fig. 2d), the shock wave at t = τ_pb is too far away from the heating region and no longer can transform into a detonation wave. After some time, the rarefaction wave overtakes the shock front and starts to reduce its amplitude. Calculations show that for the fuel density ρ₀ = 10³ρ_s, a similar flow pattern without the target ignition occurs for the beam with τ_pb = 50 ps.

The aforementioned lack of the ignition by the beam with the intensity of 10¹⁸ W/cm² for the fuel density ρ₀ = 100ρ_s ≈ 22 g/cm³, and the ignition by the beam with the intensity of 10¹⁹ W/cm² for the same density, as well as the lack of the ignition by this beam for the fuel density ρ₀ ≈ 220 g/cm³, correspond to the density dependence of the minimum beam intensity required for ignition (Atzeni, 1999) based on results of two-dimensional calculations, which gives a value of approximately 6 × 10¹⁸ W/cm² for ρ₀ = 22 g/cm³ and about 5 × 10¹⁹ W/cm² for ρ₀ = 220 g/cm³.

Results of calculations for a target with ρ₀ = 100ρ_s, R _α = 0.1 mm, and for the beam with τ_pb = 50 ps, J ₀ = 10¹⁹ W/cm² are presented below. Consider first the case of the fuel layer half-width H = 0.5 mm, which corresponds to the value of the parameter Hρ₀ ≈ 1 g/cm².

Figure 3 shows the ignition process of the target. The functions T _e(x), T _i(x), ρ(x), u(x) and the burn-up rate χ_B(x), which determines the local burn-up factor B _loc by Eq. (2.15), are given at the termination of the beam 50 ps and three subsequent times 100, 150 and 200 ps. A small discontinuity in the function χ_B(x) is caused by the discontinuity in the approximation formula for 〈σv〉_DT(T _i) from Brueckner and Jorna (1973).

Fig. 3. The temperature (a, the solid lines correspond to ions, the dashed lines—to electrons), the burn-up rate (b), the density (c) and the mass velocity (d) spatial profiles in the formation of the detonation wave at t = 50 (1), 100 (2), 150 (3), and 200 ps (4) for the target with ρ₀ = 100ρ_s, R _α = 0.1 mm, H = 0.5 mm and for the beam with τ_pb = 50 ps.

At t = 50 ps, the maximal ion temperature is about 100 MK, which is sufficient for the ignition. The shock wave has not yet formed, and the burning wave, which can be identified with behavior of the function χ_B(x), is continuous. At t = 100 ps, the shock wave is formed, and the burning wave is a little behind it and still continuous. At t = 150 ps, the burning wave has caught up with the shock wave and transformed it into the detonation one. Comparing functions T _i(x) and χ_B(x) at t = 150 and 200 ps, one can see that the temperature and the burn-up rate at the detonation wave front quickly grow with time. Before the detonation wave one can see a precursor caused by electron heat conduction and α-particles (near the wave front) and by the self-radiation of high-temperature plasma (at large distances from the front).

To conclude this section, we present the numerical results for the target with the initial density ρ₀ = 100ρ_s at different values of R _α, those limiting trajectory of α-particles, and for different models of α-particle transport: the Fokker–Plank equation and the track method. Table 1 shows the values of the burn-up factor for targets with two values of H = 0.5 (B ₀) and 2.5 mm (B ₁). These values correspond to a sufficiently late point in time, when their significant growth ceases. Zero values B ₀ and B ₁ means that their value is less than 0.01, and are interpreted as a lack of ignition.

Table 1. The burn-up factors for the targets with ρ₀ = 100ρ_s, H = 0.5 (B ₀) and 2.5 mm (B ₁), for different values of the parameter R _α and different models of α-particle transport: the Fokker–Plank equation (FP) and the track method (TM)

At R _α = 0.15 mm, the factors B ₀ and B ₁ are the same for both models, and slightly differ from respective values at R _α = 0.1 mm. At R _α = 0.05 mm, the Fokker–Plank equation still leads to the ignition of the target with somewhat smaller values of B ₀ and B ₁ compared with the case of R _α = 0.1 mm, and in the case of the track method B ₁ = 0. This means that the detonation wave is not formed for H = 2.5 mm since the shock wave at large distances is extinguished by the rarefaction wave. If the symmetry plane is located close enough to the edge of the target (H = 0.5 mm), the shock wave at the moment of reflection is strong enough to ignite the target, creating a reflected detonation wave. At R _α = 0.03 mm, the ignition is possible only if H = 0.5 mm with small values of B ₀. At R _α = 0.01 mm the ignition is absent. To estimate the energy of ignition, we set R _α = 0.1 mm, wherein the reliable ignition occurs in both models.

4. REFLECTION OF DETONATION WAVE

In the case of H = 0.5 mm, reflection of the detonation wave from the symmetry plane x = 0 is shown in Figure 4. The functions T _i(x), ρ(x), u(x) and the local burn-up factor B _loc(x) are presented in the six time points. The first two points (200 and 250 ps) meet the detonation wave before its reflection from the symmetry plane. Note the relatively small values of the local burn-up factor (about 0.05). Pay attention to the velocity profile at t = 250 ps (curve 2 in Fig. 4d). One can see that the precursor of the detonation wave transforms into a flow with a linear velocity profile between the wave and the symmetry plane. Such flows, which are characterized by the dependence

Fig. 4. The spatial profiles of the ion temperature (a), the local burn-up factor (b), the density (c) and the mass velocity (d) in reflection of the detonation wave at points of time 200 (1), 250 (2), 300 (3), 350 (4), 400 (5), and 450 ps (6) for the parameters of the problem, as shown in Figure 3.

(4.1)$$ u\lpar x\comma \; t\rpar = {\rm \varphi} \lpar t\rpar x\comma \;$$

are well known in hydrodynamics (Sedov, 1972) and combined into a large family of solutions, each member of that is determined by an arbitrary function of one argument and three arbitrary constants. For all of these solutions, the following formula takes place

(4.2)$$ {\rm \rho} \lpar s\comma \; t\rpar = {\rm \rho} \lpar s\comma \; {t_1}\rpar \exp\left({ - \int\limits_{{t_1}}^t \, {\rm \varphi} \lpar t^{\prime}\rpar {\rm d}t^{\prime}} \right)\comma \;$$

which is obtained by integrating an ordinary differential equation arising under the substitution (4.1) to the equation of discontinuity (2.3). In the case under consideration, since ϕ(t) < 0, the density in the region of the precursor and, in particular, in the symmetry plane, grows with time.

The last four points of time in Figure 4 (300, 350, 400, and 450 ps) give the flow pattern after reflection of the detonation wave. One can see that the flow with a linear velocity profile arises again between the reflected detonation wave and the symmetry plane. Its special feature is the proximity of the thermodynamic functions to the constant upon x value, which depends on time. The local burn-up factor increases significantly after reflection of the detonation wave and continues to grow because the burn-up rate χ_B(x,t) decreases with time slowly enough, remaining almost constant function of x. Note that the spatial homogeneity of the thermodynamic functions after the reflected wave front is not typical for the solutions of the equations of hydrodynamics in the case of spherical or cylindrical geometry. For example, a well-known solution on the problem of a converging spherical shock wave (Guderley, 1942; Stanyukovich, 1955), extended through time after the shock wave collapse, has zero density at the point of symmetry (Goldman, 1973).

Solution of the equations of hydrodynamics with a linear velocity profile and constant on x values of thermodynamic functions is contained in the above-mentioned family from Sedov (1972). It also can be obtained directly from the equation of motion (2.4). Putting in Eq. (2.4) ∂p/∂x = 0 and substituting the velocity in the form (4.1), we obtain an ordinary differential equation $\dot{\rm \varphi} + {\rm \varphi}^2 = 0$, whose solution has the form

(4.3)$$ {\rm \varphi} \lpar t\rpar = \displaystyle{1 \over {C + t}}\comma \;$$

where C is an arbitrary constant. Substituting (4.3) to Eq. (4.2), we obtain

(4.4)$$ {\rm \rho} \lpar t\rpar = \displaystyle{{{{\rm \rho} _0}C} \over {C + t}} = {{\rm \rho} _0}C{\rm \varphi} \lpar t\rpar \comma \;$$

where constant ρ₀ is chosen so that ρ(0) = ρ₀.

Since the density ρ is independent of the spatial coordinate x, the velocity has a linear profile along the Lagrangian coordinate s = ρx, the slope of which, as follows from Eqs. (4.1), (4.3), and (4.4), u/s = (ρ₀C)⁻¹, is independent of time. Shown in Figure 5 velocity profiles along the Lagrangian coordinate indicate that this property is satisfied with good accuracy for the problem in the region between the symmetry plane and the reflected detonation wave.

Fig. 5. Mass velocity profiles along the Lagrangian coordinate after reflection of the detonation wave for the same parameters of the problem and the time points as in Figure 4.

Using the simulation results, one can determine the approximate value of constant C in Eq. (4.3). It is convenient to substitute time $\bar t = t - t_{\ast}$ for t in Eqs. (4.3) and (4.4), where $t_{\ast}$ is the time point of appearance of the reflected detonation wave. Note that the density on the symmetry plane as a function of time ρ(0, t) increases as the incident detonation wave approaches to the symmetry plane, and begins to decrease in accordance with Eq. (4.4) after appearance of the reflected wave. It is therefore natural to define the time of formation of the reflected wave $t_{\ast}$ and constant ρ₀ in Eq. (4.4) by the condition max_t(ρ(0, t)) = ρ(0, $t_{\ast})={\rm \rho}_0$. The result is ρ₀ ≈ 80 g/cm³, $t_{\ast}$ ≈ 270 ps.

The calculation results are compared with Eqs. (4.3) and (4.4) for the four time points shown in Figure 4 after reflection of the detonation wave, t = 300, 350, 400, and 450 ps. As the slope ${\rm \varphi}$ and density ρ for each time, the values at the point x approximately two times less than the coordinate of the front of the reflected detonation wave are selected. We denote by $\bar{\rm \varphi}$ and $\bar {\rm \rho}$ the corresponding values included in Eqs. (4.3) and (4.4). At first, we put $\bar{\rm \varphi} = {\rm \varphi}$ and calculate C from Eq. (4.3), then $\bar {\rm \rho}$ by Eq. (4.4) and the relative error of calculating the density ${\rm \delta \rho} = \vert 1 - \bar {\rm \rho}/{\rm \rho}\vert $. The results of calculations are shown in the first three columns of Table 2. For the first three values of t, constant C varies within 10% and the error δρ is within 6%. At the last time point t = 450 ps, constant C changes by about 10% more, and the error δρ grows to 14%. Nevertheless, it is possible to select a single value of constant C for all four points of time and get a relatively small magnitude errors ${\rm \delta\varphi} = \vert 1 - \bar {\rm \varphi}/{\rm \varphi}\vert $ and δρ. Such a possibility is demonstrated in the last two columns of Table 2, where the values of these errors are presented for C = 40 ps. It can be seen that both errors are within 3%.

Table 2. Comparison of simulated values of φ and ρ with the analytical Eqs. (4.3) and (4.4): the value of C from Eq. (4.3) as well as the relative error δρ for the given value φ and the relative errors δρ₄₀ and δφ₄₀ for C = 40 ps at some points of time t

Let us consider now the target with the layer half-thickness H = 2.5 mm, which corresponds to the value of the parameter Hρ₀ ≈ 5 g/cm². Figure 6 shows the ion temperature, the density, the mass velocity, and the local burn-up factor as functions of x at eight points of time. The first four points give a picture of the flow before the detonation wave reflection from the symmetry plane. Detonation wave intensity increases with time, as can be seen from the profiles of the ion temperature and the mass velocity, reaching values of the order of 700 MK and 2 Mm/s at the instant of reflection. Such a high velocity in the incident detonation wave at the reflection leads to the ion temperature after reflection of about 1.3 GK (see curve 5 in Fig. 6a). Further, as in the above case H = 0.5 mm, the flow with a linear velocity profile and the weak dependence of the thermodynamic functions upon x arises between the reflected detonation wave and the symmetry plane.

Fig. 6. Spatial profiles of the ion temperature (a), the local burn-up factor (b), the density (c) an the mass velocity (d) at the time points 0.2 (1), 0.4 (2), 0.6 (3), 0.8 (4), 1 (5), 1.2 (6), 1.4 (7) and 1.8 ns (8), corresponding to the incident (1–4) and the reflected (5–8) detonation waves for H = 2.5 mm (another parameters of the problem are the same as in Fig. 3).

Note the rapid rise of the local burn-up factor after reflection of the detonation wave as close to the symmetry plane as in the periphery. At the last instant shown in Figure 6, the local burn-up factor in the region 0 ≤ x ≤ 3.5 mm, containing about 80% of the fuel mass, varies from about 0.7 to 0.6.

The burn-up factor B(t) is presented in Figure 7 for the two above values of H. The point of reflection of the detonation wave t = t _* is indicated. One can see that the expansion stage t > t _* gives the main contribution to the final value of the factor in the case of H = 0.5 mm. For H = 2.5 mm the burn-up factor increases from 0.4 to 0.7 at the expansion stage.

Fig. 7. The time dependence of the burn-up factor for H = 0.5 (1) and 2.5 mm (2) (another parameters of the problem are the same as in Fig. 3); asterisks denote the instant of the detonation wave reflection.

5. INTEGRATED CHARACTERISTICS

Generalized characteristics of the target for three values of the initial density ρ₀ = 5ρ_s, 25ρ_s, and 100ρ_s are shown in Table 3. For ρ₀ = 100ρ_s, the characteristics are obtained by the results of the above calculations with the parameter R _α = 0.1 mm. For other values of ρ₀, the parameter R _α is chosen from the condition ρ₀R _α = const. For each value of ρ₀ two values of the layer half-width H, corresponding to the two given values of the parameter Hρ₀ ≈ 1 and 5 g/cm², are chosen. Interpreting the parameter R _α as the radius of the cylindrical target, shown in Figure 1, one can define the target ignition energy E _ig = πR _α²I(∞) and the mass M = 2HπR _α²ρ₀. These values are given in the table together with the values of the gain with respect to the energy of neutrons and the burn-up factor after the burning process.

Table 3. The ignition energy for one proton beam E _ig, the target mass M, the burn-up factor B and the gain G for three values of ρ₀, the cylinder radius R _α~ρ₀⁻¹ and different values of the layer half-width H determined by two given values of the parameter Hρ₀

As in the case of the well-known approximate formula for the expansion of a spherical target (Basko, 2009), the burn-up factor B, as well as the gain G, are determined by the parameter Hρ₀. If Hρ₀ ≈ 1 g/cm², B ≈ 0.3, G ≈ 200. If Hρ₀ ≈ 5 g/cm², the burn-up factor increases by more than a factor of two, which, together with the increase in the fuel mass five times, gives the gain G > 2000.

Since the parameter R _α, determined by the condition R _αρ₀ = const, provides the target ignition, the ignition energy decreases with increasing ρ₀ under the law E _ig ~ ρ₀⁻² in accordance with the known theoretical estimate (Tabak et al., 1994). For ρ₀ = 100ρ_s ≈ 22 g/cm³, the ignition energy E _ig = 160 kJ.

6. CONCLUSIONS

Computational experiments with two proton beams of the same energy and different intensities show that the beam intensity required for the ignition increases with the initial fuel density in accordance with a known formula, obtained by results of two-dimensional calculations.

If the target is ignited, a detonation wave, which is adjacent to the front of the rarefaction wave, occurs. The formation of the detonation wave precedes a short time interval with a shock wave that goes ahead of a slow combustion wave. When the combustion front overtakes the shock wave, the latter transforms into the detonation wave.

There is no ignition if the distance between the shock front and the region of the fuel heating by protons is too large at the end of irradiation. Over time, the rarefaction wave overtakes the shock front and reduces the intensity of the shock wave.

For small values of R _α, limiting trajectories of α-particles, and small target thickness (Hρ₀ ≈ 1 g/cm²), the target ignition is possible after the shock wave reflection from the symmetry plane. With increasing H, starting from a certain value, the target ceases to ignite due to the aforementioned decrease in the intensity of the shock wave that is overtaken by the rarefaction wave.

Upon reflection of the detonation wave from the symmetry plane, the flow with the linear velocity profile along the spatial variable x and with a weak dependence of the thermodynamic functions upon x, close to one of the representatives of a certain family of solutions of hydrodynamic equations, occurs between this plane and the reflected detonation wave. Combustion efficiency at this stage is largely due to the spatial homogeneity of the burn-up rate (the product of the DT reaction rate and the triton number density) whose integral along the trajectory of a Lagrangian particle defines the local burn-up factor.

Solution of the equations of hydrodynamics with a linear velocity profile and constant in x values of thermodynamic functions is found analytically up to an arbitrary constant. The possibility to choose this constant so that the solution describes with good accuracy the change in the density with time after reflection of the detonation wave is shown.

For the parameter Hρ₀ ≈ 1 g/cm², the main contribution to the final value of the burn-up factor B ≈ 0.3 is given by the expansion stage after reflection of the detonation wave. For Hρ₀ ≈ 5 g/cm², the burn-up factor is increased by the expansion stage from 0.4 to 0.7. At Hρ₀ ≈ 1 g/cm², the gain with respect to the energy of neutrons G ≈ 200. At Hρ₀ ≈ 5 g/cm², the burn-up factor increases by more than a factor of two, which, together with the increase in the fuel mass five times, gives the gain G > 2000. The above values of the burn-up factor and the gain should be considered as the maximum possible values for the real cylindrical target of the same initial density and length.

The quasi-1D model, limiting trajectories of α-particles by a cylinder of a given radius, reproduces known theoretical dependence of the ignition energy E _ig ~ ρ₀⁻². For ρ₀ = 100ρ_s ≈ 22 g/cm³, the quasi-1D model gives E _ig = 160 kJ.