Issue 
A&A
Volume 506, Number 2, November I 2009



Page(s)  1065  1070  
Section  Numerical methods and codes  
DOI  https://doi.org/10.1051/00046361/200912483  
Published online  18 August 2009 
A&A 506, 10651070 (2009)
Twodimensional adaptive mesh refinement simulations of colliding flows
M. Niklaus^{1,3}  W. Schmidt^{2,3}  J. C. Niemeyer^{2,3}
1  Deutsches Fernerkundungsdatenzentrum, Deutsches Zentrum für Luft
und Raumfahrt, Oberpfaffenhofen, Germany
2  Institut für Astrophysik, Universität Göttingen,
FriedrichHundPlatz 1, 37077 Göttingen, Germany
3  Lehrstuhl für Astronomie, Institut für Theoretische Physik und
Astrophysik, Universität Würzburg, Am Hubland, 97074 Würzburg, Germany
Received 13 May 2009 / Accepted 5 August 2009
Abstract
Context. Colliding flows are a commonly used
scenario for the formation of molecular clouds in numerical
simulations. Turbulence is produced by cooling, because of the thermal
instability of the warm neutral medium.
Aims. We carried out a twodimensional numerical
study of colliding flows to test whether statistical properties
inferred from adaptive mesh refinement (AMR) simulations are robust
with respect to the applied refinement criteria.
Methods. We compare probability density functions of
various quantities, as well as the clump statistics and fractal
dimension of the density fields in AMR simulations to a
staticgrid simulation. The static grid with 2048^{2} cells
matches the resolution of the most refined subgrids in the
AMR simulations.
Results. The density statistics are reproduced
fairly well by AMR. Refinement criteria based on the cooling time or
the turbulence intensity appear to be superior to the standard
technique of refinement by overdensity. Nevertheless, substantial
differences in the flow structure become apparent.
Conclusions. In general, it is difficult to separate
numerical effects from genuine physical processes in
AMR simulations.
Key words: hydrodynamics  turbulence  instabilities  ISM: kinematics and dynamics  methods: numerical  ISM: clouds
1 Introduction
Computational fluid dynamics in astrophysics rely on numerical methods that are capable of covering a huge range of scales. Apart from smoothed particle hydrodynamics (Monaghan 1992), adaptive mesh refinement (AMR) has been applied to variety of problems. This method was developed by Berger & Oliger (1984) and Berger & Colella (1989). Among the widely used, publicly available AMR codes for astrophysical fluid dynamics are FLASH (Fryxell et al. 2000), Enzo (O'Shea et al. 2004) and Ramses (Teyssier 2002). Although there are comparative studies of AMR vs. SPH (for example, O'Shea et al. 2005; Agertz et al. 2007; Commerçon et al. 2008), the degree of reliance of AMR in comparison to nonadaptive methods has met only little attention so far.
Especially for turbulent flows, it is a nontrivial question whether the solutions obtained from AMR simulations agree with the correct solutions of the fluid dynamical equations at a given resolution level. For this reason, we systematically compare AMR and staticgrid simulations for a particular test problem in this article. We chose a scenario that has been investigated in the context of molecular cloud formation, namely, the frontal collision of opposing flows of warm atomic hydrogen at supersonic speed (Hennebelle et al. 2008; Walder & Folini 2000; Hennebelle & Audit 2007a; Heitsch et al. 2006; VázquezSemadeni et al. 2007). Because of the cooling instability at densities and temperatures of a few thousand Kelvin, the gas becomes highly turbulent at the collision interface. Since the instabilities develop on length scales much smaller than the integral scale, this problem is computationally extremely demanding. The twodimensional resolution study of Hennebelle & Audit (2007a) showed that the properties of the turbulent multiphase medium evolving in these simulations is highly resolutiondependent, and numerical convergence is seen only at resolutions well above 1000^{2}. In threedimensional simulations, such high resolutions are infeasible if static grids are used. Consequently, Hennebelle et al. (2008) and Banerjee et al. (2008) applied refinement by fixed density thresholds and refinement by Jeans mass, respectively, in their threedimensional highresolution AMR simulations.
In this article, we consider twodimensional colliding flows without selfgravity and magnetic fields for a systematic comparison of AMR simulations to a reference simulation on a static grid. We analyze both statistical properties and the morphology of the gas fragmentation due to the cooling instability. This work is organized as follows: in Sect. 2 the numerical methods are described and the setup of the simulations will be presented in detail. In Sect. 3, we compare the results from the different simulations. Section 4 concludes this paper with a summary of the main results and general remarks on AMR.
2 Numerical methods and simulation setup
The simulations presented in this article are accomplished using the open source code Enzo (O'Shea et al. 2004; Bryan & Norman 1997). The compressible Euler equations are solved by means of the staggered grid, finite difference method Zeus (Stone & Norman 1992a,b; Stone et al. 1992). We included the cooling function defined by Audit & Hennebelle (2005) in these equations:
The primitive variables are the mass density , the velocity and the specific total energy e of the fluid. The total energy per unit mass is given by
(4) 
where is the adiabatic exponent and the pressure P is related to the mass density and the temperature T via the ideal gas law:
(5) 
The constants , and denote the Boltzmann constant, the mean molecular weight and the mass of the hydrogen atom, respectively. The gas is assumed to be a perfect gas with and .
Figure 1: Phase diagramm of vs. . The solid black curve indicates the equilibrium curve defined by , and the straight lines (orange in online version) are the isothermals for 7000 K and 50 K. The contour shading indicates the probability density. 

Open with DEXTER 
The cooling function of Audit & Hennebelle (2005) includes the cooling by finestructure lines of CII and OI, further the cooling by H (Lyline) and the electron recombination onto positively charged grains. The heating is due to the photoelectric effect on small grains and polycyclic aromatic hydrocarbons (PAH) caused by the farultraviolet galactic background radiation. For more information about this cooling function see Bakes & Tielens (1994); Wolfire et al. (1995,2003); Spitzer (1978) and Habing (1968). The pressureequillibrium curve resulting from the cooling function is plotted as black curve in Fig. 1. For the numerical solution of the fluid dynamical equations, we used the radiative cooling routine implemented into Enzo. For each hydrodynamical time step, the state variables are iterated over several subcycles, and the resulting total energy increment for the whole time step is added.
Figure 2: Contour plot of the mass density after 5 Myr of evolution in the staticgrid simulation. The density scale is logarithmic. 

Open with DEXTER 
For our numerical study, the twodimensional setup of Audit & Hennebelle (2005) and Hennebelle & Audit (2007b) was adopted with small modifications. The initial gas content corresponds to the warm neutral material (WNM) in the ISM, i.e., the temperature is T=7100 K, the pressure is P=7 10^{13} erg cm^{3} and the number density of neutral hydrogen is n=0.71 cm^{3}. From the left and the right boundaries, warm gas with identical thermodynamic properties flows into the computational domain, where the cosineshaped inflow velocity profile is modulated with small perturbations, realised by randomly shifted phases. These phase shifts are kept constant for the different simulations, so that initial conditions are exactly the same for all runs to ensure comparability. Following Hennebelle et al. (2008), the top and bottom boundary conditions are periodic. The physical dimensions of the computational domain are 20 20 pc. The two inflows of gas collide in the middle of the domain. The supersonic collision causes a steep rise in the gas density that triggers the thermal instability, and gas undergoes a transition into the phase of the cold neutral material (CNM) in the ISM. In this phase, the gas has temperatures in the range 30100 K and number densities within 2050 cm^{3} (Ferrière 2001). The thermal instability produces highly turbulent structures (see Fig. 2) with Mach numbers up to 20. The challenge for AMR is to track these turbulent structures as accurately as possible.
Figure 3: Pdfs of the number density n ( upper panel) and the temperature T ( lower panel) in loglogscale for the different AMRcriteria (black curves) compared to the static grid simulation (red curves in the online version). 

Open with DEXTER 
A reference simulation was run with a static grid of 2048^{2} cells. Then the same setup was evolved in AMR simulations with a rootgrid resolution of 128^{2} cells and 4 levels of refinement. The resolution between adjacent refinement levels increases by a factor of 2. Hence, the effective resolution at the highest level of refinement is 2048^{2}. In these simulations, we employed three different types of refinement criteria:
 1.
 refinement by overdensity (OD);
 2.
 refinement by cooling time (CT);
 3.
 refinement by rate of compression and enstrophy (RCEN).
For comparison of the simulation results, we calculated probability density functions (pdf) of several quantities. To analyze the gas fragmentation in each simulation, we adapted the clumpfind algorithm implemented by Padoan et al. (2007) to nonisothermal problems. The algorithm identifies the smallest dense regions that fulfill the Jeans criterion for gravitationally unstable gas. Since the clump samples found on the twodimensional grids used in our simulations are insufficient for the calculation of clump mass spectra, only the total number and the mean size of the clumps are used for quantitative comparisons. In addition, we computed the fractal dimension of gas at densities higher than (corresponding to the minimum density of gas in the cold phase) by means of the boxcounting method (Federrath et al. 2009).
3 Results
Due to the gradual accumulation of gas in the simulation domain, no strict statistical equilibrium is approached. For this reason, we evolved the flow until noticeable smallscale structure has developed and the separation of the gas into two phases has emerged. As shown in Fig. 1, two distinct phases are found at time . At this time, the central flow region is in a turbulent state (a contour plot of the mass density of the gas is shown in Fig. 2). Thus, we carry out our analysis for . While the main fraction of the gas is situated in the warm phase with temperatures between 5000 and 10 000 K and low densities cm^{3}, the cold gas with temperatures between 30 K and 100 K and densities in the range 30350 cm^{3} can be found close to the equilibrium curve.
Figure 4: Clump distributions ( upper panel) and gas with number density ( lower panel) for the AMR simulations compared to the static grid simulation. The clumps (displayed in red in the online version) were identified by a clumpfind algorithm and are superimposed on a grayscale density plot. 

Open with DEXTER 
The pdfs of the mass density and the temperature obtained from different AMR simulations are plotted in Fig. 3. In principle, all refinement criteria reproduce the distributions found in the staticgrid simulation quite well, although there is a trend of slightly more cold gas at the cost of warm gas. The discrepancy is more pronounced for refinement by overdensity (OD) compared to the other criteria, and it becomes worse for the higher density thresholds (OD3 and OD4; not shown in the figure). Nevertheless, it appears that the thermodynamic properties of the gas are quite robust in AMR simulations.
The gravitationally unstable clumps of gas identified by the clumpfind algorithm in the staticgrid simulation at time are depicted in Fig. 4a. The corresponding results for the AMR runs are shown in Figs. 4ad. Table 1 lists the total number and the mean size of the clumps for each simulation. Also listed are the fractal dimensions of the gas regions with number density cm^{3}, which are plotted in Figs. 4eh.
Table 1: Number N and average size of the clumps in the CNM; fractal dimension D of gas with number density greater than 20 cm^{3}.
For refinement by OD, the fragmentation of the CNM is severely underestimated. The number of clumps is roughly half of the number in the staticgrid simulation, and the clumps are typically larger. The lower degree of cold gas fragmentation results in a smaller fractal dimension (also see Fig. 4f). If the criteria OD3 and OD4 are applied, the number of clumps decreases further, while their average size increases. In the case of criterion OD4, a slightly higher fractal dimension is obtained, because the cold phase tends to fill broad, areafilling regions. The cooling time criterion CT yields an amount of dense clumps with an average size that compares well to the reference simulation (see Fig. 4c), although the degree of fragmentation appears to be overestimated slightly. However, we found that this overestimation decreases with the further evolution of the colliding flows and, thus, appears to be transient. Refinement by RCEN also reproduces the number of clumps and the fractal dimension of dense gas very well. However, there are some anomalously big clumps, which contribute to an average clump size that is systematically too large. In the plot showing gas at density cm^{3} (see Fig. 4h), on the other hand, such anomalous structures are not visible. Although refinement by RCEN does not overproduce gas in the cold phase (as one can see from the excellent agreement of the density and temperature pdfs in Figs. 3c and f), there appears to be a bias toward bigger clumps with this refinement method.
In contrast to the phase separation and gas fragmentation, striking deviations of the turbulent flow properties in the AMR vs. static grid simulations become apparent. Generally, a lot of turbulent smallscale structure is missing in the AMR simulations. Even for the criterion RCEN, which is based on control variables related to turbulence, this is apparent from the contour plots of the squared vorticity modulus shown in Fig. 5. Basically, the perturbations of the velocity field imposed at the inflow boundaries are quickly smoothed out in AMR simulations, so that turbulence is only produced by secondary (e.g., KelvinHelmholtz) instabilities at the collision interface in the central region of the computational domain. The reason is that all AMR cirteria, including RCEN, select relatively large fluctuations, whereas smaller perturbations are suppressed. On a static grid, on the other hand, the perturbations are transported from the boundaries to the centre and actively contribute to the production of turbulence. Consequently, small eddies are present in almost the whole domain in this case. Accordingly, the probability distribution of vorticity is markedly different (see Fig. 6). In contrast, Schmidt et al. (2009) found very close agreement of the vorticity pdfs in a staticgrid and an AMR simulation with criterion RCEN for turbulence in a periodic box with largescale forcing. Our results thus indicate that the merits of different refinement schemes are nonuniversal but rather depend on the properties of individual flow structures.
Figure 5: Plots of the squared vorticity in the staticgrid simulation a) and in the AMR simulation with refinement criterion RCEN b), where the contours of the refinement levels are shown (colorcoded in the online version). 

