2-D Lagrangian studies of symmetry and stability of laser fusion targets

Computer Physics Communications 43(1986)107—124 North-Holland, Amsterdam

107

2-D LAGRANGIAN STUDIES OF SYMMETRY AND STABILITY OF LASER FUSION TARGETS Stefano ATZENI Assocazione EURA TOM— ENEA sulla Fusione, Centro Ricerche Energia Frascati, C.P. 65, 00044 Frascati, Roma, Italy

The use of 2-D Lagrangian codes for studying the symmetry and the stability of laser fusion targets is critically discussed. Physical models and their finite difference implementation are illustrated, with particular re(erence to the three-temperature code, DUED. Applications to “model problems” and to full-scale laser target experiments illustrate the potential of the approach.

1. Introduction In directly driven laser and particle-beam fusion (or, in general, in directly driven pellet Inertial Confinement Fusion (ICF)) a certain number of laser (or particle) beams irradiate a spherical target containing some amount of fusionable materials (a few mg of a D—T mixture in prospective reactor targets). As a reaction to the heating and ablation of the outer layers of the target, the non-ablated material is accelerated inward, i.e. implodes. The inner fuel is eventually compressed, and a certain portion of it is also substantially heated. If a sufficiently large and dense quantity of D—T fuel is heated at a temperature above about 4 to 5 keY, the fuel ignites (i.e. the thermonuclear reactions are self-sustained) and a thermonuclear burn-wave is propagated to the remaining fuel, which burns for a time of the order of its hydrodynamic confinement time [1—4]. In indirectly driven ICF schemes laser or partide beams are used instead to generate a suitable X-ray field, which, in turn, drives the target implosion [3,5]. It is usually stated [6—8]that currently considered ICF concepts require a quasi-spherical implosion, the actual degree of the allowed amount of asymmetry being the subject of research. Indeed, an asymmetric implosion could result in reduced compression of the fuel, and/or cause mixing of the fuel with the other materials of the target, and

could also hinder the development of the outward propagating burn-wave. Asymmetric implosions could be generated by both the macroscopic non-uniformities of the driving field (e.g. of the laser beam), and by the occurrence of plasma and fluid instabilities. Among the latter we mention the Rayleigh—Taylor instability of the ablation front and the thermal and ponderomotive self-focusing of the laser beam in the low-density plasma corona (see, for example, refs. [3,4]). The importance of studying these processes, which necessarily requires the use of 2-D or 3-D numerical codes, is further enhanced by recent trends in target design. Indeed, theoretical one-dimensional predictions [6,9,10] and recent experimental results [11] indicate the advantages of using targets consisting of very thin, multi-layered shells (with aspect ratio, i.e. radius to thickness ratio, A = 50 to 200), irradiated by moderate-intensity, short-wavelength laser pulses (with A ~ 0.5 ~im and IA2 < 1014 W cm2 ~m2). Owing to their small thickness and to the use of short laser wavelengths, these targets are likely to be very sensitive to non-uniform irradiation [3,4,6]. Indirectly driven targets also call for extensive numerical studies aimed at contributing to the understanding of the complex processes which lead to the generation of a quasi-uniform X-ray field inside a suitable laser- (or beam-) target configuration.

0010-4655/86/$03.50 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)

108

5,

A Izeni

/ Syrnrneg,-~’and stahiliii.

The numerical simulation of ICF experiments requires the description of a large number of basic processes concerning fluid-dynamics, plasma physics, dense-matter physics, Equation-Of-State (EOS), radiative transfer, nuclear reactions, etc. These processes have typical time and space scales varying over many orders of magnitude. The same applies to the mechanisms leading to the relevant target asymmetries and instabilities. For instance, while the most dangerous modes of the Rayleigh— Taylor instability have wavelengths of the order of a few times the thickness of the shell (or 10—100 ~m) [4,6], thermal self-focusing occurs over a length of the order of several mm [121. Also, thermonuclear burn should last less than 100 Ps, while the implosion may take up to several dozen ns. Most of the processes of interest for the previously quoted moderate-intensity, short-wavelength regime can be studied by means of the Lagrangian model of an unmagnetized fluid. Indeed, processes which require a microscopic or kinetic treatment (such as the generation and transport of superthermal electrons) seem to play a minor role in experiments performed in this regime [4,6,8]. The same should apply to magnetic field effects, whose treatment, included in the LASNEX Lagrangian ICF code [13] and in a few Eulerian codes [14], will not be discussed in this paper. The choice of the Lagrangian approach is supported by its ability to deal with fluid systems including regions with different typical scale-lengths, and which undergo strong compression and expansion. Lagrangian codes can also easily discriminate interfaces between different materials and between fluid and vacuum (i.e. free boundary problem) [15]. In this paper we discuss the main computational aspects of flexible and relatively fast-running 2-D Lagrangian fluid codes for ICF, and illustrate their potentialities by means of several examples. We describe physical models and numerical algorithms, which are chosen, as usual for complex multi-dimensional codes, in an attempt to attain a compromise between the accuracy of the models and the flexibility, robustness and execution speed of the numerical implementalions. In particular we discuss the consistency and conservation properties, on an arbitrarily dis-

of laser fusion tar,geo

torted, quadrilateral-zone Lagrangian mesh, of the finite-difference schemes for a couple of space operators (~‘and v-xc~) that enter the models of many different processes. The performance of the considered algorithms on vector computers is also briefly discussed. The basic equations, the mesh and the evolulion scheme of the simplest conceivable Lagrangian fluid codes are briefly illustrated in section 2. Such codes are of interest because they can be considered as a basic framework for codes including more complex physical models, as will become apparent from the material discussed later. Furthermore, despite their simplicity, they can already be used to study, in an approximate fashion, some crucial aspects of the stability of ICF targets [16]. The subsequent part of this paper is more specifically devoted to laser fusion codes, and starts, in section 3, with a description of the physical model of the DU ED three-temperature (3-T) thermonuclear code [17]. The finite-difference solution of the 3-T equations, which are parabolic partial differential equations, is treated in section 4. A numerical method for the ray-tracing treatment of laser—plasma interaction is outlined in section 5. The applications of 2-D Lagrangian ICF codes to schematized situations, i.e. to “model problems”, and to full-scale laser—target experiments are discussed in sections 6 and 7, respectively, with the help of examples of problems run by the DUED code. Some conclusions are drawn in section 8. Concluding this introductory section, it is worth while to recall that purely Lagrangian codes cannot handle strongly sheared flows, in the presence of which the mesh is substantially distorted, or may even become pathological (in the sense, for example, that opposite sides of a quadrilateral zone intersect each other [18]). However, despite this limitation, these codes can still be used to study a number of problems relevant to the symmetry and stability of laser fusion targets. If appropriately designed, they can also satisfactorily deal with situations in which some “small” highly sheared region is present [18]. The other often quoted drawback of Lagrangian codes, with respect to orthogonal-mesh Eulerian codes, is the involved treatment of some differential operators

S. Atzeni

/ Symmetry

and stability of laser fusion targets

(e.g. the diffusioin operator, (~7.XV)) on an arbitrarily distorted quadrilateral mesh. This problem has, in fact, been succesfully solved by the introduction of accurate, consistent and energy-conserving algorithms, and by the use of fast, partially vectorizable solvers for the solution of the resulting linear system (see section 4).

2. Basic features of 2-D Lagrangian codes In the Lagrangian description of a fluid, when the state (e.g. the density p and the specific energy er’) and the macroscopic velocity u are given at the time t = t0 in all points of a domain supple~,

mented by boundary conditions at t ~ t0 on the frontier of one aims at finding out, for t> t~, the position R, the velocity and the state of the generic fluid element, which at time t = t0 was at the point R0, namely one aims at finding out R = R(R0, t), u = u(R0, t), p = p(R0, t) and ~ =e~°(R0, t). ~,

R

109

2.]. The Lagrangian mesh Discretization of the fluid equations requires the introduction of a mesh (or net, or grid), defining a finite set of grid-points and of fluid elements (“zones”) on which the equations are solved with a finite-difference scheme. Several types of grid topology based on triangular, quadrilateral, or even exagonal zones, have been suggested, not to mention more complex, non-purely Lagrangian, reconnecting grids [19]. Here we only consider the quadrilateral zone grid, because of its simplicity and flexibility and the wide experience which has been accumulated in more than twenty years of routine use. In this description we shall essentially follow the original treatment of Schulz [18]. To simplify the notation and make discretization easier, we replace the pair of coordinates (R0, Z0) (the third coordinate being ignorable, whether it be Cartesian or the axial coordinate of a cylindrical system) with a pair of dimensionless Lagrangian coordinates (i, j), so that R(R0, t) is

