Modeling ablator grain structure impacts in ICF implosions

High-density carbon (HDC) shells^1,2 are used as ablators in inertial confinement fusion (ICF) experiments that aim to generate fusion by compressing and heating the deuterium–tritium (DT) fuel that is layered within.^3–6 The high density of these HDC shells (⁠ $∼ 3.5$ g/cm³) gives certain advantages in comparison with other candidate ablators, such as plastic (CH, $∼ 1$ g/cm³), which have lower density. These advantages, described in detail elsewhere,^3,6,7 relate both to various aspects of the ablator (shell) physics itself and to benefits derived from changed requirements for the (indirect) radiation drive, in particular the shorter required irradiation time.

HDC capsules (shells) were used in recent ICF experiments that reached the so-called burning plasma regime and achieved record fusion yields.^8,9 Despite these promising results, there are certain known or suspected complications with using HDC as the ablator material. Among these are a larger perturbation due to the fill-tube that supplies the DT fuel⁷ and perturbations associated with incomplete melt of the HDC when shocked.¹⁰ Additionally, HDC implosions have persistently observed compressed-fuel (DT) areal density (ρr) lower than simulation predictions,^11,12 while the HDC itself achieves the expected compression.¹³ The present work focuses on another (possibly related) complication, the suspected small-scale density nonuniformity of HDC ablators¹⁴ due to grain structures. HDC is polycrystalline diamond, consisting of diamond grains (sp3-bonded carbon) surrounded by graphitic carbon (sp2-bonded) at a lower density;^14,15 we describe this in more detail later.

These complications generate uncertainty about the degree to which HDC designs can be pushed, for example, toward greater compressions (higher ρr) and higher yield (or, uncertainty about the design accommodations that might be required to do so). Present HDC designs aim for a first shock strength >12 Mbar to mitigate the risk of perturbations associated with incomplete melt of the HDC.^10,16,17 In general, for ICF, the first shock is important in setting the initial entropy in the DT fuel and, therefore, its compressibility (see, e.g., Ref. 18). While the expected transferred shock strength into the DT fuel is below the $∼ 2$ Mbar threshold required to stay near the Fermi degenerate adiabat, the desire to stay above 12 Mbar for the first shock in HDC is nonetheless a design constraint.

The presence of small-scale density nonuniformity due to grain structures might also introduce limitations or constraints on HDC designs. As a shock passes through a nonuniform density field, the post-shock flow is perturbed; these effects can be calculated precisely in the (linear) limit of weak density nonuniformity.^19,20 At the fuel–ablator (DT-HDC) interface, these perturbations will lead to Richtmyer–Meshkov-like instability (RMI-like)^21,22 and early-time mixing of fuel and ablator. The resulting perturbed interface and associated flows then serve as the starting point for further perturbation or mixing growth, due, for example, to Bell–Plesset effects,^23,24 compression,^25–27 and various fluid instabilities (e.g., Rayleigh–Taylor^28,29) of the local density and temperature gradients during the acceleration and deceleration of the capsule during the implosion.

Such flow perturbations and mixing have the potential to degrade HDC implosions in different ways. Here we discuss effects which do not involve mix of the HDC ablator into the hot spot (“hot mix”), since we do not expect the small-scale grain perturbations alone to generate substantial hot mix (although they might contribute to it in concert with other larger perturbation seeds such as the fill-tube^7,30). Entropy increase in the DT fuel associated with fuel–ablator mixing has been proposed as one possible implosion degradation mechanism,^31,32 causing reduced DT fuel compression. Such mixing of hot ablator (HDC) can also in principle heat the DT fuel, increasing its entropy and again reducing compression. Simulations with other small-scale perturbation sources, internal and surface defects and voids, find that reduced compression is possible without substantial entropy or mixing effects, due to nonradial displacement of fuel mass from perturbed flows.³³ There is the possibility for impacts even if grain-like perturbations do not lead to reduced compression: as HDC implosions reach higher fuel burn-up fraction, it is possible that such fuel–ablator mixing may limit burn, owing to the outer portion of the DT containing too much mixed mass. While steps can be taken to eliminate certain capsule imperfections, we might regard the grain structures of HDC as ineliminable, although they do have controllable size.¹⁴ Therefore, it is important to be able to anticipate the impacts of these structures on ICF implosions.

A prior estimate of the impact of grains, for the National Ignition Facility (NIF) shot N161204,³⁴ found a 25% decrease in fusion yield from grains alone.³⁵ However, in general, the impacts of grains in implosions are quite uncertain, in large part because the very small spatial scales (∼nm) that are involved make direct resolution impractical in ICF simulations (there are also modeling uncertainties associated with treating grains, some of which we touch on later in the discussion, Sec. VI). As a result, these aforementioned simulations needed to implement a strategy to model the grains at a reduced resolution, a process with likely no perfect answer and uncertain impacts. Indeed, we will show in this work that, even within a specific grain model, there are different ways of carrying out such a de-resolving process and these lead to different results for the grain-induced perturbations and mixing.

The overarching goal of the present work is to develop grain modeling as a step toward reducing the design uncertainties associated with grains. Further, we aim to develop a model that is straightforward to utilize in hydrodynamic codes and which has a notion of convergence to a “true” grain case. First, in Sec. II, we present a grain model which has both these properties; this model also does a reasonable job at matching the limited available experimental measurements related to grains and their associated perturbations. As part of presenting this grain model, in Sec. III, we introduce two different techniques for de-resolving the grain model to greatly reduce the computational cost associated with grain modeling.

Second, in Sec. IV, we study the perturbations generated by shock interaction with grains, in a setup that allows us to simulate resolved grains (a small portion of HDC ablator, rather than a full capsule). This allows us to study the effects of de-resolving the grains in comparison with the true case. We show, in particular, that the perturbed kinetic energy deposited by the first shock transiting the grains depends on the method by which the grains are scaled (de-resolved), and we devise a method by which the deposited kinetic energy perturbation is approximately conserved as the grains are de-resolved. This method holds constant the density dispersion, $σ ρ$⁠, of the medium during scaling (as well as the mean density, $ρ ¯$⁠). As we will discuss in Sec. IV B, this is a sensible thing to do when regarded from the perspective of the aforementioned linear theory,²⁰ if the grains are considered (approximately) as a random density field: doing so means preserving a number of impacts of nonuniform density, including modifications to the post-shock mean pressure and density, as well as modifications to the shock speed. Caution is warranted in utilizing this density-dispersion-conservation method, as the geometry and lengthscales of the perturbations can matter, a point we highlight throughout.

Third, in Sec. V, we study the impact of HDC grains in an ICF implosion, and the effects of de-resolving the grains on fuel–ablator mixing. In contrast to the single-shock setup of Sec. IV, the HDC is now subjected to multiple shocks and a complex acceleration history, and there is an interface between the HDC and DT fuel. Because of the great computational expense in modeling grains, we construct a hydrodynamic-only surrogate implosion, which closely tracks the source radiation-driven implosion, as shown in Sec. V A. For this ICF implosion case, we focus, in particular, on the fuel–ablator mix width as a quantitative measure. We show that quantitative, and in some cases qualitative, mixing dynamics depend on how the grains are de-resolved. Further, we show that while not perfect (as to be expected), the strategy of conserving the density dispersion produces reasonable agreement for the mix width time-history in the tested cases, even at much lower resolution. This supports a key goal of the present work, enabling faster future estimates of grain impacts. Note that the use of the hydrodynamic-only implosion complicates our ability to study here some of the possible grain impacts introduced above; these will be projects for future work.

Section V also shows two other results worth highlighting. The result suggests that a microscale grain model (⁠ $∼ 2 × 1$ μm grains) leads to somewhat greater mixing than a nanoscale grain model (⁠ $∼ 200 × 200$ nm grains), despite the latter leading to significantly more perturbed kinetic energy deposition (a fact that highlights the importance of geometry and lengthscale). Experiments with HDC capsules have moved from microcrystalline-grain capsules to nanocrystalline-grain ones, and the highest-performing experiments to date use nanocrystalline HDC (though this could also be due to other effects, for example, changes in capsule quality or differences in the melt properties of nanocrystalline vs microcrystalline HDC). Section V also shows, qualitatively, the complex interaction between conduction-driven density gradients at the HDC–DT interface and perturbations from the grains. To the extent the conductivity of DT or HDC are uncertain in the warm dense matter regime of the capsule during the implosions, this can be expected to generate uncertainty in the mixing dynamics.

