Multi-scale initial conditions for cosmological simulations

We discuss a new algorithm to generate multi-scale initial conditions with multiple levels of refinements for cosmological "zoom-in" simulations. The method uses an adaptive convolution of Gaussian white noise with a real space transfer function kernel together with an adaptive multi-grid Poisson solver to generate displacements and velocities following first (1LPT) or second order Lagrangian perturbation theory (2LPT). The new algorithm achieves RMS relative errors of order 10^(-4) for displacements and velocities in the refinement region and thus improves in terms of errors by about two orders of magnitude over previous approaches. In addition, errors are localized at coarse-fine boundaries and do not suffer from Fourier-space induced interference ringing. An optional hybrid multi-grid and Fast Fourier Transform (FFT) based scheme is introduced which has identical Fourier space behaviour as traditional approaches. Using a suite of re-simulations of a galaxy cluster halo our real space based approach is found to reproduce correlation functions, density profiles, key halo properties and subhalo abundances with per cent level accuracy. Finally, we generalize our approach for two-component baryon and dark-matter simulations and demonstrate that the power spectrum evolution is in excellent agreement with linear perturbation theory. For initial baryon density fields, it is suggested to use the local Lagrangian approximation in order to generate a density field for mesh based codes that is consistent with Lagrangian perturbation theory instead of the current practice of using the Eulerian linearly scaled densities.Comment: 22 pages, 24 figures. MNRAS in press. Updated affiliation

Coupled atmosphere-wildland fire modeling with WRF-Fire

We describe the physical model, numerical algorithms, and software structure of WRF-Fire. WRF-Fire consists of a fire-spread model, implemented by the level-set method, coupled with the Weather Research and Forecasting model. In every time step, the fire model inputs the surface wind, which drives the fire, and outputs the heat flux from the fire into the atmosphere, which in turn influences the atmosphere. The level-set method allows submesh representation of the burning region and flexible implementation of various ignition modes. WRF-Fire is distributed as a part of WRF and it uses the WRF parallel infrastructure for parallel computing.Comment: Version 3.3, 41 pages, 2 tables, 12 figures. As published in Discussions, under review for Geoscientific Model Developmen

Bayesian Strong Gravitational-Lens Modeling on Adaptive Grids: Objective Detection of Mass Substructure in Galaxies

We introduce a new adaptive and fully Bayesian grid-based method to model strong gravitational lenses with extended images. The primary goal of this method is to quantify the level of luminous and dark-mass substructure in massive galaxies, through their effect on highly-magnified arcs and Einstein rings. The method is adaptive on the source plane, where a Delaunay tessellation is defined according to the lens mapping of a regular grid onto the source plane. The Bayesian penalty function allows us to recover the best non-linear potential-model parameters and/or a grid-based potential correction and to objectively quantify the level of regularization for both the source and the potential. In addition, we implement a Nested-Sampling technique to quantify the errors on all non-linear mass model parameters -- ... -- and allow an objective ranking of different potential models in terms of the marginalized evidence. In particular, we are interested in comparing very smooth lens mass models with ones that contain mass-substructures. The algorithm has been tested on a range of simulated data sets, created from a model of a realistic lens system. One of the lens systems is characterized by a smooth potential with a power-law density profile, twelve include a NFW dark-matter substructure of different masses and at different positions and one contains two NFW dark substructures with the same mass but with different positions. Reconstruction of the source and of the lens potential for all of these systems shows the method is able, in a realistic scenario, to identify perturbations with masses >=10^7 solar mass when located on the Einstein ring. For positions both inside and outside of the ring, masses of at least 10^9 solar mass are required (...).Comment: 21 pages, 15 figures, 4 tables; accepted for publication in MNRA

Performance comparison of point and spatial access methods