R.

~ r________

‘‘~1rmTn11~~~—GHOST—ZONES

~~fl~:ILASER

i=0

...,

Z

i=J--l

(a)

(b)

Fig. 1. Mesh initialization: (a) for a cylindrical target irradiated by a laser beam (cylindrical coordinates); (b) for a hollow cylinder imploded by an outer pressure profile with a periodicity ~ = i~/6(Cartesian coordinates); an analogous mesh can be used for a sphere. Dashed lines respresent “ghost zones”.

110

S. A t:eni

/

Symmetry and stability of laser fusion targets

replaced by R(i, j, t) or, in shorthand notation, R.~(t).The physical (or Eulerian) coordinates are related to the logical (or Lagrangian) coordinates

specific power dissipated by the artificial viscosity. In tensor notation one has (ae~~/at)a = (au,,,/aR1)/p. The tensorial form for the artificial

through the Jacobian J = ~Pi(R,Z)/3(i, j). For typical ICF targets the definition of an initially orthogonal (or quasi-orthogonal) mesh is then immediate, as shown in, for example, fig. 1. The domain is then mapped onto the rectangle 1
stress was introduced to deal with shock waves propagagting obliquely to the mesh [18]. Eq. (4) is usually replaced by a temperature equation

~.

ç~

= —

(B + p)~— +

(~

a~ = (a~/aT)0 is

)

at

where the specific heat, B = (a~,/aP-’)T and p is now given as p(p, T). The standard finite-difference (FD) sequence of advancement of eqs. (1 )—(3) and (5) from a generic time level ~‘ to the next one, ~ui 4-1 = ‘ + ~ t~4-1/2, is obtained in a straightforward manner if the velocities are defined n+ 2= (t”at+ intermediate t’~’)/2,andlevels all other 1/2 (with t’~’~ variables at integer levels n, n + 1, etc.) ‘Employing time-centered FD, we may then write ~,

u”

1/2 = ~n

+

1/2

~

[—

( c’p)”

1 2

/ +p’~],

2.2. Simple finite-difference fluid model

1 =

To discuss the general features of 2-D Lagrangian fluid codes, we first consider the simpie system of the fluid equations governing the evolution of a compressible, inviscid fluid, subjected to gravity only. We then write

R”

=

T”

+

aR/at = a ~—(pJR)

=

0,

a

1

(3) /

a~~

= —p—-~---— + ~-~----)

,

(4)

where a/at indicate partial time derivatives in the Lagrangian frame (i.e. with i and j kept fixed), R = 1 in Cartesian coordinates and R = R in cylindrical coordinates, g is the gravity, p =p(p, e~) is the fluid pressure, ~a is the so-called artificial viscosity stress tensor (needed to treat the discontinuities introduced by shocks in the framework of a continuous model [15,18]), and (as/at)” is the

—(B +p)’~4-1z”

C, —1

x ~

(2)

~,

(8)

,

,,+1/2

(6)

(7)

z~t”4-1~2~~4-1/2,

Ri R”4- ‘i”4-’

~y’4-i =

T”4-’

+

—1

(1)

(5)

/ +

~

“4-

I—I at

)

“

1/2~

(9)

1

where space indexes have been omitted, it being understood that each equation is solved in all grid-points (eqs. (6) and (7)) or on all~zones (eqs. (8) and (9)) before going to the next equation. The Jacobian and the EOS parameters are computed after solving eqs. (7) and (9), respectively. The temperature equation and the EOS are usually iterated to time-center the values of p. C,, and B in eq. (9), thus improving the accuracy of the solution. It should be pointed out that the continuity

S. A tzeni

/

Symmetry and stability of laser fusion targets

equation in the form of eq. (8), with J computed as a zone-centered quantity, gives meaningful values of the density also for pathological zones [18]. The alternative form p”4-’ = pOAoR0/A~ik~~~1, where A is the zone area (computed from the elementary geometry expression) could instead lead to negative densities. Numerical stability of the solution of eqs. (6)—(9) is ensured if the time step t~tsatisfies, in any zone, the 2-D Courant—Friedrichs—Lewy condition [18], t~t~ ~L/2C (C is the sound speed and L~Lis a typical zone length) and a further condition, related to the artificial stress, which essentially imposes a maximum (around 10%) on the relative variation of the density of a zone in a time step [18]. 2.3. Solution of the momentum equation As far as the spatial operators are concerned, the momentum equation is first differenced by writing gradient and divergence operators in terms of the Lagrangian variables (e.g. V J’[(aR/ ai)a/aj (aR/aj)a/ai], with R us (Z, —R)), and is then discretized in the i—j space [18]. In this context it is worth discussing how to get through a characteristic difficulty of Lagrangian FD schemes for the momentum equation. To show this, let us consider the contribution of the pres—

R ~

/