Fourth, following Sec. V, we further discuss the results, assumptions, and caveats in Sec. VI, before briefly concluding in Sec. VII. While we do not discuss it here, we note that other ablator materials, for example, beryllium, can also have grain structures or small-scale density nonuniformity. A different type of small-scale density nonuniformity is present in wetted-foam ignition designs (e.g., Refs. 36 and 37) where liquid DT fuel is held within a (carbon-based) foam lattice. This lattice structure leads, as in the present grain case, to very small-scale density nonuniformity, leading to perturbed shocks and post-shock flow.^20,38–43 Measured in terms of the normalized density dispersion, $σ ρ / ρ ¯$⁠, the wetted foams are more nonuniform than the grains presently considered, which will lead to larger fractional perturbations. Prior work has found a few percent change in post-shock mean pressure and density, and corresponding changes in shock velocity, when a shock is propagated through a wetted foam with resolved foam lattice, compared with the case when the lattice is ignored (homogeneous case).⁴¹ As in the grain cases, the foam lattice is generally impractical to resolve in integrated ICF simulations. Aspects of the approaches developed here may be applicable in these other scenarios.

We now describe the model for grains used in the present work. This model is inspired in part by measurements made during the course of work on fabricating high-density carbon (HDC) materials with controllable grain sizes.¹⁴ Such a model has been used in simulations to compare against the experimentally observed spectrum of velocity perturbations on a shock that has passed through an HDC ablator,⁴⁴ comparing with data similar to that in Ali et al.¹⁰ The comparison is favorable, but for the particular HDC samples used (and unlike in capsule implosion cases), the simulations suggest that surface roughness made the dominant contribution to the shock velocity perturbations. Thus, while these data do place some constraints on the grain modeling, it would be desirable in the future to place more stringent experimental constraints on grain models. While our grain model is based on observations, it is still an idealized and approximate picture; a variety of caveats are discussed later in Sec. VI.

In this model, the grain structure of the (ablator) material is taken to be made up of two separate components, the grain, at a density ρ_g, and the interstitial, at a density ρ_i. Here the grain is (diamond-like) HDC, at a density $ρ g = 3.5006$⁠, and the interstitial is taken to have density $ρ i = 2.15$ (graphite). Given a total volume of material, $V = V g + V i$⁠, made up a volume V_g at ρ_g and V_i at ρ_i, we have an average density $ρ ¯ = ( ρ g V g + ρ i V i ) / V$⁠.

The interstitial material is confined to thin layers which separate the grains. At present, we work in 2D and represent the grains as a repeated pattern (in two orthogonal directions) of an individual unit; a diagram of this unit is shown in Fig. 1. It consists of an $α g × g$ rectangular grain portion (density ρ_g), surrounded on all sides by an interstitial (density ρ_i) of width $i / 2$⁠. Thus, α is a parameter that controls the aspect ratio of the grain unit. When repeatedly tiled, the (interstitial) spacing between grains is then i. In practice for our simulations (described in more detail later), we generally randomly perturb the “corners” of all grain units (the black circles labeled with x, y positions in Fig. 1), so that the resulting grain pattern has grains of irregular (but quadrilateral) shape and some range of size (which is still on order of $α g × g$⁠). However, for computing various quantities, we find it is often sufficient to take the non-perturbed representation, thereby greatly simplifying computation (and various practical aspects of setting up simulations).

In 2D, grains (and interstitials) are columns; they have an infinite extent in the direction orthogonal to the plane. This extent is identical for the grain and interstitial portions, so we simply drop it everywhere when writing volumes. We have then that

Using these, we may calculate an average ablator density

Figure 2 shows a comparison of Eq. (3) against a set of measurement pairs¹⁴ of mean HDC density and grain size. In the present work, we take i = 4 nm to be the “base truth” (Fallon and Brown¹⁵ measure carbon grain boundary thickness in polycrystalline diamond $i ≲ 5$ nm). Taking this value for i means the 2D model does a reasonable job at reproducing the experimental mean density for a range of $α ≥ 1$ without changing i. In HDC for spherical ablators, it is often the case that grains are elongated in the radial direction, at least in the micrometer-scale grain cases where grain shape is most easily observed.

At present, we will focus on cases with g > 100 nm, in which case all the curves in Fig. 2 are reasonably close. Of course, grains are 3D in reality; for α = 1, we find $i ∼ 1.7$ nm yields a curve in Fig. 2 that is comparable with the $α = 1 , i = 2.5$ nm and $α = 4 , i = 4$ nm cases. The present approximations to the grain geometry are discussed further in Sec. VI.

As discussed in Sec. I, it is, in general, extremely computationally expensive, or computationally infeasible, to simulate an ICF implosion with well-resolved $∼ 4$-nm interstitials, even in 2D. This motivates the use of an alternate grain representation, which aims to replicate the impact of (fine-scale) grains at a coarser resolution. Once we are forced (by resolution considerations) to move away from our conceptualization of the true grains, it is not so clear how we should choose to represent grains. We may choose a representation that is rather different from our present conceptualization of the underlying true grains (e.g., the microstructure of Ref. 35).

Here we consider a class of simple scaling strategies, whereby the interstitial width i in the grain model described in Sec. II is (artificially) widened (⁠ $i → i ̃ , i ̃ > i$⁠). Scaling to a wider interstitial $i ̃$ reduces the resolution requirements, since resolving the thin interstitial width is the most stringent resolution requirement. A nice feature of the present scaling strategies is that they are convergent to our true grain representation as the simulation resolution (and therefore interstitial width $i ̃$⁠) is refined toward i.

We now define the class of scaling strategies, working from the grain representation of Sec. II. A given true grain pattern is determined by a combination of the parameters i, g, α, ρ_g, and ρ_i, plus a prescription for randomly perturbing the corners of the grains. A simple way to get a convergent representation is to fix the number and sizes of grains, which corresponds to fixing the corners of each grain (with the “lower left” corner labeled by x_j, y_k, Fig. 1), that is, we fix $( x j , y k )$ for all j and k. This means that the lengths $L x = α g + i = x j + 1 − x j$ and $L y = g + i = y k + 1 − y k$ are fixed, such that on rescaling $i → i ̃$⁠, we have $L ̃ x = L x$ and likewise $L ̃ y = L y$⁠. Throughout, we will specify the parameters in the scaled system with a tilde (as in $i ̃$⁠).

Then, the conditions $L ̃ x = L x$ and $L ̃ y = L y$ give, simply, that

We are now left to specify the grain and interstitial densities of the scaled system, $ρ ̃ g$ and $ρ ̃ i$⁠. In addition to those quantities already kept fixed, we would also like the mean density of the material (ablator), Eq. (3), to be conserved under changes in $i ̃$⁠, that is,

With this equation, and two unknowns, there are then an infinite number of possible pairs (⁠ $ρ ̃ g , ρ ̃ i$⁠) for a given $i ̃$⁠. This is the class of scaling strategies we draw from here.

Figure 3 shows an example of the space of scaling strategies for an underlying true case with i = 4 nm, g = 1 μm, and α = 2 (with $ρ g = 3.5006$ and $ρ i = 2.15$⁠). The black dashed lines represent, for different given $i ̃$⁠, the relation Eq. (6), and therefore any scaled grain representation that conserves mean density should, when scaled to a specific $i ̃$⁠, fall on the associated curve.

Also shown in Fig. 3 are the densities selected by two different prescriptions for choosing (⁠ $ρ ̃ g , ρ ̃ i$⁠). We dub these strategies $i Δ ρ$⁠, which aims to conserve the mass deficit in each interstitial, and $A ¯$⁠, which aims to conserve an average grain Atwood number. Sections III A and III B describe these two strategies.

Conceptually, this prescription for choosing $ρ ̃ g , ρ ̃ i$ aims to keep the mass deficit (per unit length) represented by the interstitial fixed as it is scaled. For the base case, this mass deficit is defined $i ( ρ g − ρ i ) = i Δ ρ$⁠. Since the deficit is calculated with regard to ρ_g, this approach logically regards ρ_g (⁠ $= 3.5006$⁠) as the true density in the limit of vanishing interstitials, whereby the ablator would be uniform. Fixing $ρ ̃ g = ρ g$⁠, we require that the quantity $i Δ ρ$ be fixed as the interstitials are scaled (⁠ $i Δ ρ = i ̃ Δ ρ ̃$⁠), which means that

