Description

The development of radiation hydrodynamical methods that are able to follow gas dynamics and radiative transfer (RT) self-consistently is key to the solution of many problems in numerical astrophysics. Such fluid flows are highly complex, rarely

Information

Category:
## Spelling

Publish on:

Views: 0 | Pages: 36

Extension: PDF | Download: 0

Share

Transcript

a r X i v : 0 9 0 5 . 2 9 2 0 v 1 [ a s t r o - p h . C O ] 1 8 M a y 2 0 0 9
Mon. Not. R. Astron. Soc.
000
, 1 (2008) Printed 18 May 2009 (MN L
A
TEX style ﬁle v2.2)
Cosmological Radiative Transfer Comparison Project II: TheRadiation-Hydrodynamic Tests
Ilian T. Iliev
1
,
2
,
3
⋆
, Daniel Whalen
4
, Garrelt Mellema
5
, Kyungjin Ahn
6
,
7
, Sunghye Baek
8
,Nickolay Y. Gnedin
9
, Andrey V. Kravtsov
10
, Michael Norman
11
, Milan Raicevic
12
,Daniel R. Reynolds
13
, Daisuke Sato
14
, Paul R. Shapiro
6
, Benoit Semelin
7
, Joseph Smidt
15
,Hajime Susa
16
, Tom Theuns
12
,
17
, Masayuki Umemura
14
1
Astronomy Centre, Department of Physics & Astronomy, Pevensey II Building, University of Sussex, Falmer, Brighton BN1 9QH, United Kingdom
2
Universit¨ at Z¨ urich, Institut f¨ ur Theoretische Physik, Winterthurerstrasse 190, CH-8057 Z¨ urich, Switzerland
3
Canadian Institute for Theoretical Astrophysics, University of Toronto, 60 St. George Street, Toronto, ON M5S 3H8, Canada
4
T-2 Nuclear and Particle Physics, Astrophysics, and Cosmology, Los Alamos National Laboratory, Los Alamos, NM 87545, U.S.A.
5
Dept. of Astronomy and Oskar Klein Centre, AlbaNova, Stockholm University, SE-10691 Stockholm, Sweden
6
Department of Earth Science Education, Chosun University, Gwangju 501-759, Korea
7
Department of Astronomy, University of Texas, Austin, TX 78712-1083, U.S.A.
8
LERMA, Observatoire de Paris, 77 av Denfert Rochereau, 75014 Paris, France
9
Fermilab, MS209, P.O. 500, Batavia, IL 60510, U.S.A.
10
Dept. of Astronomy and Astrophysics, Center for Cosmological Physics, The University of Chicago, Chicago, IL 60637, U.S.A.
11
Center for Astrophysics and Space Sciences, University of California, San Diego, 9500 Gilman Drive, La Jolla, CA 92093-0424, U.S.A.
12
Institute for Computational Cosmology, Durham University, Durham, United Kingdom
13
Department of Mathematics, 208 Clements Hall, Southern Methodist University, Dallas, TX 75275, USA
14
Center for Computational Sciences, University of Tsukuba, Tsukuba, Ibaraki 305-8577, Japan
15
Department of Physics and Astronomy, 4129 Frederick Reines Hall, UC Irvine, Irvine, CA 84602, U.S.A.
16
Department of Physics, Konan University, Kobe, Japan
17
Department of Physics, University of Antwerp, Campus Groenenborger, Groenenborgerlaan B-171, B2020 Antwerp, Belgium
18 May 2009
ABSTRACT
The development of radiation hydrodynamical methods that are able to follow gas dynamicsand radiative transfer self-consistently is key to the solution of many problems in numeri-cal astrophysics. Such ﬂuid ﬂows are highly complex, rarely allowing even for approximateanalytical solutions against which numerical codes can be tested. An alternative validationprocedure is to compare different methods against each other on common problems, in or-der to assess the robustness of the results and establish a range of validity for the methods.Previously, we presented such a comparison for a set of pure radiative transfer tests (i.e. forﬁxed, non-evolving density ﬁelds). This is the second paper of the Cosmological RadiativeTransfer (RT) Comparison Project, in which we compare 9 independent RT codes directlycoupled to gasdynamics on 3 relatively simple astrophysical hydrodynamics problems: (5)the expansion of an H II region in a uniform medium; (6) an ionization front (I-front) in a
1
/r
2
density proﬁle with a ﬂat core, and (7), the photoevaporationof a uniform dense clump.Results show a broad agreement between the differentmethods and no big failures, indicatingthat the participating codes have reached a certain level of maturity and reliability. However,many details still do differ, and virtually every code has showed some shortcomings and hasdisagreed, in one respect or another, with the majority of the results. This underscores the factthat no method is universal and all require careful testing of the particular features which aremost relevant to the speciﬁc problem at hand.
Key words:
H II regions—galaxies:high-redshift—intergalactic medium—cosmology:theory—radiative transfer— methods: numerical
⋆
e-mail: I.T.Iliev@sussex.ac.uk
1 INTRODUCTION
The transfer of ionizing radiation through optically-thick media isakey process in many astrophysical phenomena. Some examples in-
2
I. T. Iliev, et al.
cludecosmological reionization (
e.g.
Gnedin 2000; Nakamoto et al.2001; Razoumov et al. 2002; Sokasian et al. 2003; Ciardi et al.2003; Iliev et al. 2006b; Kohler et al. 2007), star formation (
e.g.
Hosokawa & Inutsuka 2005; Iliev et al. 2006a; Razoumov et al.2006; Susa & Umemura 2006; Ahn & Shapiro 2007; Whalen &Norman 2008b), radiative feedback in molecular clouds (Mellemaet al. 2006a; Mac Low et al. 2007; Dale et al. 2007a,b; Krumholzet al. 2007; Gritschneder et al. 2009), and planetary nebulae (
e.g.
Mellema et al. 1998; Lim & Mellema 2003). In some of theseproblems, fast, R-type I-fronts predominate. Those fronts propa-gate faster than the hydrodynamic response of the gas, so gas mo-tions do not affect the I-front evolution. In these cases the radiativetransfer could be done on a ﬁxed density ﬁeld (or a succession of such ﬁelds), and dynamic coupling to the gas is generally not re-quired. However, the majority of astrophysical and cosmologicalapplications involve slow, D-type I-fronts (or a combination of R-type and D-type, as we describe in detail in section 3.1), so theradiative transfer and gasdynamics should be directly coupled andevolved simultaneously. Until recently, self-consistent radiationhy-drodynamical codes for radiative transport have been rare, but thisunsatisfactory situation is now rapidly changing due to the devel-opment of a number of such codes using a variety of numericalapproaches.A number of radiative transfer methods have been developedin recent years, both stand-alone and coupled to hydrodynamics.High computational costs necessitate the usage of various approx-imations. Thus, it is of prime importance to validate the numericalmethods developed and to evaluate their reliability and accuracy.Tests with either exact or good approximate analytical solutionsshould always be the ﬁrst choice for code testing. Extensive testsuites of radiation hydrodynamical I-front transport in a variety of stratiﬁedmediawithgood approximate analytical solutions doexist(Francoetal. 1990; Whalen&Norman2006) and arestringent testsof coupling schemes between radiation, gas, and chemistry. How-ever, an alternative and complementary approach is to compare avariety of methods on a set of well-deﬁned problems in astrophys-ical settings. This is the approach we have taken in this project.Our aim is to determine the type of problems the codes are(un)able to solve, to understand the srcin of any differences in-evitably found in the results, to stimulate improvements and fur-ther developments of the existing codes and, ﬁnally, to serve as abenchmark for testing future algorithms. All test descriptions, pa-rameters, and results can be found at the project website:
http
:
//www.cita.utoronto.ca/
∼
iliev/rtwiki/doku.php
.The ﬁrst paper of this comparison project discussed the re-sults from ﬁxed density ﬁeld tests (Iliev & et al. 2006, hereafterPaper I), i.e. without any gas evolution. We found that all partic-ipating codes are able to track I-fronts quite well, within
∼
10%of each other. Some important differences also emerged, especiallyin the derived temperatures and spectral hardening. We found thatsome of these differences were due to variations in microphysics(chemical reaction rates, heating/cooling rates and photoionizationcross-sections), whileothers were due to the method itself,
e.g.
howthe energy equation is solved, how many frequency bins are usedfor the spectral evolution, etc. We concluded that the tested radia-tive transfer methods are producing reliable results overall, but thatnot all methods are equally appropriate for any given problem, es-pecially in cases when obtaining precise temperatures and spectralfeatures is important.We now extend our previous work by considering aset of radi-ation hydrodynamical tests. In the spirit of Paper I, we have chosena set of test problems which are relatively simple, so as to be most
Figure 1.
Legend for the line plots.
inclusive given the current limitations of the available codes (
e.g.
1-D or 2-D vs. 3-D codes). At the same time, our tests considerproblems of astrophysical importance, and cover a wide variety of situations that test the attributes of each method, including its ra-diative and hydrodynamic components and their coupling.The efﬁciency, optimization and performance of the codesare very important, especially for the most complex andcomputationally-intensive problems. However, there are a numberof complications, which we discussed in Paper I, preventing usfrom doing such testing in a meaningful way at present. We there-fore leave it for future work.All test results for this study had to be supplied on a regularCartesian grid of
128
3
computational cells. This relatively modestresolution was chosen in the interests of inclusivity, so that evencodes which are not yet fully optimized in terms of either compu-tations or memory can participate in the comparison. We note thatproduction runs at present are typically run at
256
3
or better resolu-tion. Codes which utilize Adaptive Mesh Reﬁnement (AMR) gridsor particles have been requested to run the problem at the resolu-tion which approximates as closely as possible the ﬁxed-grid onefor fair comparison. Their results have then been interpolated ontoa regular grid for submission.
2 THE CODES
In this section we brieﬂy describe the nine radiative transfer codesparticipating in this stage of the comparison project, with refer-ences to more detailed method papers if available. Details of thecodes and their basic features and methods are summarized in Ta-ble 1. Figure 1 provides a legend allowing the reader to identifywhich line corresponds to which code in the ﬁgures throughout thepaper. The images we present are identiﬁed in the correspondingﬁgure caption.
2.1 Capreole
+
C
2
-Ray and TVD
+
C
2
-Ray (G. Mellema, I.Iliev, P. Shapiro, M. Alvarez)
C
2
-Ray (Mellema et al. 2006b) is a grid-based short characteristics(
e.g.
Raga et al. 1999) ray-tracing code which is photon-conservingand causally traces the rays away from the ionizing sources up to
Cosmological Radiative Transfer Comparison II
3
Table 1.
Participating codes and their current features.Code Grid Parallelization hydro method rad. transfer methodCapreole+
C
2
-Ray ﬁxed shared/distributed Eulerian, Riemann solver short-characteristics ray-tracingTVD+
C
2
-Ray ﬁxed shared/distributed Eulerian, TVD solver short-characteristics ray-tracingHART AMR shared/distributed Eulerian, Riemann solver Eddington tensor momentRSPH particle-based distributed SPH long-characteristics ray-tracingZEUS-MP ﬁxed distributed Eulerian 3-D ray-tracingRH1D sph. Lagrangian no Lagrangian 1-D ray-tracingCoral AMR no Eulerian, ﬂux-vector splitting short-characteristics ray tracingLICORICE AMR shared SPH Monte-Carlo ray-tracingFlash-HC AMR distributed Eulerian, PPM Hybrid characteristics ray-tracingEnzo-RT ﬁxed distributed Eulerian, PPM Flux-limited diffusion
each cell.Explicitphoton-conservation isassured bytaking aﬁnite-volume approach when calculating the photoionization rates, andby using time-averaged optical depths. The latter property allowsfor integration time steps that are much larger than the ionizationtime scale, which results in a considerable speed-up of the calcula-tion and facilitates the coupling of the code to gasdynamic evolu-tion. The code and the various tests performed during its develop-ment are described in detail in Mellema et al. (2006b).The frequency dependence of the photoionization rates andphotoionization heating rates are dealt with by using frequency-integrated rates, stored as functions of the optical depth at the ion-ization threshold. In its current version the code includes only hy-drogen and no helium, although it could be added in a relativelystraightforward way.The transfer calculation is done using short characteristics,where the optical depth is calculated by interpolating values of grid cells lying along the line-of-sight towards the source. Becauseof the causal nature of the ray-tracing, the calculation cannot eas-ily be parallelized through domain decomposition. However, us-ing OpenMP and MPI the code is efﬁciently parallelized over thesources and grid octants (Iliev et al. 2008b). The code has beenapplied for large-scale simulations of cosmic reionization and itsobservability (Iliev et al. 2006b; Mellema et al. 2006c; Iliev et al.2007a,b; Holder et al. 2007; Dor´e et al. 2007; Iliev et al. 2008a) ongrid sizes up to
406
3
and up to
∼
10
6
ionizing sources, running onup to 10,240 computing cores.There are 1D, 2D and 3D versions of the code that are avail-able. It was developed to be directly coupled with hydrodynamicscalculations. The large time steps allowed for the radiative trans-fer enable the use of the hydrodynamic time step for evolving thecombined system. The
C
2
-Ray radiative transfer and nonequilib-rium chemistry code has been coupled to several different gasdy-namics codes, utilizing both ﬁxed and adaptive grids. The tests inthis project were mostly performed with the version coupled to thehydrodynamics code Capreole developed by Garrelt Mellema andbased on Roe’s approximate Riemann solver. The ﬁrst gasdynamicapplication of our code is presented in Mellema et al. (2006a). Ad-ditionally, one of the tests has also been run with C
2
-Ray coupledto a different hydro solver, namely the TVD method of Trac & Pen(2004) (see Test 6 below).
2.2 Hydrodynamic Adaptive Reﬁnement Tree (HART) (N.Gnedin, A. Kravtsov)
The Hydrodynamic Adaptive Reﬁnement Tree (HART) code is animplementation of the AMR technique and uses a combination of multi-level particle-mesh and shock-capturing Eulerian methodsfor simulating the evolution of the dark matter particles and gas,respectively. High dynamic range is achieved by applying adaptivemesh reﬁnement to both gas dynamics and gravity calculations.The code performs reﬁnements locally on individual cells, andcells are organized in reﬁnement trees (Khokhlov 1998). The datastructureisdesigned both to reduce the memory overhead for main-taining atree and to fully eliminatethe neighbor search required forﬁnite-difference operations. All operations, including tree modiﬁ-cations and adaptive mesh reﬁnement, can be performed inparallel.The advantage of the tree-based AMR is its ability to control thecomputational mesh on the level of individual cells. This results ina very efﬁcient and ﬂexible (and thus highly adaptive) reﬁnementmesh which can be easily built and modiﬁed and, therefore, effec-tively match the complex geometry of cosmologically interestingregions: ﬁlaments, sheets, and clumps. Several reﬁnement criteriacan be combined with different weights allowing for a ﬂexible re-ﬁnement strategy that can be tuned to the needs of each particularsimulation. The adaptive reﬁnement in space is accompanied by atemporal reﬁnement (smaller time steps on meshes of higher reso-lutions).The ART code was initially developed by A. Kravtsov in col-laboration with A. A. Klypin and A. M. Khokhlov (Kravtsov et al.1997; Kravtsov 1999; Kravtsov et al. 2002). N. Gnedin joinedthe HART code development team in the spring of 2003 and hasadopted the OTVET algorithm for modeling 3D radiative transferfor the ART mesh structure and implemented a non-equilibriumchemical network and cooling (e.g. Gnedin et al. 2008).
2.3 RSPH (H. Susa, M. Umemura, D. Sato)
The Radiation-SPH (RSPH) scheme is speciﬁcally designed to in-vestigate the formation and evolution of ﬁrst-generation objects at
z
>
∼
10
, where the radiative feedback from various sources playsimportant roles. The code can compute the fraction of chemicalspecies e
−
, H
+
, H, H
−
, H
2
, and H
+2
by fully implicit time integra-tion. It also can deal with multiple sources of ionizing radiation, aswell as with Lyman-Werner band photons.Hydrodynamics is calculated by the smoothed particle hydro-dynamics (SPH) method. It uses the version of SPH by Umemura
4
I. T. Iliev, et al.
(1993) with the modiﬁcation according to Steinmetz & Mueller(1993), and adopts the particle resizing formalism by Thacker et al.(2000). The present version does not use the so-called entropy for-malism (Springel & Hernquist 2002). The non-equilibrium chem-istry and radiative cooling for primordial gas are calculated usingthe code developed by Susa & Kitayama (2000), where H
2
coolingand reaction rates are taken from Galli & Palla (1998a).As for the photoionization process, the on-the-spot approxi-mation is employed (Spitzer 1978), meaning that the transfer of ionizing photons directlyfrom thesource issolved, but diffuse pho-tons are not transported. Instead, it is assumed that recombinationphotons are absorbed in the same zone from which they are emit-ted. Due to the absence of the source term in this approximation,the radiation transfer equation becomes very simple. Solving thetransfer equation reduces to the easier problem of assessing the op-tical depth from the source to every SPH particle.The optical depth is integrated utilizing the neighbour listsof SPH particles. It is similar to the code described in Susa& Umemura (2004), but can now also deal with multiple pointsources. In the new scheme fewer grid points are created on thelight raythan initspredecessor. Instead, just onegrid point per SPHparticle is created in the particle’s neighborhood. The ’upstream’particle for each SPH particle on its line of sight to the source isthen found. Then the optical depth from the source to the SPH par-ticle is obtained by summing up the optical depth at the ’upstream’particle and the differential optical depth between the two particles.The code is parallelized with the MPI library. The computa-tional domain is divided by the Orthogonal Recursive Bisectionmethod. The parallelization method for radiation transfer is sim-ilar to the Multiple Wave Front method developed by Nakamotoet al. (2001) and Heinemann et al. (2006), but it is adapted to ﬁt theSPH code as described in (Susa 2006).The code computes self-gravity using a Barnes-Hut tree,which is parallelized as well. A Tree-GRAPE version of the codehas also been developed. This code has been applied to radiativefeedback in primordial star formation (Susa & Umemura 2006;Susa 2007; Hasegawa et al. 2009), as well as the regulation of star formation in forming galaxies by ultraviolet background (Susa2008).
2.4 ZEUS-MP (D. Whalen, J. Smidt, M. Norman)
ZEUS-MP solves explicit ﬁnite-difference approximations to Eu-ler’s equations of ﬂuid dynamics self-consistently with a 9-speciesprimordial gasreactionnetwork (
H
,
H
+
,
He
,
He
+
,
He
++
,
H
−
,
H
2
,
H
+2
and
e
) and ray-tracing radiative transfer, which is used to com-pute the radiative rate coefﬁcients required by the network and thegas energy equation. Our method is described in detail elsewhere(Whalen & Norman 2006, 2008b); here, we review multifrequencyupgrades to the radiative transfer and improvements to the subcy-cling scheme (Whalen & Norman 2008a).The ZEUS-MPRTmodule evaluates radiative ratecoefﬁcientsby solving the static equation of transfer in ﬂux form. To obtain thetotal rate coefﬁcient
k
for a zone we sum the
k
ν
computed for agiven binned photon emission rate over all energies by looping thesolution to the transfer equation over them. In tests spanning 40to 2000 energy bins, good convergence is found with 120 bins, 40bins spaced evenly in energy from 0.755 eV to 13.6 eV and 80 binsthat are logarithmically-spaced from 13.6 eV to 90 eV.Successive updates to the reaction network and gas energy areperformed over the minimum of the chemical time
t
chem
= 0
.
1
n
e
+ 0
.
001
n
H
˙
n
e
.
(1)and the photoheating/cooling time
t
hc
= 0
.
1
e
gas
˙
e
ht/cool
(2)until the larger of these two times has been crossed, at which pointfull hydrodynamical updates of gas densities, energies, and veloc-ities are performed. These times are global minima for the entiregrid. Chemical times are deﬁned in terms of electron ﬂow to ac-commodate all chemical processes rather than just ionizations orrecombinations. Adopting the minimum of the two times for chem-istry and gas energy updates enforces accuracy in the reaction net-work when
t
chem
becomes greater than
t
hc
(in relic H II regions,for example).ZEUS-MP is now fully parallelized for three-dimensional ap-plications. We have updated the H and He recombination and cool-ing rates responsible for some minor departures between ZEUS-MP and the other codes in Paper I in the temperature structureof H II regions, and now use the most recent data from Hummer(1994) and Hummer & Storey (1998). Our code has been validatedwith stringent tests of R-type and D-type I-fronts in a variety of stratiﬁed media (Franco et al. 1990; Whalen & Norman 2006) andapplied to both cosmological and astrophysical problems, such asthe breakout of UV radiation from primordial star-forming clouds(Whalen et al. 2004), the formation of dynamical instabilities ingalactic H II regions (Whalen & Norman 2008b), the circumstel-lar environments of gamma-ray bursts (Whalen et al. 2008b), thephotoevaporation of cosmological minihalos by nearby primordialstars (Whalen et al. 2008a), and Pop III supernovae explosions incosmological H II regions (Whalen et al. 2008c).
2.5 RH1D (K. Ahn, P. Shapiro)
RH1D is a 1D, Lagrangian, spherically-symmetric, radiation-hydrodynamics code for a two-component gas of baryons and col-lisionless dark matter coupled by gravity (Ahn & Shapiro 2007).For the baryonic component, the Euler equations and the equationof state are solved, together with multi-frequency, multi-speciesradiative transfer equations and a reaction network with nine pri-mordial species (
H
,
H
+
,
He
,
He
+
,
He
++
,
H
−
,
H
2
,
H
+2
and
e
).Dark matter dynamics, governed by the collisionless Boltzmannequations, takes a simpliﬁed form in spherical symmetry. The codesolves an effective set of Euler equations for a dark matter ﬂuid,based upon the “ﬂuid approximation” of dark matter dynamics fora spherically symmetric system with an isotropic velocity disper-sion, derived and justiﬁed elsewhere (Ahn & Shapiro 2005). Theseeffective Euler equations are identical to those for an inviscid, idealgas with a ratio of speciﬁc heats
γ
= 5
/
3
.The Euler equations are solved using the so-called “leap-frog”method, where the Lagrangian position (radius) and velocity (ra-dial velocity) are staggered in time to achieve a second-order accu-racy in time steps, both for baryonic and dark matter ﬂuid. Theusual artiﬁcial viscosity scheme is used to capture shocks. Wetypically adopt a few thousand uniformly spaced bins in radius.Non-equilibrium rate equations for the nine primordial species aresolved using the backward differencing scheme of Anninos et al.(1997). For
H
−
and
H
+2
, due to their relatively fast reaction rates,the equilibrium values may be used.Radiative transfer is performed by ray-tracing, taking accountof the optical depth to bound-free opacity of H I, He I, He II,
H
−
,and
H
2
, as well as bound-free and dissociation opacity of
H
+2
. The
Cosmological Radiative Transfer Comparison II
5
optical depth to the Lyman-Werner band photons of
H
2
, which arecapable of dissociating
H
2
, is treated using a pre-calculated self-shielding function by Draine & Bertoldi (1996), which is deter-mined by the
H
2
column density and gas temperature. Diffuse ﬂuxisnot explicitly calculated, but isaccounted for implicitlyby adopt-ing case B recombination rates. The radiative reaction rates arecalculated using a photon-conserving scheme, which enables thecode to treat optically-thick shells (
e.g.
Razoumov & Scott 1999;Abel et al. 1999a). A wide range of radiation frequency (energy),
hν
∼
[0
.
7
−
7000]eV
, is covered by a few hundred, logarithmi-callyspaced bins,together withadditive,linearlyspacedbinswhereradiative cross sections change rapidly as frequency changes. Foreach frequency and species, the corresponding radiative reactionrate is calculated, then summed over frequency to obtain the netradiative reaction rate.The radiative transfer scheme is able to treat 1) an internalpoint source, 2) an external, radially-directed source, and 3) an ex-ternal, isotropic background. The transfer for (1) and (2) is 1D, per-formed along the radial direction only. For (3), the transfer is 2D innature, and at each point the mean intensity is required to calculatethe radiative rates, which involves an angle integration. The radia-tivetransfer calculation isperformed for each pre-selected angle (
θ
,measured from the radial direction), and then the angle integral iscalculated using the Gaussian quadrature method.The code adopts a very stringent time step criterion for accu-racy. The minimum of dynamical, sound-crossing, cooling/heating,and species change time scales, which is multiplied by a coefﬁcientsmaller than unity (
∼
0
.
1
), is chosen as the time step. All the Eulerequations and rate equations are solved with this time step, whichmakes the whole calculation self-consistent. This code has beentested extensively and used to study the radiative feedback effectsby the ﬁrst stars on their nearby minihalos (Ahn & Shapiro 2007).
2.6 Coral (I. Iliev, A. Raga, G. Mellema, P. Shapiro)
CORAL is a 2-D, axisymmetric Eulerian ﬂuid dynamics AMRcode (see Mellema et al. 1998; Shapiro et al. 2004, and refer-ences therein for detailed description). It solves the Euler equationsin their conservative ﬁnite-volume form using the second-ordermethod of van Leer ﬂux-splitting, which allows for correct and pre-cise treatment of shocks. The grid reﬁnement and de-reﬁnementcriteria are based on the gradients of all code variables. When thegradient of any variable is larger than a pre-deﬁned value the cell isreﬁned, and when the criterion for reﬁnement is not met the cell isde-reﬁned.The code follows, by a semi-implicit method, the non-equilibrium chemistry of multiple species (H, He, C II-VI, N I-VI, O I-VI, Ne I-VI, and S II-VI) and the corresponding cooling(Raga et al. 1997; Mellema et al. 1998), as well as Compton cool-ing. The photoheating rate is the sum of the photoionization heat-ing rates for H I, He I and He II. For computational efﬁciency allheating and cooling rates are pre-computed and stored in tables.The microphysical processes – chemical reactions, radiative pro-cesses, transfer of radiation, heating and cooling – are implementedthough the standard approach of operator-splitting (i.e. solved ateach time-step, side-by-side with the hydrodynamics and coupledto it through the energy equation). The latest versions of the codealso include the effects of an external gravity force.Currently the code uses a black-body or power-law ionizingsource spectrum, although any other spectrum can be accommo-dated. Radiativetransfer of theionizing photons istreatedexplicitlyby taking into account the bound-free opacity of H and He in thephotoionization and photoheating rates. The photoionization andphotoheating rates of H I, He I and He II are pre-computed for thegiven spectrum and stored in tables vs. the optical depths at the ion-izing thresholds of these species, which are then used to obtain thetotal optical depths. The code correctly tracks both fast (by evolv-ing on an ionization timestep,
∆
t
∼
˙
n
H
/n
H
) and slow I-fronts.The code has been tested extensively and has been applied tomany astrophysical problems,
e.g.
photoevaporation of clumps inplanetary nebulae (Mellema et al. 1998), cosmological minihalophotoevaporation during reionization (Shapiro et al. 2004; Ilievet al. 2005), and studies of the radiative feedback from propagat-ing ionization fronts on dense clumps in damped Lyman-
α
systems(Iliev et al. 2006a).
2.7 LICORICE: LIne COntinuum Radiative tranferIntegrated Computing Engine (S. Baek, B. Semelin, F.Combes)
The LICORICE code has three main components: TreeSPH tocompute gravity and hydrodynamics, continuum radiative trans-fer with hydrogen and helium ionization physics, and Lyman-alphaline transfer. The latter is not relevant to this comparison and hasbeen described elsewhere. The ionizing continuum transfer hasbeen described in details in Baek et al. (2009).The current version of LICORICE does not include H
2
for-mation, or diffuse radiation from recombinations, but they will beincorporated in the future. LICORICE uses SPH particles for thegas dynamics and an adaptive grid for the radiative transfer. Physi-cal quantities are interpolated from one to the other as required.The ﬂuid dynamics are followed using a TreeSPH method.The implementation is described in detail in Semelin & Combes(2002) and Semelin & Combes (2005). Since there are many va-rieties of SPH, we summarize the main features of our algorithmhere. We use a spherically-symmetric spline-smoothing kernel and50 neighbours to compute the SPH quantities using an arithmeticaverage between the neighbours of the smoothing length
h
and thesimple viscosity scheme by Monaghan (1992).For the tests in this paper we implemented transmissiveboundary conditions. This was achieved as follows: for each SPHparticle within a distance of the simulation box boundary smallerthan its smoothing length
h
, we create a symmetrical ’ghost’ par-ticle on the other side of the boundary. All physical quantities forthis ghost particle are equal to those of the initial particle, includingthe velocity. The ghost particles are used as neighbours to computethe SPH quantities of real particles. The ghost particles are erasedand recreated at each time step.The continuum radiative transfer is solved using a MonteCarlo approach similar to the one employed in the CRASH code(Maselli et al. 2003). Here we summarize only the differences be-tween LICORICE and CRASH. We compute the gas density ateach particle’s position with the SPH smoothing kernel, and phys-ical quantities such as ionization fraction and temperature are up-dated according to these particle densities. The density ﬁeld is gen-erally smooth, but may sometimes show spurious ﬂuctuations if theparticle number density changes sharply. This is a well known butunavoidable problem with SPH.The radiation ﬁeldisdiscretized into photon packets and prop-agated through cells along directions chosen at random. The cellsform an adaptive grid which is derived from the tree structure of the particle distribution. Our adaptive grid is built to keep the num-ber of particles in each cell within a given range (1 to 8 and 1 to 1ranges have been used). This yields greater resolution in the denser

Related Search

Similar documents

We Need Your Support

Thank you for visiting our website and your interest in our free products and services. We are nonprofit website to share and download documents. To the running of this website, we need your help to support us.

Thanks to everyone for your continued support.

No, Thanks