Open with DEXTER 
Figure 6: Pdfs of the vorticity modulus for AMRcriterion RCEN (black curve) compared to the static grid simulation (red curve in the online version). 

Open with DEXTER 
4 Conclusions
We performed twodimensional simulations of colliding flows of warm atomic hydrogen with a radiative cooling function as source term in the energy equation. The goal of our study was the systematic comparison of AMR simulations, where different criteria for refinement were applied, to a reference simulation on a static grid.
While the probability distributions of mass density and temperature are well reproduced in AMR simulations, regardless of the refinement technique, differences become apparent in the fragmentation properties of the cold gas phase. As indicators, we used the total number of clumps and their average size. The clumps were identified by a clumpfind algorithm. In addition, we calculated the fractal dimension of dense gas, assuming a number density threshold of 20 cm^{3}. Remarkably, the largest deviations from the clump statistics and fractal dimension extracted from the staticgrid simulation, were encountered for refinement by overdensity, which is a commonly used refinement criterion in astrophysical AMR simulations. The deviations increase with the chosen density threshold. In this regard, it is important to note that Hennebelle et al. (2008) applied a densitybased refinement criterion, where the thresholds were chosen even higher than those considered in our study. Good agreement, on the other hand, was obtained if the cooling time or enstrophy in combination with the rate of compression (the negative rate of change of the velocity divergence) were applied.
Table 2: CPU time for the simulation runs.
Substantial problems with AMR became apparent with regard to turbulent flow properties. Basically, none of our AMR runs were able to reproduce even remotely the smallscale structure of turbulence and the probability distributions of turbulent flow variables such as the vorticity modulus. This defficiency can be attributed to the selection effects introduced by adaptive techniques. The definition of thresholds for triggering refinement either selects strong local fluctuations (for example, large shear that gives rise to KelvinHelmholtz instabilities) or largescale perturbations such as accumulations of mass that become Jeansunstable in selfgravitating gas. In this respect, the test problem we investigated in this work is particularly tough, because turbulence stems from smallscale instabilities that are seeded by weak initial perturbations. The varying grid resolution in AMR simulations inevitably modulate the growth of these instabilities and, as a consequence, the production of turbulence is suppressed. This defficiency might be overcome by the application of a subgrid scale model, which transports turbulent energy contained in small eddies that are resolved on finer grids across coarser grid regions (see Maier et al. 2009).
The key point of using AMR is the computational cost for a given effective resolution. Indeed, Table 2 demonstrates that a substantial reduction of computation time is achieved with AMR, especially if refinement by overdensity is applied. So, AMR is essentially a tradeoff, where fast computation has to be weighted carefully against the question whether the essential physics of the specific problem is captured. Regarding threedimensional simulations with high spatial resolution the specific impacts of the numerics are hard to predict, since the direct comparison to a high resolution computation on a static grid is not feasible and thus not available so far. Nevertheless, the physical results could be affected as much as in our twodimensional simulations. So, besides the enormeous reduction in computational time especially in 3D, caution is to be recommended when discussing the results gained from AMR simultions and the conclusions drawn from them.
AcknowledgementsWe thank Patrick Hennebelle and Edouard Audit for providing the cooling function that was used for this numerical study. We also thank Paolo Padoan for sharing his clumpfind tool. Computations described in this work were performed using the Enzo code developed by the Laboratory for Computational Astrophysics at the University of California in San Diego (http://lca.ucsd.edu).
References
 Agertz, O., Moore, B., Stadel, J., et al. 2007, MNRAS, 380, 963 [NASA ADS] [CrossRef]
 Audit, E., & Hennebelle, P. 2005, A&A, 433, 1 [NASA ADS] [CrossRef] [EDP Sciences]
 Bakes, E. L. O., & Tielens, A. G. G. M. 1994, ApJ, 427, 822 [NASA ADS] [CrossRef]
 Banerjee, R., VazquezSemadeni, E., Hennebelle, P., & Klessen, R. 2008 [arXiv:0808.0986]
 Berger, M. J., & Colella, P. 1989, J. Comp. Phys., 82, 64 [NASA ADS] [CrossRef]
 Berger, M. J., & Oliger, J. 1984, J. Comp. Phys., 53, 484 [NASA ADS] [CrossRef]
 Bryan, G. L., & Norman, M. L. 1997, in Computational Astrophysics; 12th Kingston Meeting on Theoretical Astrophysics, ed. D. A. Clarke, & M. J. West, ASP Conf. Ser., 123, 363
 Commerçon, B., Hennebelle, P., Audit, E., Chabrier, G., & Teyssier, R. 2008, A&A, 482, 371 [NASA ADS] [CrossRef] [EDP Sciences]
 Federrath, C., Klessen, R. S., & Schmidt, W. 2009, ApJ, 692, 364 [NASA ADS] [CrossRef]
 Ferrière, K. M. 2001, Rev. Mod. Phys., 73, 1031 [NASA ADS] [CrossRef]
 Fryxell, B., Olson, K., Ricker, P., et al. 2000, ApJS, 131, 273 [NASA ADS] [CrossRef]
 Habing, H. J. 1968, Bull. Astron. Inst. Netherlands, 19, 421 [NASA ADS]
 Heitsch, F., Slyz, A. D., Devriendt, J. E. G., Hartmann, L. W., & Burkert, A. 2006, ApJ, 648, 1052 [NASA ADS] [CrossRef]
 Hennebelle, P., & Audit, E. 2007a, A&A, 465, 431 [NASA ADS] [CrossRef] [EDP Sciences]
 Hennebelle, P., & Audit, E. 2007b, A&A, 465, 431 [NASA ADS] [CrossRef] [EDP Sciences]
 Hennebelle, P., Banerjee, R., VázquezSemadeni, E., Klessen, R. S., & Audit, E. 2008, A&A, 486, L43 [NASA ADS] [CrossRef] [EDP Sciences]
 Maier, A., Iapichino, L., Schmidt, W., & Niemeyer, J. 2009, ApJ, submitted
 Monaghan, J. J. 1992, ARA&A, 30, 543 [NASA ADS] [CrossRef]
 O'Shea, B. W., Bryan, G., Bordner, J., et al. 2004, arXiv Astrophysics eprints
 O'Shea, B. W., Nagamine, K., Springel, V., Hernquist, L., & Norman, M. L. 2005, ApJS, 160, 1 [NASA ADS] [CrossRef]
 Padoan, P., Nordlund, Å., Kritsuk, A. G., Norman, M. L., & Li, P. S. 2007, ApJ, 661, 972 [NASA ADS] [CrossRef]
 Schmidt, W., Federrath, C., Hupp, M., Kern, S., & Niemeyer, J. C. 2009, A&A, 494, 127 [NASA ADS] [CrossRef] [EDP Sciences]
 Spitzer, L. 1978, Physical Processes in the Interstellar Medium, 1st edn. (WileyInterscience)
 Stone, J. M., & Norman, M. L. 1992a, ApJS, 80, 753 [NASA ADS] [CrossRef]
 Stone, J. M., & Norman, M. L. 1992b, ApJS, 80, 791 [NASA ADS] [CrossRef]
 Stone, J. M., Mihalas, D., & Norman, M. L. 1992, ApJS, 80, 819 [NASA ADS] [CrossRef]
 Teyssier, R. 2002, A&A, 385, 337 [NASA ADS] [CrossRef] [EDP Sciences]
 VázquezSemadeni, E., Gómez, G. C., Jappsen, A. K., et al. 2007, ApJ, 657, 870 [NASA ADS] [CrossRef]
 Walder, R., & Folini, D. 2000, Ap&SS, 274, 343 [NASA ADS] [CrossRef]
 Wolfire, M. G., Hollenbach, D., McKee, C. F., Tielens, A. G. G. M., & Bakes, E. L. O. 1995, ApJ, 443, 152 [NASA ADS] [CrossRef]
 Wolfire, M. G., McKee, C. F., Hollenbach, D., & Tielens, A. G. G. M. 2003, ApJ, 587, 278 [NASA ADS] [CrossRef]
All Tables
Table 1: Number N and average size of the clumps in the CNM; fractal dimension D of gas with number density greater than 20 cm^{3}.
Table 2: CPU time for the simulation runs.
All Figures
Figure 1: Phase diagramm of vs. . The solid black curve indicates the equilibrium curve defined by , and the straight lines (orange in online version) are the isothermals for 7000 K and 50 K. The contour shading indicates the probability density. 

Open with DEXTER  
In the text 
Figure 2: Contour plot of the mass density after 5 Myr of evolution in the staticgrid simulation. The density scale is logarithmic. 

Open with DEXTER  
In the text 
Figure 3: Pdfs of the number density n ( upper panel) and the temperature T ( lower panel) in loglogscale for the different AMRcriteria (black curves) compared to the static grid simulation (red curves in the online version). 

Open with DEXTER 
In the text
Figure 4: Clump distributions ( upper panel) and gas with number density ( lower panel) for the AMR simulations compared to the static grid simulation. The clumps (displayed in red in the online version) were identified by a clumpfind algorithm and are superimposed on a grayscale density plot. 

Open with DEXTER  
In the text 
Figure 5: Plots of the squared vorticity in the staticgrid simulation a) and in the AMR simulation with refinement criterion RCEN b), where the contours of the refinement levels are shown (colorcoded in the online version). 

Open with DEXTER 
In the text
Figure 6: Pdfs of the vorticity modulus for AMRcriterion RCEN (black curve) compared to the static grid simulation (red curve in the online version). 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.