Strictly speaking, once we fix $ρ ̃ g = ρ g$⁠, Eq. (6) should be solved to give $ρ ̃ i$⁠, rather than Eq. (7). However, the latter comes from an intuitive definition, and one can show that the error incurred in using Eq. (7) instead of Eq. (6) is $i ( i ̃ − i ) ( ρ g − ρ i ) / [ ( g + i ) ( α g + i ) ]$⁠, a small quantity for $i ≪ g$⁠.

Indeed, the example in Fig. 3 uses Eq. (7) for the $i Δ ρ$ curve, but Eq. (6) for the dotted black curves. We see very good agreement in this case. This $i Δ ρ$ strategy results in raising the interstitial density to be closer to the fixed grain density as $i ̃$ grows.

We can define an average Atwood number, relevant for shock interaction with a medium of nonuniform density,^45,46 as

The density dispersion $σ ρ = ( ρ 2 ¯ − ρ ¯ 2 ) 1 / 2$⁠, with the overline representing a spatial average here. In the present grain model, with two regions of volume V_g and V_i, at densities ρ_g and ρ_i, we have

where we have assumed that $ρ g ≥ ρ i$⁠. Equation (9) can then be expressed in terms of the quantities g, i, and α using Eqs. (1) and (2).

In the $A ¯$ strategy, we seek to simultaneously preserve this average Atwood number and the mean density, solving the system $A ¯ ̃= A ¯, ρ ¯ ̃ = ρ ¯$ for $ρ ̃ g , ρ ̃ i$⁠. This system can be solved analytically for the present grain model, but in practice, here we solve numerically to find the new densities for a given $i ̃$ and base case.

Figure 3 shows, again for an example case, the resulting density pairs for various values of $i ̃$⁠. We see that preserving $A ¯$ requires increasing both ρ_g and ρ_i as $i ̃$ increases; this is in contrast to the $i Δ ρ$ strategy when ρ_g is fixed.

While it is not feasible to well resolve grains with $i ∼ 4$ nm interstitials in an ICF implosion simulation, we can do so in a reduced problem, which is the focus of this section. This reduced problem studies the effects on post-shock perturbations when a shock is passed through a field of grains which are scaled according to the strategies outlined in Secs. III A and III B. First, we describe the setup and then we display the results.

Panel (a) in Figs. 4 and 5 shows the basic setup of our single-shock test problem, for a microscale grain case and a nanoscale grain case, respectively. Here cases with $i ̃ = 32$ are shown to make the features more visible than in the i = 4 case, while still showing the full domain. The $A ¯$ scaling strategy has been used. The base (true interstitial width) microscale case has g = 1 μm, i = 4 nm, and α = 2, as in the example shown in Fig. 3 illustrating the grain and interstitial densities under scaling of i. The base nanoscale case has g = 200 nm, i = 4 nm, and α = 1. We use these two base cases both in Sec. IV and later in ICF implosion simulations in Sec. V, so that these are the base cases when we refer to the “microscale case” or “nanoscale case.”

We simulate this setup in the radiation-hydrodynamics code HYDRA.^47,48 A $∼ 10$ (x) × $∼ 5$ μm² (y) region of grains is initialized, followed by a $∼ 1$ μm shocked region at a density $∼ 2 ρ g ≈ 7$ g/cm³. This density jump is chosen to correspond to a shock somewhat exceeding 12 Mbar, following the current strategy for reducing the possible impact of incomplete HDC melt.^10,16,17 The (post-shock) velocity and temperature are chosen then to match from the shock Hugoniot for HDC; here we use equation of state data for HDC from LLNL's LEOS library.^49,50 The pre-shock material has a uniform (cryogenic) temperature $∼ 12$ K.

The boundary conditions at the right-hand-side (RHS), $x = x max$⁠, of the domain are set to maintain these initial post-shock conditions (i.e., mass inflow at the appropriate velocity, density, and temperature). The boundary conditions at the LHS (⁠ $x = x min$⁠) of the domain do not come into play, since the material at that boundary is stationary through to the final simulation time. The top (⁠ $y = y max$⁠) and bottom (⁠ $y = y min$⁠) boundaries are slip (zero normal velocity).

Conceptually, the grain corners are first initialized as a uniform grid, as in Fig. 1, so that a given domain length $L x$ contains $N x = L x / L x$ grains, and similarly for the y direction. This grid of x_j, y_k grain corners is then perturbed by random uniform values $δ x j , δ y k$⁠, with $δ x j$ a random number drawn uniformly from the interval $[ − 0.4 L x , 0.4 L x ]$⁠. The $δ y k$ are likewise drawn from the interval $[ − 0.4 L y , 0.4 L y ]$⁠. Grain corners on the edges of the domain, or at the edge of the grain region, are only perturbed in the direction that is edge-parallel, so that they remain on the edge.

With the grain corners set in this fashion, the density in the interstitials, which are taken to be rectangles of width i (⁠ $i ̃$⁠) and extent from grain corner to grain corner, is then set using HYDRA's shape library.⁵¹ When a simulation cell is fully within the interstitial, its density is set to ρ_i (⁠ $ρ ̃ i$⁠). When a simulation cell is cut by the interstitial edge, the cell density is set to the appropriate volume-weighted average density between ρ_i and ρ_g (⁠ $ρ ̃ i , ρ ̃ g$ in the scaled cases). See Fig. 6 for an illustrative example of the results of this from one of the simulation cases. This feature assists with directly utilizing the analytic expressions of Secs. II and III to set grain and interstitial densities which approximately achieve the desired aim (conserving $ρ ¯$ and either $i Δ ρ$ or $A ¯$⁠).

For the single-shock problem discussed in this section, we run purely hydrodynamic simulations (no radiation), with an Eulerian (fixed), square grid. We note though that we have analyzed some similar simulations using a (moving) arbitrary Lagrange Eulerian (ALE) mesh, without obvious differences; the ALE simulations typically require a fair amount of relaxation to prevent mesh tangling from the grain-induced flows. The resolution of a given simulation run is specified in terms of the number of “cells per interstitial,” for example, if i = 4 and there are four cells per interstitial, then $d x = d y = i / 4 = 1$ nm. In all cases, we run with thermal conduction following Lee and More.⁵² The simulations use artificial viscosity to capture the shock, and in the results displayed here, we neglect physical viscosity. We have also run simulations including physical viscosity, which show only small differences from the present results. See Sec. VI for further discussion on these points. We note that thermal conduction does make a substantial difference to the perturbations generated by the grains. Its inclusion both somewhat decreases the post-shock velocity perturbations, and, more substantially, smooths the post-shock density; these impacts will be discussed in Sec. IV B.

With the simulations thus initialized, the shock progresses from right to left (in the negative x direction). After 0.37 ns, the shock has passed through nearly the entirety of the $∼ 10$ μm grain length. At this point, we analyze quantities in a narrow post-shock region of x extent $∼ 410$ nm. This region is indicated by vertical dashed lines in panels (b), (c), and (d) in Figs. 4 and 5, which show post-shock density and velocity fields for, respectively, the microscale and nanoscale grain cases scaled to $i ̃ = 32$ nm, as in panel (a).

In the post-shock, this averaging region averages over a few grains in the nanoscale grain case, but only a partial grain in the microscale case. While, due to time evolution and (numerical) dissipation in the downstream region, we should expect quantitative results of any analysis quantities to possibly depend (with very modest sensitivity) on the placement (and width) of this post-shock averaging region, our interest here is primarily in relative changes of quantities between different scalings of i. Further, it is the behavior in the immediate post-shock that sets the initial conditions for the downstream (post-shock) evolution. Note that different cases are all analyzed after an identical propagation time (and, to a very good approximation, propagation distance) for the shock; the apparent differences in Figs. 4 and 5 are due to slight differences in the total length of the grain region.

In the present work, we focus on two closely related quantities. First, we study the post-shock velocity dispersion σ_v. With the velocity deviations $v x , y ′ = v x , y − v ¯ x , y$⁠, this velocity dispersion is defined

with v_x and v_y the x and y velocity. Here it is understood that spatial averages (overline) are taken over the post-shock averaging region just described.

Second, we study the (average) perturbed kinetic energy density, which we call $E ′$⁠, in the frame moving with the average post-shock flow. As will be discussed in Secs. IV B and VI, in the present cases, $E ′$ and σ_v capture nearly identical information owing to the relative uniformity of the post-shock density field. We define