1 [aR pJ j~8i

ap

8j

—

aR

a~

3j

3z J~J

Here the values of quantities such as [(aR/a~)(a~/aJ)/~.J],~must be obtained as ~1~i—1/2,j +

=

(1

(ij) -

/ —

________

j

=

i.e. as some average of the values at (i— 1/2, ~1) and (i + 1/2, j), which are easily available as [18] + 1/2

— —

j

(~R/~i)I+l/2.J.,/s+ (~R/~i),+,/2,1+,/2 (Jp)1±l/2,1~~l/2

+

(Jp)~±l/2,f÷l/2

2.j + 1/2 P 2,j — 1/2)~ (11) x (p, l/ approach 1 [18] l/ The standard consists in choosing —

±

±

—‘

it

-

________

—

j+l

j —a-

i-i

z

Fig. 2 Evaluation of 4’ =(p’rp) in point (i, J) 4’jJ ~ obtained as an average value over 4’,— l/2,~ and 4’, + 1/2.j’ which refer to the upper and lower couple of zones, respectively (see main text).

the values of the coefficient ~jin eq. (10) in such a way as to ensure (approximate) energy conservation, which is attained if each contribution to is weighted proportionally to the average size of the two zones on which the pressure derivative (ap/aj in eq. (11)) is computed, in the Lagrangian direction orthogonal to that of differentiation (i in the case of eq. (11)). For instance, in the simple case shown in fig. 2, where for simplicity we consider rectangular zones, one would have ~ = (R 1 R,_,,1)/(R,+,,~ R,,~). This means that the contribution to from the most distant of the two points (i.e. point (i 1/2, j) has a larger —

—

weight, giving a poor approximation of the acceleration in (i, j)). In more formal terms, this is an inconsistent FD representation of the differential

(10)

— ~)41+1/2,J’

+1

-

—

sure gradient to the z-component of the acceleration of a generic grid-point (i, ])~which is

111

momentum Actually, we have instead found that itequation. is preferable to ensure a good approximation of the local acceleration, even if this is attained at the expense of local energy conservation. In the case of fig. 2 this is done by setting in eq. (10) ~ = 1 [(R,1 R,,1)/(R,~,1 R.,~)]. Indeed, when adjacent zone have corn—

—

—

parable mass and very different size (which is often the case around the ablation front of a laser fusion target), the latter procedure is needed to obtain qualitatively correct local results (i.e. proper movement which is crucial in the description of of grid-points, processes where physical “singularities” occur). The error in the global energy con-

112

S. A tzen,

/

Symmetrj’ and stability of laser fusion targets

ap’

servation introduced by this procedure seems to be practically negligible.

Cye

liT

3. A three-temperature laser fusion model

aT = Cvr~~

A relatively simple, but at the same time rather complete, fluid model for laser fusion targets is

where ~ Bfl=(a~/aP’)J.,~ is the specific energy of the fluid $ (with /3 = i, e, r), p1, is the pressure, f~is the thermal flux,

at

—~=

(13)

‘Qer+SFC+SL,

that implemented by the DUED code [17]. This essentially is a 3-T model [20], and also takes into account the fusion reactions, the fuel burn-up, the diffusion of charged fusion products, the refraction and the absorption of the laser-light, as well as the properties of the real materials (through an appropriate EOS). For many processes more accurate treatments than here exist 2-D (some being incorporated those in the given the few existing Lagrangian codes for laser fusion [13,21,22]); in particular we may quote the multigroup treatment of radiative transfer and the non-local thermodynamic equilibrium (non-LTE) treatment of the plasma corona. Nevertheless, the field of application of the present model is very wide, and allows most aspects of moderate-intensity laser fusion experiments to be studied. Furthermore, this model also serves to illustrate, in a general fashion, the typical computational aspects of any ICF Lagrangian model.

—

(Br + Pr)

~ap’

—

V

F~+ Qer’

(14)

is the specific power exchanged by the power fluids /3Q~ and y, SFe (SF~) is the fusion-specific delivered to the electrons (ions), and SL is the specific power delivered by the laser to the elec-

trons. The radiation temperature equation, eq. (14), is obtained from the radiative transfer equation on the assumption that the radiation field has4/c.a Planckian spectrum, for which P~r= 4a51 where ~s is Stefan’s constant and c is the speed of light. The radiative flux is approximated by using a flux-limited diffusive model, i.e. by writing Fr = Furthermore, the same energy exchange term is formally written in the form as that used for electron—ion energy exchange, i.e. Qer ~ ~

(To— TheI~)/Ter. three thermal fluxes are then written as ~

=

~(Xett)pV7~,

(15)

with (Xeff)$C1$(9is)Xti.

3.1. The 3-Tmodel

where In a 3-T model [20] the ions and the electrons of the plasma, as well as the radiation (i.e. the thermal photons), are considered as three separate

(16)

x~is the appropriate thermal conductivity.

a~(O~) is a factor which takes flux limitation into account, and = I x~ V7~/( F 1~) (17) ~,

fluids, each being in LTE, and therefore characterized by its own temperature, T1, 7~and 7’, respectively. Each fluid conducts heat, and the electrons exchange energy with both the ions and the radiation. The single-temperature equation of section 2 is then replaced by the following three temperature equations

) ap

aT. +p C~5—~ = -(B.

at

,

/ ~~\“

setting

—~-——+~-~_)

—V~I~+ Qei+

SF1,

i.e. the ratio of the flux given by the Fourier law to the relevant limiting flux. The introduction of the factor a~ is motivated by the fact that the Fourier law is not valid for steep temperature gradients. Indeed, in any case, causality requires I <(F,~)~= ~ where u13 is the typical velocity of the particles of the species /3. For the electrons such a limitation is then enforced by

(12)

1)