In the past few years a large number of multidimensional point access methods, also called multiattribute index structures, has been suggested, all of them claiming good performance. Since no performance comparison of these structures under arbitrary (strongly correlated nonuniform, short "ugly") data distributions and under various types of queries has been performed, database researchers and designers were hesitant to use any of these new point access methods. As shown in a recent paper, such point access methods are not only important in traditional database applications. In new applications such as CAD/CIM and geographic or environmental information systems, access methods for spatial objects are needed. As recently shown such access methods are based on point access methods in terms of functionality and performance. Our performance comparison naturally consists of two parts. In part I we w i l l compare multidimensional point access methods, whereas in part I I spatial access methods for rectangles will be compared. In part I we present a survey and classification of existing point access methods. Then we carefully select the following four methods for implementation and performance comparison under seven different data files (distributions) and various types of queries: the 2-level grid file, the BANG file, the hB-tree and a new scheme, called the BUDDY hash tree. We were surprised to see one method to be the clear winner which was the BUDDY hash tree. It exhibits an at least 20 % better average performance than its competitors and is robust under ugly data and queries. In part I I we compare spatial access methods for rectangles. After presenting a survey and classification of existing spatial access methods we carefully selected the following four methods for implementation and performance comparison under six different data files (distributions) and various types of queries: the R-tree, the BANG file, PLOP hashing and the BUDDY hash tree. The result presented two winners: the BANG file and the BUDDY hash tree. This comparison is a first step towards a standardized testbed or benchmark. We offer our data and query files to each designer of a new point or spatial access method such that he can run his implementation in our testbed

Latent Gaussian modeling and INLA: A review with focus on space-time applications

Bayesian hierarchical models with latent Gaussian layers have proven very flexible in capturing complex stochastic behavior and hierarchical structures in high-dimensional spatial and spatio-temporal data. Whereas simulation-based Bayesian inference through Markov Chain Monte Carlo may be hampered by slow convergence and numerical instabilities, the inferential framework of Integrated Nested Laplace Approximation (INLA) is capable to provide accurate and relatively fast analytical approximations to posterior quantities of interest. It heavily relies on the use of Gauss-Markov dependence structures to avoid the numerical bottleneck of high-dimensional nonsparse matrix computations. With a view towards space-time applications, we here review the principal theoretical concepts, model classes and inference tools within the INLA framework. Important elements to construct space-time models are certain spatial Mat\'ern-like Gauss-Markov random fields, obtained as approximate solutions to a stochastic partial differential equation. Efficient implementation of statistical inference tools for a large variety of models is available through the INLA package of the R software. To showcase the practical use of R-INLA and to illustrate its principal commands and syntax, a comprehensive simulation experiment is presented using simulated non Gaussian space-time count data with a first-order autoregressive dependence structure in time

Block Structured Adaptive Mesh and Time Refinement for Hybrid, Hyperbolic + N-body Systems

We present a new numerical algorithm for the solution of coupled collisional and collisionless systems, based on the block structured adaptive mesh and time refinement strategy (AMR). We describe the issues associated with the discretization of the system equations and the synchronization of the numerical solution on the hierarchy of grid levels. We implement a code based on a higher order, conservative and directionally unsplit Godunov's method for hydrodynamics; a symmetric, time centered modified symplectic scheme for collisionless component; and a multilevel, multigrid relaxation algorithm for the elliptic equation coupling the two components. Numerical results that illustrate the accuracy of the code and the relative merit of various implemented schemes are also presented.Comment: 40 pages, 10 figures, JPC in press. Extended the code test section, new convergence tests, several typos corrected. Full resolution version available at http://www.exp-astro.phys.ethz.ch/miniati/charm.pd

Object-oriented construction of a multigrid electronic-structure code with Fortran 90

We describe the object-oriented implementation of a higher-order finite-difference density-functional code in Fortran 90. Object-oriented models of grid and related objects are constructed and employed for the implementation of an efficient one-way multigrid method we have recently proposed for the density-functional electronic-structure calculations. Detailed analysis of performance and strategy of the one-way multigrid scheme will be presented.Comment: 24 pages, 6 figures, to appear in Comput. Phys. Com