Before covering the analysis of the post-shock velocity dispersion σ_v, we examine the convergence of this quantity in the post-shock as a function of resolution. Figure 7 shows σ_v in a post-shock analysis region as a function of resolution, measured by cells per interstitial. These simulations use a setup identical to those described in Sec. IV A and discussed in Sec. IV B but a smaller domain. This is done because the smaller domain eases the use of many cells per interstitial (high resolution).

The case analyzed is the nanoscale grain case (g = 200 nm, α = 1), and the x extent of the post-shock averaging region is identical to that used elsewhere. Since the domain is smaller, fewer grains are averaged over in the y direction. Also, the shock propagates through a shorter length (⁠ $∼ 1$ μm) of grains before analysis.

Our main purpose here is to show that post-shock quantities like σ_v can converge somewhat slowly in resolution, but at four cells per interstitial, σ_v is nonetheless within about 75% of its estimated converged value (in this example). Here we only consider convergence within the hydrodynamic model presently employed; since the scales considered are rather small, especially at high cells per interstitial, physics outside of the present model could, in reality, be important. See also Sec. VI.

In all simulations that follow, we use four cells per interstitial. Thus, a simulation with, say, i = 256, has $256 / 4 = 64$ times fewer cells in each direction than the base i = 4 case, which we expect to yield roughly a $64 3 ∼ 10 5$ reduction in computational work in 2D when factoring in relaxed time step requirements (assuming the time step is determined by the hydrodynamics).

As an additional convergence-related comment, we note that the quantitative results presented here and later in Figs. 8 and 9 in Sec. IV B will vary somewhat if the random seed used to set the grain field is changed. This is because the analysis region (indicated in Figs. 4 and 5) gives us a finite sampling of the post-shock perturbations; we should expect larger seed-to-seed variation in the microscale case than the nanoscale case, since fewer grain-interstitial boundaries are contained in the analysis region for the former. As a concrete example, we have tested three different random seedings for a $i ̃ = 64$ nm, $A ¯$ scaling case analyzed shortly in Sec. IV B, for both nanoscale and microscale grains. In the nanoscale case, we find the difference between the maximum and minimum σ_v over the three seeds is $∼ 1.4$% of the mean σ_v, while this difference is $∼ 14$% for the microscale case.

The pre-shock density nonuniformity caused by the grains, and visible in panel (a) for the example cases shown in Figs. 4 and 5, leads to post-shock velocity perturbations, shown in panels (c) and (d). Most of these post-shock velocity perturbations are due to vortical motions, which are concentrated along the grain-interstitial boundaries. Some perturbed velocity is also associated with sound waves, which are also apparent in these panels and which travel largely in the direction opposite the shock propagation. These sound waves originate from shock reflection at grain-interstitial boundaries. In the case of sufficiently small pre-shock density perturbations, this problem can be analyzed using linear theory;¹⁹ this is the subject of separate work⁵³ (the present cases are not in the linear regime). We note that the impact of sound waves generated from certain types of internal shell defects has been studied previously,⁵⁴ and the present reflected sound waves may interact at the ablation front, an aspect which we do not study here.

The post-shock density, shown in panel (b) of these same figures, is, in general, much more uniform than the pre-shock density. For example, in the true i = 4 nm case with nanoscale grains (g = 200 nm, α = 1), $σ ρ / ρ ¯ ∼ 0.07$⁠, while in the post shock (labeled with the subscript 1), we find $σ ρ , 1 / ρ ¯ 1 ∼ 0.007$⁠, approximately a factor of 10 smaller. We observe similar magnitude drops (from similar or smaller initial magnitude) in all other cases shown in this work. This relative uniformity of the post-shock density field is driven in part by thermal conduction; Fig. 10 shows, for a small region of the aforementioned i = 4 nm case, a comparison of the post-shock density with and without thermal conduction enabled. While on the topic of thermal conduction, we also note that the post-shock velocity perturbations are modestly reduced from their “no-conduction” values when thermal conduction is included; this is shown in Fig. 7 where the purple plus point (no conduction) has a $∼ 12$% higher average velocity dispersion than the comparable simulation with conduction (blue star).

Thus, to a good approximation, the post-shock velocity dispersion σ_v is directly related to the post-shock kinetic energy perturbation deposited by the shock on traveling through the grains, $E ′ ∼ ρ ¯ 1 σ v 2$ (the impact of any density–velocity correlations being small owing to the relative uniformity of the post-shock density variations compared with the mean density). Here we (again) define this perturbed energy in the post-shock frame, and $ρ ¯ 1$ is the post-shock mean density. Previous work has shown that, for different cases with large pre-shock density perturbation (⁠ $σ ρ / ρ ¯ ∼ 1$⁠), the post-shock velocity dispersion is proportional to the average pre-shock Atwood number

in line with the expectations from analysis of the Richtmyer–Meshkov instability.^21,45,46 Here C is a coefficient that depends on the shock velocity and that may depend on various geometric aspects of the pre-shock density perturbations. For small $σ ρ / ρ ¯ , A ¯ ∼ σ ρ / ρ ¯$ and the proportionality between σ_v and $σ ρ / ρ ¯$ implied by Eq. (12) will also come out of other linear analyses of shock interaction with nonuniform density,^20,55 with these approaches also providing predictions for the shock velocity dependence and numerical coefficient, assuming a random density field and an ideal gas (and linear dynamics without thermal conduction).

The four panels of Fig. 8 show that this proportionality also holds to a reasonable approximation for the current cases. Panel (a) shows a plot of σ_v vs $A ¯$ for various $i ̃$ for the nanoscale grain case. Here the densities $ρ ̃ g , ρ ̃ i$ are selected using the fixed $i Δ ρ$ strategy described in Sec. III A. We see that $A ¯$ decreases with increasing i in this strategy, and that the resulting σ_v fall approximately on a line, as expected from Eq. (12). The i = 4 case falls slightly lower in σ_v than expected from the trend; we associate this primarily with a faster evolution timescale for the perturbations associated with these narrow interstitials, which then can (numerically) dissipate somewhat on the timescale of the analysis. Note that the decrease in $A ¯$ in the $i Δ ρ$ strategy can be solved for using the expression for $ρ ̃ i$ in this case, Eq. (7), in $A ¯$⁠, Eq. (8).

Panel (b) of Fig. 8 shows, also for the nanoscale grain case, the σ_v that result from following the fixed $A ¯$ strategy described in Sec. III B. First, we observe that solving the system $A ¯ ̃= A ¯, ρ ¯ ̃ = ρ ¯$ analytically for the (unpertubed) grain model is quite successful at preserving $A ¯$ from the i = 4 case as i is scaled. Second, we see that σ_v is, therefore, largely preserved as expected from the relation Eq. (12), although again the i = 4 case is slightly below the expectation.

Panels (c) and (d) show the same things, respectively, as panels (a) and (b), but for the microscale grain case. Since this case has a lower spatial density of interstitials, the values of $A ¯$ are smaller for a given interstitial width i (as expected from the analytic grain representation). Otherwise we see the same features: the linear relation of Eq. (12) holds to a good degree, with the i = 4 case falling slightly below the trend, and solving the system $A ¯ ̃= A ¯, ρ ¯ ̃ = ρ ¯$ nearly preserves $A ¯$ from the i = 4 case. For the cases shown in panels (a)–(d), the post-shock mean velocity is $v ¯ x ∼ 13.3$ km/s, so that $σ v / v ¯ x ≲ 4 %$⁠.

The maximum useful scaling of the interstitial width i, from the perspective of reducing the resolution requirements (at fixed cells per interstitial), is $i ̃ = ( g + i ) / 2$⁠, since this yields $g ̃ = ( g + i ) / 2$⁠. That is, the scaled grain and interstitial have the same linear width and require an equal number of cells (in the shorter dimension, assuming $α ≥ 1$⁠).

The coefficient C which encodes geometric information is apparently quite similar in the nanoscale and microscale cases and relatively insensitive to the scaling of the interstitial $i ̃$ within a given case (since the σ_v vs $A ¯$ show a quite linear relation). Although it is not easy to distinguish when comparing panel (a) with panel (c) in Fig. 8, a comparison shows C_micro > C_nano (modestly) for the present cases. Because larger velocity perturbations are associated with interstitials that are perpendicular to the shock front rather than parallel to it, we suspect this difference is a result of the fact that $α micro > α nano$⁠, which, for a given $A ¯$⁠, leads to a higher proportion of horizontal interstitial and therefore (modestly) larger σ_v.