a~=min(1, f9

S. A tzeni

/ Symmetry and stability

or, alternatively, ae = (1 + Oe/f) where f is a numerical coefficient whose value should be in the range 0.03—0.6 [3,4]. For the radiation one has 4, and we set [23] (F,im)r=4as~ ~,

ar(Or)

of laser fusion targets

temperatures below about 30 keY. We then write aEa =

_____________________ 6+ (9/16)Or

6

=

+ (3/4)~r +

113

—~(~ +~

+ s~ (19)

+V•(Xeff)VE~

(9/16)O~

Tae~

where p,,, is the a-particle pressure, s,, is the Notice that in any numerical scheme the value of the factor a~ at a given time level is computed using either the value of O~ at the previous time level, or at the previous iteration at the same time level. The energy exchange terms are written as

(18)

~

a-particle power density source, (Xeff)a is the effective diffusion coefficient, ‘T,,e is the typical ato-electron energy relaxation time [28], and 1~is the fraction of energy delivered to the electrons (for which we write .F~= 33/(33 + 1~(keV))[20]). The specific powers 5F1 and SFe which enter eqs. (12) and (13), respectively, are E 1—F E~, 5Fe —

where C~ 0= (3/2)N0kZ*, Tp~ is the energy relaxation time, and Z” is the average ion charge. In the DUED code standard expressions are used for Xe’ and ~ [24]. For Xr and ~r simple ~

PT,,e

F

e;

The diffusion coefficient is X~Xae(1 +

2b1(v0) ~‘ be ( v~)j

analytical expressions are employed which take into account free—free transitions, Thompson scattering and, in a rough way, free—bound transitions. The quantities Cv,, B,, p, (and Cve, Be, Pe) are obtained from a tabulated EOS [25], as functions of the couple (p, T~)(or (p, 7~));the average ion charge state Z’~’is also given by the EOS as = Z*(p, 7)

where X~e= (1/9)V~Tae is the diffusion coefficient of the earlier model [28], b1 and be are the drag coefficients due to interactions with the ions and the electrons [29],respectively, and v0 is the initial velocity of the a-particles. The limiting flux is (Fiim)a = v0E~/3.

3.2. Thermonuclear reactions

3.3. Laser—plasma interaction

The model of DUED takes into account the D—D and D—T reactions, the burn-up of deuterium and tritium, and the transport and the energy deposition of the 3.5 MeV a-particles produced by the D—T reactions. Charged fusion products of the D—D reactions instead deposit their energy locally and instantaneously, while neutron interaction with the target is neglected. a-particle transport can be handled by rather accurate multigroup models [20,26,27]. However, for simplicity, we use a much simpler single-group model for the diffusion of the energy density Ea of the non-thermal a-particles. This is a heuristic extension [27] of a previously developed model [28], which was found to be appropriate for plasma

In moderate-intensity, short-wavelength laser fusion experiments most of the laser light is absorbed collisionally, i.e. by inverse bremsstrahlung, in the underdense plasma corona. Typical time and space scales of these plasmas allow a ray treatment of laser-plasma interaction. The ray equations are, in general, [30] x = aca/aK, (20a)

a’a/ax, (20b) where w and K are the frequency and the local wave-vector of the wave field. The ray equations for the laser light propagating in an unmagnetized plasma are found by substituting in eqs. (20a, b) E=

—

114

S. A tzeni

/

Symmetry and ,itahility of laser fusion targets

the dispersion relation K2_~_

-i

~P

c2

+

c2

—

4.1. The splitting method

(i ‘ + i— ~a 1

=

0,

where w.~,is the electron plasma frequency, and v is the electron—ion collision frequency. Setting the term in brackets equal to unity (indeed / ~ 1), we finally get d2x d(ct)2

—

i(P)

(21)

Pc

=

where ~ is the critical density at which the laser frequency equals the plasma frequency (notice that p/p~=(~a~/~)2). As we shall show in the next section, we need not solve the equations for the ray amplitudes to compute power deposition on the Lagrangian mesh. We only solve, along the rays, the power attenuation equation dW/dx = WK,B, (22) where K,B = 2 Im(k(~c)) is the inverse bremsstrahlung absorption coefficient. (The drawback of this model is that it does not yield the field intensities, and therefore ponderomotive effects cannot be taken into account.) —

4. Numerical solution of the 3-T equations The numerical solution of the system of the 3-T equations (eqs. (12)—(14)) can be performed in various ways. We consider two of them, offering alternative advantages. The first one makes use of splitting of the right-hand side (r.h.s.) operators of eqs. (12)—(14) in tems of the different physical processes. The second consists in the simultaneous, fully implicit solution of the three equations over the whole mesh (i.e. of 3P equations, where P denotes the total number of zones). The latter approach is somewhat more complex and requires about four times as much computer memory as the former, but is more accurate and ensures more rapid convergence when both some conductivities and some exchange rates are simultaneously high, while the specific heats are of different orders of magnitude.

4.1.1. TheGeneralities splitting method used here essentially consists in computing separately the temperature increments produced by different groups of physical processes. For instance, inspection of eqs. (12)— (14), (15) and (18) shows that the terms at the r.h.s., of eqs. (12)—(14) can be divided into four different groups, namely, i) ii) iii) iv)

hydrodynamical terms; outer energy sources; diffusive terms; energy exchange terms proportional to (7~— Tv).

We then apply the splitting method by considering in a first splitting (splitting I) the terms of groups i) and ii), then in splitting lithe diffusive terms, and, eventually, in splitting III, the energy exchange terms. Explicit difference schemes are used in splitting I, while implicit schemes are used in splittings II and III. The solution procedure is iterated until some convergence criterion is fulfilled. This allows time centering of the various coefficients which depend on the temperatures, in an attempt to improve both the accuracy and the energy conservation of the numerical solution. Notice that in any of the three splittings, time centering of a coefficient f~ ‘~“~(T)is done by using the values at time level n and at time level n + 1, i.e. at the end of splitting III of the previous iteration.

4.1.2.

Splitting I

In splitting I, which is formally analogous to the solution of eq. (9), we obtain, by straightforward, explicit computations, some intermediate temperatures, which we call T’~. Splitting II In this splitting we solve the diffusive part of each of the three eqs. (12)—(14) by a fully implicit scheme. Indeed, use of an explicit scheme is ruled out, for most problems, by the stability condition L~t~ pC~(/.~L)2/2X [15], which is usually considerably more stringent that the Courant condition 4.1.3.

S. Atzeni

/

Symmetry and stability of laser fusion targets

seen earlier. We then write 2

~

~~/

115

separately. For instance, we write

(T~,—T~ 1)= Evx~/2vT*1

k,I

(23)

(XP~)..÷,/2~[(Fhm)p(~)lJ+i/2]

~

where ~R1,f÷,/2 is the distance between the centres of the zones (i + 1/2,

with k = 1,2,..., I—i, and 1=1,2,..., J— 1. In eq. (23), for the sake of simplicity, we have replaced the zone indices (i + 1/2, j + 1/2) with (k, 1) (this means that zone (k, 1) has its lower left corner in the grid-point (i, J)). The diffusion operator can be differenced by using the 9-point scheme developed by Kershaw [31]. This scheme assures consistency of the finite-difference representation and local energy conservation. Furthermore, it leads to a linear system, with a symmetric matrix, which is particularly suitable for the ICCG solution [32]. According to such a scheme, the term in square brackets in eq. (23) can be written as

j + 1/2) and (i

—

1/2, j+ 1/2). To solve the system of temperature equations, we order the unknown temperatures Tk*, by columns in a vector T*, with elements T1,*,with p = 1, 2,..., P = (I 1)(J 1), 50 that p = (I 1)(1 1) + k. Using eq. (24), we can now write the set of eqs. (23) in matrix form, as —

—

—

—

Ai~* E,

(25)

=

.

-

-

with the b elements of A and B being given respective y, y pVC.,, Apq

=

~

Vpt~pq

—

Apq

and

[v

Xn

/2VT*] k,1 =

~

k’!’

~

-

(24)