This preceding discussion highlights a feature of the present scaling approach to relaxing the grain resolution requirements: it (to a degree) preserves certain geometric features of the grains and their associated perturbations. Of course, the full grain geometry is not preserved, and we should expect, in general, that, even in the case where an integral quantity like σ_v is approximately conserved, other quantities, even integral ones (e.g., the components $σ v x , σ v y$⁠), may change under scaling.

While the scaling approach for our present simplified grain model keeps certain geometric features, the scaling strategies (⁠ $i Δ ρ$ and $A ¯$⁠) are in principle agnostic of the underlying grain conception (or more generally, density perturbation conception). Assuming two pre-shock densities, ρ_g and ρ_i, arranged in some layout with volumes V_g and V_i, the $i Δ ρ$ strategy fixes $ρ ̃ g = ρ g$ and scales $ρ i → ρ ̃ i$ so that $ρ ¯$ is preserved as the relative volumes are changed. The $A ¯$ strategy instead solves for $ρ ̃ g$ and $ρ ̃ i$ to preserve $ρ ¯$ and $A ¯$ as the relative volumes are changed. Extension of this $A ¯$ strategy to other cases is straightforward, though additional constraints may be needed if extra variables are introduced.

Considered in this (grain agnostic) context, it is interesting to draw from the linear theory for a shock passing through pre-shock nonuniformity.^19,20,55,56 Suppose that the pre-shock density nonuniformity is a random field that is locally homogeneous and isotropic, and that the density perturbations are small, $σ ρ / ρ ¯ ≪ 1$ (in which case $A ¯ ≈ σ ρ / ρ ¯$⁠), as in Velikovich, Huete, and Wouchuk.²⁰ Then, as shown in that work, it is possible to solve (to lowest non-vanishing order) the shock interaction problem, providing not only the downstream perturbations but also corrections to the mean jumps in the shock Hugoniot relation (see Davidovits et al.⁴⁶ for an approach to approximate jump corrections in a large pre-shock density perturbation case). In this weak random field case, these corrections depend only on $σ ρ / ρ ¯$ and complex functions of the shock strength (Mach number M) and adiabatic exponent γ. Thus, if the grains were such a random field, and while scaling the interstitials we preserved $σ ρ / ρ ¯$ (as in the $A ¯$ strategy), then we would expect not only that σ_v is preserved but also other key quantities like corrections to the shock jumps and shock velocity. It is beyond the scope of the present work to study the degree to which this is the case for the (not truly random) grains and the scaling strategies discussed here. Rather, we move to study the impact of grains and scaling in an integrated ICF implosion context, where the grain-containing material (ablator) is subject to multiple shocks and a complex acceleration history.

We now turn to studying the impact of grains (in our simplified grain model) on an ICF implosion, with a particular focus on the qualitative and quantitative mixing of the HDC ablator into the deuterium–tritium (DT) fuel. As just highlighted at the end of Sec. IV B, the presence of grains could in principle perturb the implosion in other ways (e.g., shock velocity, initial HDC compression). Since resolving small interstitials (⁠ $i ∼ 4$ nm) in a capsule simulation is extremely computationally challenging, we resort to scaling the interstitials wider. Thus, although we ultimately have an interest in grain impacts on ICF implosions beyond the mix measures presented here, the present work aims to examine the impact of this grain scaling on the mixing dynamics and thereby give a reference point and assessment for the uncertainty involved in such grain modeling. However, we will also highlight complex physics features of the mixing dynamics with general relevance for ICF implosions. Before showing these results in Sec. V B, we first construct a hydrodynamic-only surrogate for an implosion (Sec. V A), which greatly reduces the expense of grain simulations.

Figure 11 gives an overview of the construction of a 2D hydrodynamic simulation which is matched to a radiation-hydrodynamic simulation. Here we use as our source the NIF experiment N170601,⁵⁷ a cryogenic layered shot with an HDC shell. The initial inner radius of the ice as well as the thicknesses of the DT ice and HDC ablator are marked in Fig. 11. We simulate the radiation-driven implosion of N170601 in 1D using HYDRA.^47,48 Using information from this simulation, we construct a simulation on a 2D domain that can reasonably replicate all key features, but that can be run without any radiation physics. Note that for our purposes, it is enough to have such a case reproduce the main features of a general ICF implosion, rather than exactly reproduce the source implosion (N170601).

From the full (radiation-driven) simulation, we extract the velocity vs time information for the point in the HDC that is initially $∼ 5.4$ μm from the HDC–DT interface. This position forms the outer boundary of our reduced-domain hydrodynamic-only implosion, r_max. We also extract velocity vs time information for the point in the DT that is initially $∼ 47.7$ μm from the HDC–DT interface. This forms the inner boundary condition of our hydrodynamic-only implosion, r_min. Figure 11 shows this reduced simulation domain, while Fig. 12 shows $r max ( t ) , r min ( t )$ and the radial velocities of these points, $v r , max = r ̇ max , v r , min = r ̇ min$⁠. Here the overdot indicates a time derivative.

The 2D hydrodynamic simulations are conducted on a $δ θ = 1$ degree wedge extending from the waist of the capsule and make use of HYDRA's arbitrary Lagrange Eulerian (ALE) features. This allows the radii of the outer and inner boundaries to track exactly the radiation-driven case once their velocities are specified as functions of time. At the left and right boundaries (constant θ boundaries), we use slip boundary conditions (zero normal flow) as in the simulations of Sec. IV. We also specify the density and temperature values at the inner and outer boundaries as functions of time, using values from the radiation-driven case.

Figure 13 compares the interface dynamics between the unperturbed (no grains) hydrodynamic-only simulation and the full 1D radiation-driven simulation. Panel (a) compares the HDC–DT interface radius vs time, $r int ( t )$⁠, while panel (b) compares the interface radial velocity vs time, $v r , int ( t ) = r ̇ int$⁠. We see that the hydrodynamic-only simulation matches the behavior of the radiation-including simulation very well.

For the perturbation and mixing dynamics at the HDC–DT interface, it is not sufficient to only match the interface time history. We must also have density and temperature profiles around the interface that follow the full simulation (or a general ICF implosion). To achieve this, we extract the electron-radiation coupling (joules/gram) as a function of time from the full simulations and add corresponding energy sources (or sinks, depending on the sign at each time) in the HDC and in the DT. At present, we find one source for the DT and one source for the HDC are sufficient. Although during stagnation the DT-radiation coupling is important, for the times we examine here, the primary radiation effect is to volumetrically heat the HDC. Figure 14 compares, for the hydrodynamic and full simulations, the density and temperature profiles vs radius at four times throughout the implosion. We see the hydrodynamic simulation tracks the full profiles quite well, though some more discrepancy emerges later in time (⁠ $t ≳ 7.47$ ns). Note the emergence, in panels (e) and (g), of steep local density gradients near the interface, which are due to thermal conduction heating the DT near the interface (lowering its density) and cooling the HDC near the interface (raising its density). These will be discussed further in Sec. V B.

Figure 15 shows a region of the initial density and mesh for one of the ICF implosion cases analyzed shortly in Sec. V B. As noted in Sec. IV A 1, we use an initial resolution of four cells across an interstitial, so that this $i ̃ = 32$ case has an (azimuthal) resolution of $32 / 4 = 8$ nm. For the ICF simulations, the initial radial zoning is finer than the azimuthal zoning, which is necessary for stability of the hydrodynamics with the ALE mesh; thus, there can be greater than four cells across an interstitial initially, depending on its orientation.

In Sec. IV, we examined the impact of grains (and scaling of grains) in terms of the post-shock velocity dispersion σ_v, which, when combined with the post-shock mean density, gives approximately the kinetic energy perturbation deposited when a shock passes through the field of grains. Here we consider, in the ICF implosion scenario described in Sec. V A, a rather different measure of the impact of grains. In particular, we consider the time evolution of a mix width measuring the interpenetration of the grain-containing HDC ablator with the DT fuel.

Let us now define the mix width used here. All fields in our 2D simulations are functions of $r , θ$⁠. Divide the radial extent into a set of (radial) bins, with edges, $[ R 1 , … , R j , … , R N ]$⁠. Let $r j = ( R j + R j + 1 ) / 2$⁠. The DT mass as a function of r_j, $m D T ( r j )$⁠, is computed by summing the mass within cells (c) within each bin, $m D T ( r j ) = ∑ r c ∈ [ R j , R j + 1 ] m D T , c$⁠, and likewise for the carbon (HDC) mass, $m C ( r j )$⁠. In practice to construct the bins we take a local average (in θ) radial resolution (the cells can be somewhat irregular owing to ALE), so that we calculate $m D T ( r j )$ and $m C ( r j )$ at (approximately) the grid resolution.

From these masses, we define the DT fraction as

Define the radial position, $r ξ$⁠, of a given DT mass fraction ξ to satisfy $f D T ( r ξ ) = ξ$ (where $ξ ∈ [ 0 , 1 ]$⁠). Then, the “5–95” mix width, $W 5 − 95$⁠, is defined

In principle, a given f_DT profile may be non-monotonic, such that are multiple points where $f D T ≈ 0.05$ or $f D T ≈ 0.95$⁠. For this reason, we take the maximum spatial difference of any such points. Since $W 5 − 95$ is a threshold-dependent measure, there can be some noise in this measure, but in practice, we find this is slight for the cases here. We have also studied a “full” mix width, $W 0 − 100$⁠, defined as in Eq. (14) but using $r > 0$ and $r < 1$⁠. In principle, $W 0 − 100$ captures the full extent of the mixing, while $W 5 − 95$ tracks more with the interface perturbation. We find the two mostly track each other here, with modest differences particularly early in time and so do not show $W 0 − 100$⁠.

Figure 16 shows $W 5 − 95$ vs time for the 2D hydrodynamic-only N170601 implosion, for different grain cases. Panels (a) and (b) show, for the nanoscale grain case, the result of scaling with the $i Δ ρ$ strategy and the $A ¯$ strategy, respectively. In both cases, the interstitials are refined from the maximum scaled value $i ̃ = 102$ nm, down to $i ̃ = 32$ nm. The mix width plots begin shortly before t = 4 ns, sometime after the first shock passes through the HDC–DT interface (see Fig. 13), and an RMI-like instability, seeded by the grains, has developed.

We observe that the mixing layer is re-compressed (decreases in width) three times after this initial shock; the first two recompressions come from the shocks clearly visible as jumps in $v r , int$ in Fig. 13, between 4 and 5 ns. The last re-compression is caused by more gradual acceleration of the interface after the last shock. It is important to note that during the early acceleration phase (here prior to $∼ 6$ ns), we initially have a stable configuration, with denser HDC pushing on less dense DT. The stability of this configuration during acceleration leads to this last mixing-layer re-compression.

Of particular note in panels (a) and (b) in Fig. 16 is that the $A ¯$ scaling strategy leads to $W 5 − 95$ experiencing only very modest changes as $i ̃$ is refined toward the true value of i = 4. As we will discuss further below, the mixing dynamics over the implosion are complex, in principle depending on the full time history of mixing. As such, we cannot rule out some substantial quantitative or qualitative changes as $i ̃$ is further refined. Nonetheless, in addition to possessing the potentially attractive properties highlighted in Sec. IV B, the $A ¯$ strategy also appears promising for estimating mixing impacts at lower computational cost.

It is also worth noting that the $i Δ ρ$ and $A ¯$ scaling strategies, in general, give different mix widths as a function of time for a given resolution. That is, even within this particular grain model, the mixing behavior varies depending on the de-resolution procedure. Recall that the initial conditions for the two strategies are identical at i = 4. For the cases examined here, at a given $i ̃$⁠, the $A ¯$ scaling leads to larger mix widths at all times than the $i Δ ρ$ strategy.

To further test these scalings, we also run 2D implosion simulations for the microscale grain case. The resulting mix widths for the $i Δ ρ$ and $A ¯$ strategies are shown in panels (c) and (d), respectively, of Fig. 16. In this microscale grain case, a wider range of $i ̃$ can be used (and an even greater speedup achieved). As a point of comparison, the $i ̃ = 32$ nm cases in the present work consume $O ( 10 5 )$ cpu-hours, while the $i ̃ = 502$ nm cases consume $O ( 10 2 )$ cpu-hours. Nonetheless, we again see relatively good agreement between the $i ̃ = 502$ case and $i ̃ = 32$ case in panel (d), using the $A ¯$ strategy. The same cannot be said for the $i Δ ρ$ case. Note that for wider (⁠ $i ̃ > 64$ nm) $i Δ ρ$ cases, there is no measurable mixing at the 5–95 level at early times; plots for a given $i ̃$ start when there is measurable mixing.

Figure 17 shows the ratio of $W 5 − 95$ computed for the least resolved case to $W 5 − 95$ computed for the most resolved case, for a given grain size (nanoscale or microscale) and scaling strategy (⁠ $i Δ ρ$ or $A ¯$⁠). That is, for the nanoscale grain cases, panels (a) and (b), it shows $W 5 − 95 ( t ; i ̃ = 102 ) / W 5 − 95 ( t ; i ̃ = 32 )$⁠, while for the microscale grain cases, it shows $W 5 − 95 ( t ; i ̃ = 502 ) / W 5 − 95 ( t ; i ̃ = 32 )$⁠. In panels (b) and (d), we see that the $A ¯$ strategy leads to the mix width ratio being within about $∼ [ 0.5 , 1.5 ]$ over the time history, while for the $i Δ ρ$ strategy, this range is $∼ [ 0.2 , 2 ]$⁠.

To get a sense of the degree of mixing implied by the 5–95 mix width values in Fig. 16, we can compute the penetration of the 95% DT (i.e., the 5% carbon) radius into the fuel, relative to the unperturbed (1D) fuel–ablator radius (r_int in Fig. 12). In particular, we can compute the fractional radial extent of DT fuel that is mixed at least at this 5% level. Working with the $i ̃ = 32$ microscale $A ¯$ case, we find that, through the shock and early acceleration phase (prior to $∼ 6$ ns), $∼ 10 %$ of the simulated DT radial extent is mixed at this level. This fraction grows through the acceleration phase, reaching 50% around t = 7.47 ns, when the hydro-only vs radiation-driven profile agreement is still quite good (Fig. 14) and continues to rise afterward. Since only $∼ 85 %$ of the initial DT fuel radial extent is included in the simulations (Fig. 11), these figures represent a (modest) overestimate.

We now discuss some of the qualitative features observable in the mixing dynamics, which have an impact on the time history of the quantitative mix measures. The mixing experiences a significant qualitative change during the acceleration-phase growth, observable in the density plots of Figs. 18–21. In particular, between the density snapshots of the second row in these figures (t = 5.67 ns) and the third row (t = 6.57 ns), there is instability associated with the development of thermal-conduction-driven local density gradients in both the DT and the HDC; these local density gradients are visible in unperturbed (1D) simulation profiles of Fig. 14, panels (e) and (g).

Since the early acceleration phase is stable, leading to the final observed decrease in the mix layer width, these local density gradients play a role especially in the onset of the final extended growth of the mixing layer from $t ∼ 6$ ns. Without these local density gradients (e.g., using a spatial filter to average them out, to mimic lower resolution), the profiles in panel (e) of Fig. 14, at t = 6.57 ns, would be marginally stable (interface Atwood number $∼ 0$⁠). However, the presence of these local density gradients makes for essentially two unstable regions, one in the DT and one in the HDC. In Fig. 16, we see the mix layer grows prior to t = 6.57 ns (this is true even controlling for Bell–Plesset effects^23,24). In the third row density plots in Figs. 18–21, the effect of the unstable density gradient in the ice is clearly visible. The physical arrangement on either side of the fuel–ablator interface during this part of the acceleration phase has similarities also to the typical Rayleigh–Benard setup,⁵⁸ particularly to that modified to use radiation.⁵⁹ The DT experiences heating from “below” (the HDC side), and there is a reservoir of colder DT above; a similar scenario occurs for the HDC, but the temperature-gradient scale length is much shorter.

The importance of thermal conduction for fine-scale fuel–ablator instability has been noted before, see, e.g., Hammel et al.⁶⁰ The grains give the potential for complex interaction between the time-history of mix and these local density gradients. For example, larger mixing driven by the shock-grain interaction (RMI-like mixing) at early times could lead to less build up of local density gradients due to thermal conduction, since the transition between HDC and DT is less sharp. This in turn could lead to reduced later-time RT-like growth, compared with a case with less early-time mixing. Since, during the acceleration phase, there will be a tendency to mix by RT until a stable density arrangement is reached, it is not so surprising to see a degree of “catch-up” in the mix widths of the $i Δ ρ$ cases to the $A ¯$ cases in Fig. 16. The former, after experiencing less shock-induced mixing during the initial three shocks, often experience steeper growth in $W$ during the acceleration phase, although this steeper growth is not generally sufficient to fully catch up to the more-mixed $A ¯$ cases before peak velocity (just before $∼ 8$ ns).