where k’ = k — 1, k, k + 1, 1’ = 1— 1, 1, 1 + 1, and Vk/ is the volume of zone (k, 1). The expressions of the elements of the symmetric matrix A can be found in the original paper by Kershaw (we caution the reader about the different notation used here and in ref. [31]). Flux limitation, not dealt with in Kershaw’s treatment, can be introduced in the following approximate way. We exploit the fact that Kershaw’s scheme uses a sort of average value of the thermal conductivity for each of the two Lagrangian directions (more correctly, one should actually say suitably averaged values on zone sides), which, in the notation appear in directions, the terms ~2 and 2 for theofkref. and[31], I Lagrangian respecAtively. Thus an effective flux-limited conductivity can be obtained by multiplying the appropriate Kershaw’s term by the correcting factors (ap ) and (a~)~ for the i and j Lagrangian directions, respectively. These, in turn, are obtained by modifying eqs. (16) and (17), the i and the j components of the temperature gradient being taken

PpCv

B~=—~—~J’;7~. Eq. (25) can be solved by an ICCG method [32], which can be partially vectorized, by exploiting the particular sparsity pattern of the matrix A. Indeed, A is a 9-diagonal symmetric matrix and therefore only 5 diagonals need to be stored (e.g. those with q =p, q =p 1, q =p (I 1) + 1, q =p — (I 1) and q =p — (I 1)— 1). The incomplete Cholesky (IC) decomposition and the successive Conjugate Gradient (CG) iterations are made simpler if each non-zero diagonal of A is stored in a column of a (P X 5) matrix [33]. As shown elsewhere in detail [17], this allows vectorization on a Cray computer of all but two substeps of each CG iteration, and very simple, albeit nonvectorial, programming of thefeatures, IC decomposition. As a consequence of these our ICCG package runs very fast, requiring about 2.7 p.s per unknown for the IC decomposition and 1.6 p.s per unknown for each CG iteration on a CRAYXMP/12 computer. This means that for a typical laser fusion problem with 1000 zones the solution of the 3-T equations, with 2 iterations per step and 10 CG internal iterations, takes about 120 ms per —

—

—

—

—

116

S. AIzen,

/

Symmetry and stability of laser fusion targets

step, which is a fairly reasonable time. Unfortunately, Kershaw’s scheme does not assure the positivity of the solution. In practice this means that on a highly distorted mesh some negative temperatures may appear just ahead of a thermal wave. A remedy to this can simply consist in the re-setting to zero of any negative temperature [31], or, perhaps more effectively, in identifying undershoots and overshoots of the temperatures in the proximity of the front of a thermal wave, and in eliminating them, “a posteriori”, after the solution of eq. (25). To be sure, an approach to the diffusion equations exists [34] which does not suffer from this problem, but is considerably more complicated than that considered in this subsection.

solution of this equation is a straightforward extension [17] of that of eq. (25), but the required storage grows from 14P locations [33] to 48P locations. 4.3. Solution of the a-particle energy equation In the a-particle energy equation (19), the diffusive term is differenced as discussed in subsection 4.1.3; the energy deposition term is also differenced fully implicitly. The resulting linear system, analogous to eq. (23), is again solved by the ICCG method. ~ Numerical treatment of laser—plasma interaction 5.1. Ray tracing on the Lagrangian mesh /35]

4.1.4. Splitting III

Light rays are described by eq. (21), which

In the final splitting we simultaneously solve

shows that principal rays, i.e. rays lying in the

the three temperature equations, considering only the terms coupling the three fluids. To guarantee unconditional stability, we difference them fully implicitly, and thus write, e.g. for the ions, T~~1 T*

R—Z plane (which are those we consider in our 2-D treatment), have parabolic paths in regions of constant V(p/p~). It is therefore expedient to approximate the p/pa distribution within each of the two triangular sub-zones into which a zone is divided by one of its diagonals (e.g. that joining point (i, j) to point (i + 1, j + 1)), so as to pro-

—

1

k.I

T’~’— T.”4-’ =

~

e

duce a constant V(p/p~)inside the considered sub-zone. The V(p/p~)can be computed from the

n+i/2

T~ 1

grid-point values of ~

k.I

For each zone we then have to solve a linear system of three equations, from which we obtain the updated temperatures 7’,n+i 7~’~and 7’~4-’. 4.2. Simultaneous solution of the 3-T equations In this case the three temperature equations are solved simultaneously, with the diffusive terms differenced as in the previous splitting II, and the energy exchange terms also differenced fully implicitly. We then have to solve the matrix equation AT = B, where the unknown vector T now has 3P elements, ordered in such a way that i~,= (T.) T~+~ = (~)~‘and ~s±2P = (~)~), for s = 1, 2 P, and A is a (3P x 3P) matrix, with 11 non-vanishing diagonals (the 9 considered earlier plus two diagonals with q = p R P). By analogy with the case discussed earlier we only need to store 6 diagonals in a (3P X 6) matrix. The ICCG ,

which are obtained, in

turn, by interpolation over the zone-centered values of P/PC~computed by the code [35]. Tracing a ray therefore consists in first finding its entrance in the simulation domain and then following the sequence of exits from (and entrances into) the various sub-zones it passes through. For each sub-zone which is entered by the ray, the code computes the parameters of the parabolic ray path, and from these the coordinates of the points of exits from the sub-zone and the indices of the next sub-zone crossed by the ray. Unfortunately, this conceptually simple process requires a complex programming logic and does not allow vectorization. 5.2.

Laser

energy absorption

Integration of the equations for the ray amplitudes, or even the simpler consideration of finite

S. A tzeni

/

Symmetry and stability of laser fusion targets

117

6. Applications to model problems The main goal of 2-D simulations for ICF

BEAMLET RAY

both for present experiments and for the prospec-

_______________

---~

BEAM

—

— ——

LASER

Fig. 3. Schematization of the laser beam in the ray-tracing package of the DUED code. Notice that rays are not centered in the beamlets, but are placed stochastically. At each step energy is delivered only to zones crossed by at least one ray (dashed in figure).

transverse-size rays [12], would be extremely cornplicated on a Lagrangian mesh, but can fortunately be avoided by the use of the following method. The laser beam is divided into M (10 to 100) beamlets, each defining an “energy flux tube”, and described, at a given time level, by a single ray carrying all the power of the beamlet. The ray is then traced through the plasma as if the whole path occurred in a single time step ~t. At each time step the position of the ray in the beamlet before entrance into the simulation domain is determined stochastically (see fig. 3) in order to reproduce, in the limit of a large number of time steps, the correct spatial energy distribution of the laser beam. In this way each energy flux tube is filled stochastically and the energy deposition into the mesh can be computed simply by integrating the energy attenuation equation (eq. (22)) along the path of each of the M rays. At a given time level, and for each beamlet, power is delivered only to the zone crossed by the rays, the finite transverse size of the beamlet being neglected. We have found that even when a rather small number of rays are used, after a reasonable number of time steps the numerical noise introduced by this method is, in most case, practically negligible. Note that the integration of eq. (22) is further simplified by approximating each parabolic arch by several straight segments, their actual number depending on the ray curvature.

consists assistingand, the in design of optimal processesininvolved particular, their targets differ-

tive space ent future reactor. timeHowever, scales, often the require complexity a twofold of the implosion approach. then used and to symmetry, Full-scale study such the simulations large-scale absorption ofaspects efficiency, a target as are hydrodynamic efficiency and the thermonuclear gain. The analysis of other aspects is better performed instead by means of simulations dealing with schematized situations, which we call “model problems”, in which a particular process or class of processes is studied in detail. This approach not only saves computing resources, but also allows a better understanding of the basic features of a single phenomenon, which can be addressed with wide parametric studies. Furthermore, model problems are most suitable for the study of fluid instabilities, usually performed by introducing some ad hoc perturbation into a certain equilibrium configuration. As examples of ICF-related model problems which can be addressed by a 2-D code, we briefly describe three applications of our DUED code. 6.1. Rayleigh— Taylor instability of a fluid layer We consider the Rayleigh—Taylor instability (RTI) of fluids described by the simple model of of section 2, and by an ideal-gas EOS. This can be taken as a preliminary step towards the study of the RTI in the complicated environment of the ablation front of a laser target. We have shown elsewhere [16] that the application of DUED to the single-mode RTI of a thin fluid layer produced remarkable results, showing the theoretically predicted [36] spike—antispike behaviour, previously revealed only by ad hoc codes and by a Lagrangian code employing a very fine mesh [21]. Instead we used a rough 20 x 12 zone mesh (with 20 zones per half-wavelength in the horizontal direction); the problem was run on a CRAYXMP/12 computer in less than 20 s. Here we present a new case in which the initial

118

5. A tzeni / Symmetry and stahiliti’ of laser fusion targets

Z

a ________________________________________

—

0

bubbles with a large radius of curvature. At the instant represented in fig. 4c, the mesh is strongly deformed and the run halts. Notice that we have again used a coarse mesh (200 X 12). —I—

Z (cm)

I

I

I

a

05 I

b

0

I

I

I

I I

b

C

________

0.5

0

0

5

-

R

‘

10