Figures 18–21 show, through density plots, various comparisons of the scaled microscale and nanoscale grain cases of Fig. 16. Note that, in all these figures, the first row is “zoomed in” and only displays a partial extent of the domain in both x and $y ∼ r$⁠, while the second and third rows show the full extent of the domain in x (and partial extent in y). This is done to improve visibility of the initial grain structures and the interface dynamics. We call attention to this because it may be tempting to compare wavelengths of features between the frames. In reality, there are $∼ 16$ grains in x for the microscale case and $∼ 80$ grains in x for the nanoscale case. In general, the wavelengths of the apparent dominant features are longer scale than the grains. This is particularly apparent in Fig. 21, which shows a comparison of the most-resolved nanoscale and microscale grain cases.

Also apparent from Fig. 21 (and observable quantitatively from Fig. 16), is that, for the $A ¯$ scaling cases at the same interstitial width, the microscale grain case leads to more mixing, including after the shock-dominated early phase, compare (b) and (e). Recall, from Sec. IV B, that the microscale case has less deposited perturbed kinetic energy. The larger mixing by t = 5.67 in the microscale grains is likely due to the long, coherent interstitials (likely making larger, more coherent local perturbations⁵³) as well as perhaps, the larger outer scale of the grains. For the most-resolved case shown (⁠ $i ̃ = 32$⁠), the microscale grains mix width is typically $∼ 20 − 60 %$ larger than the nanoscale case. However, comparing $W 5 − 95$ for the $i ̃ = 32$ microscale and nanoscale cases under the $i Δ ρ$ scaling, we find similar or somewhat lower mixing for the microscale case. Thus, while we anticipate the i = 4 (true) results would fall closer to the $A ¯$ results than the $i Δ ρ$ results, we do not regard as definitive the present suggestion of larger mixing in our microscale case.

Figure 18 shows a comparison of the density field for all the $A ¯$⁠, nanoscale grains, cases (⁠ $i ̃ = 32 , i ̃ = 64 , i ̃ = 102$⁠). Overall, the density behaviors are qualitatively and quantitatively quite similar, as expected from the quantitative comparison of $W$ in Fig. 16, though there is somewhat more apparent interface perturbation in panel (i), the $i ̃ = 102$ case.

Figure 19 shows a comparison, again of the density field, of the most and least resolved $A ¯$ microscale grain cases (⁠ $i ̃ = 32$ vs $i ̃ = 502$⁠). Despite being greatly de-resolved, the i = 502 case shows quite similar behavior at the peak mix width in the shock-driven phase, panel (e), compared with the i = 32 case, panel (b). As in the nanoscale case, the early acceleration-phase behavior (t = 6.57 ns), in panels (c) and (f), shows some observable difference in the interface behavior, but the quantitative difference is not so large as evidenced by Fig. 16.

Finally, in Fig. 20, we compare the $i ̃ = 502$ microscale grain case for the $i Δ ρ$ strategy to the $A ¯$ strategy. With this large scaling, $ρ ̃ i$ and $ρ ̃ g$⁠, while different for the $i Δ ρ$ case, do not yield a readily apparent grain pattern in panel (a); the grain pattern is the same as panel (d) despite not being visible. The perturbations induced in the passage of successive shocks are now small enough that no ablator (HDC) pieces get injected into the DT (or vice versa), after the three shocks, panel (b). The perturbations are now small enough during the early acceleration phase that the instability of the density gradient in the DT is greatly reduced, panel (c). This is also observable in panel (c) of Fig. 16, where the thin purple curve experiences relatively little growth around t = 6.57 in comparison with the more resolved cases. Nonetheless, the density profile will continue to grow more unstable with time and the interface will eventually mix more.

For microscale grains, the present $A ¯$ strategy gives a marginally feasible starting point to include a grain model in highly resolved 2D simulations. Here, at four cells per interstitial, the initial resolution of the $i ̃ = 502$ nm microscale grain case (2 × 1 μm grains) is $i ̃ / 4 ≈ 1 / 8$ μm. This resolution, for example, is comparable with the maximum refinement level used in various 2D implosions simulations utilizing adaptive mesh refinement (AMR),^33,35 though in the grain case, this fine resolution would be employed throughout the ablator (but one could perhaps include grains only near various interfaces as an additional approximation).

More generally though, the present grain approach is likely to be useful in wedge (partial sphere) simulations, to provide a faster benchmark by which one or another mix model can be calibrated for a given implosion. This calibrated mix model would then be deployed in further simulations (1D, 2D, or 3D) that study the given implosion. It is our hope that this grain scaling approach makes comprehensive studies of the impacts of grains in ablators more feasible, enabling both suites of 2D simulations and estimations of the impacts of grains in 3D (for other perturbation sources, it is known that 3D perturbations can reduce yield more than their 2D counterparts, e.g., Refs. 61–63).

The present approach to scaling grains is limited in its possible speedup by the grain size. A natural question is whether one might “superscale” grains, by relaxing the constraint on conservation of grain number (and size), while still preserving $ρ ¯$ (and, perhaps $A ¯ , i Δ ρ$⁠, or some other quantity). That is, one could consider combining grains together, reducing the number of interstitials and grains but increasing the size of both. We suspect that such an approach is possible (although perhaps even more approximate), but that it would require a different strategy for setting $ρ ̃ g$ and $ρ ̃ i$ than conserving $A ¯$ from the true case.

In considering such “superscaling,” it is useful to compare the nanoscale and microscale grains of the present work (Figs. 8, 16, and 21). The microscale case has smaller $A ¯$⁠, leading to less deposited kinetic energy but leads to similar or larger mix widths. If we got to a microscale grain case by combining nanoscale grains together, but preserved $A ¯$ and $ρ ¯$⁠, we should expect to end up with higher mix than the underlying true nanoscale case. Physically, in combining grains together, we expect to both increase the spatial scale of perturbations (that is, to deposit energy at comparatively longer lengthscales) and to increase the length of the interstitials, both of which may increase the level of mixing caused by the grains. Since there is substantial freedom in choosing $ρ ̃ g$ and $ρ ̃ i$ consistent with conserving $ρ ¯$ (Fig. 3), one could simply superscale by choosing these densities so that the resulting mixing behavior is closest to the underlying case. In other words, one could tune the superscaled case to match more resolved test cases; this might be particularly useful for developing a usable model for small (<100 nm) grains.

Prior simulation work, which included an examination of the impact of (grains) microstructure on the NIF experiment N161204,³⁵ suggests an approximately 25% decrement to the yield in that experiment from microstructure alone. Because, as already discussed, it is not practical to accurately resolve grain structures in such 2D simulations, this estimate used a strategy to represent grain-like perturbations at larger scale while preserving $ρ ¯$⁠. In general, we expect that moving perturbations to larger scale should create more disruptive perturbations. In contrast to the present scaling approaches, we also suspect this strategy results in increased $A ¯$ vs the base case, since regions with densities well below 2.15 (the nominal minimum density in the HDC assuming only grain-like perturbations) occur while $ρ ¯$ is preserved. This increased density contrast should tend to raise $σ ρ$⁠, while $ρ ¯$ was fixed, thereby raising $A ¯$ and likely increasing the kinetic energy perturbation deposited during shock passage. Zero-density regions (so-called voids) can, however, occur in HDC shells as a result of issues in the manufacturing process; these are a separate perturbation source from grains.

Ultimately, there are many uncertainties associated with such grain modeling, the present grain model included. Certain approximations or simplifications taken for the present study could be relaxed in a straightforward manner, while others would require substantial additional work or make the model computationally intractable. We discuss both now.

The present simulations take a particular approach to constructing an initial field of grains. It is clear that this approach is simplified in a number of ways. In general, it is possible for grains to have a range of sizes that spans different orders of magnitudes within the same sample.¹⁴ The present grains, even after being perturbed, are all on the same order of magnitude for size; this is particularly useful to help study whether and to what degree the (mean) grain size has an impact. The particular limits of our grain-corner perturbation (i.e., $δ ∈ [ − 0.4 L , 0.4 L ]$⁠, Sec. IV A), are somewhat arbitrary, apart from giving a spread of grain sizes while maintaining a similar size scale.