Fig. 4. Rayleigh—Taylor instability of a thin layer of perfect fluid with y = 5/3, triggered by a two-wavelength velocity 2, perturbation (See main text); the gravity is g = — 10i4 cm/s and the equilibrium pressure (at z = 0) is P 2. (a) t = 30 ns: (b) t = 50 ns; ~c) t 65.4 0ns 1013 (in the dyn/cm figure R = R/(314 am) and Z= Z/(314 ~sm).

_____________________________________

0.5

hydrostatic equilibrium of a fluid layer (with adiabatic index y = 5/3 and thickness ~z = 100 p.m) is perturbed by a two-wavelength perturbation u, 0 = u0

~

e~”~sin kmr,

_____

_______________________________

m=1,2

u~0= u0

~

e_~~z cos kmr,

0

m—1,2

with k~ 100 cm~ k2 =of90thecm_i. Fig. 4 shows, as=expected, theand growth perturbation with the shortest wavelength, modulated by a wave with z~k= k 5 k2. The non-linear evolution shows the presence of thin, falling spikes, and of rising —

0.5 R (cni)

3) by an azimuthally Fig. 5. Compression of a perfect-gas cylinder (with initial radius R0 =1 pressure cm and density p0 text; =1 g/cm varying outer (see main here P 2). In the figure constant-j Lagrangian curves are shown at selected 0 = io’~dyn/cm times; (a)

t=

150 ns; (b) t = 200 ns; (c) t = 217 ns. At the mesh spacing was uniform.

t

=

0

S. Atzeni

/

Symmetry and stability of laser fusion targets

119

6.2. Implosion of a cylindrical shock wave driven by a non-uniform pressure [37] This problem is potentially relevant to the behaviour of non-uniformly irradiated thick-shell targets, as well as of non-uniformily accelerated central ignitors of complex targets [3,6]. We replace the actual spherical target with a cylindrical target, which is studied in Cartesian geometry, because this configuration allows easier treatment of periodic boundary conditions. We then study the evolution of an initially uniform cylinder, subjected to a constant outer pressure, with an azimuthal profile P=P(~)=P0(1+8cos4n~s~). The simplicity of this problem allowed a wide parametric study which has revealed the existence of different types of behaviour of the imploding shock wave, depending on the values of 6, nx, and of a parameter characterizing the EOS of the material (e.g. y in the case of the ideal-gas EOS). For instance, if the point (8, n~,y) is in a certain region of the parameter space, energetic jets of imploding matter are produced which destroy the cylindrical symmetry of the implosion [37], as shown in fig. 5, which refers to a case with 6 = 0.1, = 4, and y = 5/3 (of course only a sector 0 ~ ~ ii/8 was simulated, with periodic 8). Forboundary suitably conditions at i~ = 0 and i~ = Tr/ chosen sets of (6, n~,y), instead, these jets do not appear and the shock wave implodes symmetrically[37]. 6.3. Ignition of a D— T sphere with an asymmetric hot-spot The implosion of an ICF reactor-size shell should produce a quasi-isobaric plasma, with a central hot-spot, surrounded by a denser, relatively cold fuel [9,38]. The hot-spot should then self-heat and trigger the development of an outward propagating burn-wave [2,38,39]. The 1-D phenomenology of this process has been studied in detail and a closed-form ignition condition has been derived [39]. As a first approach to the 2-D study of this problem we have studied, with the complete 3-T

....

a \5 200

.

\ — .~

.

_

~

4

\~

~••

,

¶i’t

I’t

R

b

(ppm) 200 .

“N ~ 5N

N

~\4 \

it

I

I

I

200 ~

\

3”~\

____________________________ 0

200

Z (i1~m)

Fig. 6. Ignition of a D—T sphere, triggered by an ellipsoidal not-spot (see main text). (a) iso-temperature T~= 5 key; (b) iso-temperature i = ~o keV; (c) front of the disassembling (rarefaction) wave, shown by the outer iso-density p = 450 g/cm3. (1) t = 0; (2) 1 = 75 ps; (3) t = i25 ps; (4) 1=175 ps; (5) t 225 ps; (6) t = 250 p5.

version of DUED, the case in which an ellipsoidal hot-spot with axes br and b~,drives the burn of a sphere of cold fuel, with initial radius R~.For instance, in the case shown in fig. 6, initially the

120

S. A tzeni

/

Symmetry and stability of laser fusion targets

hot-spot has b~= 62.5 p.m, b~= 125 p.m, T= 6 keV and p = 30 g cm3, and is immersed in a sphere of D—T with R~= 154 p.m, T—~350 eV and p = 500 g cm3. The total mass of D-T is 2.8 mg, the mass of the hot-spot is 0.125 mg. Fig. 6a and b show the evolution of the burn wave fig 6c that of the inward moving disassembling wave. Studies are now in progress to find out how the asymmetric configuration modifies the 1 D ignition condition [39] and the burn-up calculations [2—4].

,

06

R

B

--

(mm)

-

-

0

—~

-~.~‘-

-

ii —---_

~

—

0.6

7. Laser—target studies I

In this section we illustrate a few applications of DUED to the study of laser-irradiated foils and very thin shell. From the computational point of view, the study of these targets poses formidable problems because of the great difference between the axial and the transverse dimensions of the initial target configuration In both cases the use of purely Eulerian codes is ruled out the use of “sliding-Eulerian” codes also appears troublesome.

0.6

I

I

b

= --

.—.. --

~ 0 -

— j:_ ~

~.

0.6

7.1. Laser accelerated thin foils As an example of a thin-foil acceleration study we show in fig. 7 the application of the 2-T version of DUED to the case of an experiment carried out at Frascati [40] in which a very thin plastic target (1 p.m thick) is irradiated by a Nd laser delivering an energy of 4 J in a 6 ns pulse with a peak intensity of 2. The focal W cm spot is about 480 p.m in diameter and the radial profile of the beam is assumed to be trapezoidal. The figure shows that the absorption of the laser light produces a partially transparent plasma corona which expands towards the laser, while the non-ablated part of the foil is accelerated in the opposite direction. Agreement with the experiment is remarkable; even the bumps appearing at the edge of the accelerated part of the foil (see fig. 7b, around R = 200 p.m and Z = —300 Ism) are very similar to those observed in the experiment. One feature, however, that cannot be reproduced is the state of the accelerated foil, which, in fact, is heated by the X-rays, which are not taken into

I

0

1 Z

(mm) Fig. 7. Thin foil target accelerated ablatively by a low-intensity . laser pulse (see main text for the parameters). The mesh and a set of representative rays are shown at t = 2 ns (a) and t = 5 ns (h). Initially the foil was at Z = 0, and the mesh was uniform in both

Z and R.

10i2

account in the model. It should be observed that despite the relatively simple physical phenomenology, this problem, with 600 zones, was run in more than 5000 steps because of the severe limitation posed by the Courant condition on the very thin zones of the accelerated part of the foil. 7.2.

Laser-imploded thin shells

The study of the implosion of a very thin shell is more complex than the previous case, because the imploding matter is now substantially cornpressed. Furthermore, most irradiation schemes

S. Aizeni

/ Symmetry

and stability of laser fusion targets

R (tim)

121 -

.‘l.-...

>1

--

__~~

——

200~

0

200

200

Z (

R (~sm’)

200

1sm)

b

a

f-~.

~

~

..-

~T:~ /1

T~7

/ I --

1

I

4---~-~

‘/: ;

.~

____ 0

200

C

d

200

Z (~sm)

Fig. 8. Implosion of a very thin shell irradiated by a two-beam laser (with the parameters given in the main text). Constant-j mesh curves (solid) and light rays (dashed) are shown at selected times; (a) t = 0; (b) t = 2.25 ns; (c) t = 4.04 ns; (d) t = 4.58 ns.

employ a large number of non-axisymmetric beams, thus being intrinsically three-dimensional. Here we illustrate a case with simple, axisymmetric 2-beam irradiation. The two beams are tangential to the initial target configuration, and are focused at R = 0 and z = R0 and R0, respectively (here R0 is the initial target radius). Both beams have a flat spatial power profile; accordingly, the intensity at the poles of the target is 2.2 times larger than at the equator. Each beam delivers about 110 J at wavelength A = 0.53 p.m in a 6 ns triangular pulse, rising at

V~

—

peak power in 1 ns. The target is a glass shell with R0 = 250 p.m and thickness L1R0 = 2.5 rim, and is filled with 3 atm (STP) of D—T. Figs. 8 and 9 show that the non-uniform irradiation produces a bean-shaped compressed configuration~ At t = 4.58 ns the mesh has suffered such a large distortion that the program is stopped. Fig. 10 shows the presence of some mesh pathologies on the outer surface of the imploded shell (i.e. on the ablation front) as well as in the filling gas (the latter types being visible in the lower right-hand part of the figure). The continuation of this run

122

S. At:emn

R 2

~ --•~/f=~(

Svnimetrv and stability of laser Ju.sion targets

R)lsm( 0

‘~“N~

~

~

/

~

// -200

L

-200

0

200 --

Z )~m)

~

/

____________

—~

Fig. 9. Same as fig. 8d. hut on a different scale.

0

100

~

requires the use of a re-zoning technique, discussion of which would exceed the scope of this paper. It is worthy observing, however, that the study of more uniformily irradiated thin targets, as well as of thicker targets (regardless of the irradiation uniformity) is considerably easier and less expensive. For instance, in fig. 11 we show the configuration at the instant of the compression maxi-

Fig. 11. Mesh lines showing the imploded configuration of a shell with initial aspect ratio A = iO. irradiated by a tsso-heani . the parameters given . 1/2 laser, with in the main text.

mum of a thicker plastic target (R() = 100 ~sm, z~R11= 10 p.m) irradiated by a 2-beam f/2 laser delivering about 200 J at A = 0.53 p.m in I ns. The shape of the central region is qualitatively very similar to the X-ray images taken in experiments with similar parameters. 7.3. Laser-irradiated cavities

100

1

N /

/ R

(is

/

/

m)

__

I

/

/

Experiments have been performed for a couple of years in which a hollow sphere with an “entrance hole” is irradiated on the inner surface by a laser [41] beam(see focused at the of most the entrance hole fig. 12). In centre this case of the energy emitted by the target is confined in the cavity. If the appropriate conditions are met, this radiation can drive an almost spherically symmetnc ablative radiation wave [42,43]. Accurate study of this problem would be extremely complex, requiring the non-LTE treatment of the properties

N

Z

(i~m)

100

Fig. 10. Enlarged view of the mesh of the imploded target of fig. 8d at t 4.58 ns. The arrows indicate the boundaries of the glass shell.

of the matter,However, and some aaccurate radiativeapproach transfer technique. zeroth-order could be devised as follows, at least for the socalled LTE-AHW regime [42]. We start from the observation that i) the radia-

S. A tzeni

/ Symmetry and stability of laser fusion targets

123

concepts. At the same time these results show that progress is needed in at least two directions. Indeed, accurate and, possibly, flexible and auto-

200 R ~