Real grains also have a varying number of sides (not always four as here), and grain boundaries (interstitials) that are not only straight lines (rectangles) as approximated here.¹⁴ It is possible that more complex grain and interstitial geometries would influence the comparison of mixing caused by microscale vs nanoscale grains. Real grains structures are, of course, also 3D. All of these features could be more closely matched in simulations. As we have already highlighted at the end of Sec. IV B, the present scaling strategies are agnostic to these underlying grain features. Maintenance of $A ¯$ can be applied generally, although one must keep in mind scale (size) effects as just highlighted in the discussion above. Prior work, in a different context, has also shown that the relation between σ_v and $A ¯$ is not necessarily strongly influenced by a move from 2D to 3D.⁴⁶ 

At a more detailed level, our choices to match the mean density, Fig. 2, introduce approximations. The density reduction with grain size is believed to be due not only to the presence of graphitic carbon but also to hydrogen.¹⁴ While this hydrogen also has a tendency to be associated with grain boundaries,⁶⁴ we have here assigned all density reduction to graphitic interstitials. In principle, one could imagine a grain model containing more than only two possible densities, which would also introduce additional free parameters once the grains are de-resolved (requiring additional constraints or assumptions).

In addition, we treat transitions between interstitial and grain as arbitrarily sharp, occurring over a distance determined by the cell sizing. In practice, we find that the total vorticity deposited along these transitions depends to a very good approximation only on the total density contrast between grain and interstitial, while the peak vorticity values do depend on the gradient as one might expect. At a resolution of four cells per interstitial, for 4-nm interstitials, the 1 nm cell size is a fair bit larger than, but approaching, the expected interatomic spacing (⁠ $∼ 0.15$ nm for C–C in diamond). At the highest resolution used during convergence testing in Fig. 7, the resolution is similar to the interatomic spacing. It is possible that a kinetic treatment would find corrections to the post-shock perturbations calculated by the present fluid simulations, motivating future work to investigate whether this is the case.

There are many number of ways various material properties could play a role in grains perturbations; we give a few examples now. While we include grains here as a nonuniformity of the density, we neglect any property heterogeneities or anisotropies (e.g., anisotropic sound speeds) associated with the HDC material itself due to its crystal structure. While current HDC designs aim for a first shock >12 Mbar to try to fully melt the HDC, so that these impacts are mitigated, there are uncertainties associated with refreeze of HDC material after the first shock. One observation from the single shock simulations of Sec. IV is that, while the post-shock pressure exceeds 12 Mbar in the unperturbed case (no grains), and on average in the perturbed case (with grains), the local pressure may be lower (e.g., $∼ 10$ Mbar) in the grain case owing to overlapping shock reflections from interstitials.

As mentioned in Sec. IV, the thermal conductivity of the HDC has an impact on the size of post-shock velocity perturbations calculated for a given (simplified) grain density field. Additionally, the thermal conductivity of HDC and DT throughout the implosion has an influence in setting the near-interface density gradients discussed in Sec. V, which in turn influence the fuel–ablator mixing dynamics. To the extent this conductivity is uncertain in these regimes, this introduces modeling uncertainty for the grains problem.

The present simulations do not include (mass) diffusion. From the density and temperature profiles in Fig. 14, the DT during the shock phase and early acceleration phase is generally $1 – 8$ g/cm³ and from a few eV to $∼ 20$ eV. While we expect relatively modest diffusive effects in this regime on the nanosecond timescale, strictly speaking diffusion should be included given the very high resolutions utilized in our simulations and the associated small feature sizes produced.

For the purposes of testing the grain modeling, we have used here a hydrodynamic-only implosion, as described in Sec. V A. As noted there, we implement an energy source in the HDC to match the radiative heating observed in simulations including radiation. Recall also that we focus on the innermost HDC, near the DT fuel layer; for the source implosion (N170601), this innermost HDC layer is surrounded by a W-doped HDC layer. In this source calculation, we find that radiative heating of the HDC ahead of the first shock is quite small, so that the pre-shock HDC is $≲ 0.2$ eV (usually a fair bit less than). This upper value is below the anticipated carbon melt threshold. As a matter of course, our simulations use a so-called quiet start routine to prevent premature motion of the grains structures, which might otherwise occur owing to the fact that the interstitials and grains have the same (cryogenic) starting temperature but different densities (so that they would not be in pressure equilibrium). This routine uses a temperature threshold, here set similar to the melt temperature, which is well below the post-shock temperature after the first shock (∼2 eV), but above the pre-heat temperature; thus we assume no change in the initial grain structure from radiative pre-heat in the present simulations, which we believe is reasonable given the simulated pre-heat temperatures. However, the possibility for such effects should be kept in mind, particularly depending on the HDC use case and operating regime.

As a result of the high computational cost associated with the grain-including ICF simulations of Sec. V, we have used a one-degree angular extent for our wedge simulations (Fig. 11). Through the times displayed in the density plots of Figs. 18–21, this angular extent is sufficient to capture multiple of the largest visible perturbation wavelengths. However, it is possible that the mix widths shown in Fig. 16 are influenced by this restricted angular extent, in particular at later times approaching or during stagnation (⁠ $∼ 8$ ns). As previously noted, the agreement of the hydrodynamic-only density and temperature profiles with the radiation-driven (true) profiles (Fig. 14) is less good at later times (⁠ $t ≳ 7.47$ ns), and any such profile drift may influence the mix widths at these later times. However, our focus here is on comparing scaling strategies and examining the qualitative features of grain-induced mixing through the shock phase and early acceleration phase.

Here we have developed a simplified model for grains, following a strategy that reproduces (approximately) the experimentally observed relation between HDC grain size and density. We have proposed two different strategies for implementing these grains at reduced resolution, both of which conserve the material mean density. As we show, even for our simplified model, this density constraint alone leaves significant freedom in the de-resolving (scaling) process, and the choices made with this freedom impact both the velocity dispersion (∼kinetic energy) deposited during shock interaction and the mixing dynamics at the fuel–ablator interface in an ICF implosion. At least for the present kinetic energy and mix metrics, we find that one of the two strategies leads to reasonable agreement even when de-resolving by over an order of magnitude (⁠ $i ̃ = 32 → i ̃ = 502$⁠), corresponding to a speedup in our 2D simulations of a factor $∼ 10 3$⁠. We must, however, still be cognizant that this most-resolved case is not our true base case (i = 4), and that the present grain model is simplified in a number of ways as discussed. Additionally, we show that grain-induced mixing depends on a complex interaction history between grain-driven perturbations and conduction-driven density gradients at the fuel–ablator interface. We observe (partial) catch-up of cases with smaller initial grain-induced mixing to cases with larger initial mixing, owing to less smoothing of local density gradients. Instability of such density gradients in the DT fuel, which is triggered early in the acceleration phase in most but not all cases by the grains perturbations, leads to significant penetration of carbon ablator into the fuel at this early time. While the setup used for the ICF simulations in this work limits our ability to study the spectrum of possible grains impacts, a comparison of our microscale and nanoscale grain cases suggests the microscale case may lead to more mixing. The reduced-resolution grain model developed here should enable additional grain modeling to further test this observation, as well to examine the full spectrum of grains impacts in integrated ICF simulations. Thus, we hope this work begins the process of reducing uncertainties around the impacts of grains on ICF designs using HDC ablators.

This work was performed under the auspices of the U.S. Department of Energy by the Lawrence Livermore National Laboratory under Contract No. DE-AC52-07NA27344. S.D. was supported by the LLNL-LDRD Program under Project No. 20-ERD-058. S.D. would like to acknowledge helpful discussions with Grace Li about linear analysis of the grains problem.

This document was prepared as an account of work sponsored by an agency of the U.S. government. Neither the U.S. government nor Lawrence Livermore National Security, LLC, nor any of their employees makes any warranty, expressed or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that its use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the U.S. government or Lawrence Livermore National Security, LLC. The views and opinions of authors expressed herein do not necessarily state or reflect those of the U.S. government or Lawrence Livermore National Security, LLC, and shall not be used for advertising or product endorsement purposes.

The authors have no conflicts to disclose.

Seth Davidovits: Conceptualization (lead); Formal analysis (lead); Investigation (lead); Methodology (lead); Software (lead); Writing – original draft (lead); Writing – review & editing (equal). Christopher Weber: Conceptualization (supporting); Methodology (supporting); Software (supporting); Supervision (equal); Writing – review & editing (equal). Daniel S. Clark: Conceptualization (supporting); Methodology (supporting); Supervision (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.