(pMl(

perform quantitative evaluations of the effects of asymmetries on techniques very thin are shells. Also, matic re-zoning needed in methods order to should be developed for the treatment of radiative transfer in cavity-type geometric configurations that confine most of the radiation emitted by the target. It is also worth observing that in some cases the 2-D geometry can only be assumed as a rather crude approximation of the actual 3-D geometry.

-_ ~.

200 -200

200 Z (tam)

Fig. 12. Laser-irradiated cavity. Radiative transfer is not ineluded in the simulation. The picture therefore only serves to show the geometry of a typical experiment.

References [11 J. Nuckolls, L. Woods, A. Thiessen and G. Zimmerman, Nature 239 (1972) 139.

tion field should, in this case, be quasi-isotropic; ii) the radiation energy in the cavity is only a small fraction of the internal energy of the target. One could then simulate these features by imposing a no-flux boundary condition at the inner surface of the hollow sphere and by setting the radiative conductivity on the same inner surface to a very large value. Such a model for the radiation cannot account for the isotropization of the radiation field in the cavity and, in a sense, could be considered as 1-D. However, it should be useful to estimate the relative importance of the directly laser-driven ablation wave and of the radiation heat wave in a more realistic geometry than that considered by 1-D codes. We are presently studying this problem; the results will be reported elsewhere,

[21K.A.

Brueckner and S. Jorna, Rev. Mod. Phys. 46 (1974) 325. [3] J.J. Duderstadt and G.A. Moses, Inertial Confinement Fusion (Wiley, New York, 1982).

[4] T.H. Johnson, Proc. IEEE 72 (1984) 548. [5] C. Yamanaka et al., in: Plasma Physics and Controlled Fusion Research 1984, vol. 3 (IAEA, Vienna, 1985) p. 3. [61 S. Bodner, J. Fusion Energy 1 (1981) 3. [71R.G. Evans et al., in ref. [51,p. 25. [81R.L. McCrory et al., in ref. [51,p. 37. [9] Yu.V. Afanasev et al., JETP Lett. 21(1975) 68. [10] S. Atzeni, NucI. Fusion 24 (1984) 1220. [11] C. Yamanaka et al., Phys. Rev.Rev. Lett.Lett. 56 (1986) 1575. MC. Richardson et al., Phys. 56 (1986) 2048. [121 R.5. Craxton and R.L. McCrory, J. Appl. Phys. 56 (1984) 108. [131 GB. Zimmermann and W.L. Kruer, Comments Plasma Phys. Cont. Fusion 2 (1975) 51 (see also the 1973-1983 issues of the Lawrence Livermore Laboratory Laser ProAnnual/Semiannual [14] gram J.P. Christiansen and N.K. Report). Winsor, J. Comput. Phys. 35 (1980) 291; Comput. Phys. Commun. 17 (1979) 397. G.I. Pert, J. Comput. Phys. 43 (1981) 111.

8. Conclusions In this paper we have discussed the main computational aspects of 2-D Lagrangian codes for laser fusion. The potential of these codes has been illustrated by means of several examples of application of the DUED code. It is apparent that such 2-D codes are suitable for the investigation of many crucial issues relevant to the symmetry and stability of the most promising laser fusion target

[151 D.M. Richtmeyer and K.W. Morton, Difference Methods for Initial Value Problems, 2nd ed. (Interscience, New York, 1967). [161S. Atzeni, Comments Plasma Phys. Cont. Fus. 10 (1986) in press. [17] S. Atzeni, to be published; also Report ENEA/RT/ FUS8S/9 (Frascati, Italy, 1985). [18] W.D. Schulz, in: Methods in Computational Physics, vol. 3, eds. B. Alder et al. (Academic Press, New York, 1964) ~ 1. [19] M.J. Fritts, in: Finite Difference Techniques for Vectorized Fluid Dynamics Calculations, ed. DL. Book (Springer,

124

S. A t:eni

/

Symmetry amid stability of laser fusion targets

New York, 1981) p. 98. ID. Sofronov Ct al., in: Numerical Methods in Fluid Dynamics, edo. N.N. lanenko et al. (Mir, Moscow, 1984) p. 82. [201 G.S. Fraley, E.J. Linnebur. R.J, Mason and R.L. Morse, Phys. Fluids 17 (i974) 474. [21] C.P. Vernon ci al., Phys. Fluids 25 (1982) 1653. [22] K. Mima ci al., in ref. [5], p. 113. P.C. Thomson and P.D. Roberts, Las. Part. Beams 2 (1984) 13. E.G. Gamalii ci al., Soy. Phys. JETP 52 (i980) 230. [23] C.D. Levermore and G.C. Pomraning, Astrophys. J. 248 (i981) 321. [24] L. Spitzer, Physics of Fully Ionized Gases (Wiley, New York, 1962). [25] S. Atzeni, A. Caruso and V.A. Pais. Laser Part. Beams, 4 (1986) 393. [26] E.G. Corman et al., NucI. Fusion 15 (1975) 377. [27] S. Atzeni, Report ENEA 84.1/cc (Frascati, Italy. 1984). [28] S. Atzeni and A. Caruso. Nuovo Cimento 64B (1981) 383. [29] DV. Sivukhin, in: Reviews of Plasma Physics, vol. 4, ed. A.M. Leontovich (Consultants Bureau. New York. 1966) p. ioi. [30] M. Born and E. Wolf, Principle of Optics (Pergamon. Oxford. 1964) chap. 3. lB. Bernstein. Phys. Fluids 18 (1975) 320.

[3i] D.S. Kershaw, J. Comput. Phvs. 39 (1981) 375. [32] D.S. Kershaw, J. Comput. Phys. 26 (i978) 43. 1331 D.V. Anderson and At. Shestakov, Comput. Phys. C~ornmun. 30 (1983) 37. [34] (il. Pert. J. Comput. Phys. 42 (198i) 20. [35] D.S. Bailey. Lawrence Livermore Laboratory Report UCRL 75881 (Livermore, i974). Sec also A. Friedman, in: Laser Fusion Anual Report 1983. Rep. UCRL-50021-83 (Livermore, 1984). [36] E. Ott, Phys. Rev. Lett. 29 (i972) 1429. [37] S. Atzeni. A. Caruso and V.A. Pais. Phvs. Lett. (i986) submitted. [38] S.Yu. Guskov, ON. Krokhin and V.0. Rozanoy, Nuel. Fusion 16 (1976) 957. 39] S. Atzeni and A. Caruso, Phys. Lett. ISA (i98i) 345: Nuovo Cimento 80B (1984) 71: NucI. Fusion 23 (1983) 1092. [40] A. Caruso. in: Advances in Inertial Confinement Fusion Research (ILE, Osaka, 1984) p. 121. [41] lB. Foldes Ct al.. Europhys. Lett. 2 (1986) 221. [42] R. Pakula and R. Sigel, Z. Naturforsch. 4iA (1986) 463. [43] J. Meyer-ter-Vehn, R. Pakula. R. Sigel and K. Unterseer. Phys. Lett. A 104 (1984) 410.