HyPhy: Deep Generative Conditional Posterior Mapping of Hydrodynamical Physics

Generating large volume hydrodynamical simulations for cosmological observables is a computationally demanding task necessary for next generation observations. In this work, we construct a novel fully convolutional variational auto-encoder (VAE) to synthesize hydrodynamic fields conditioned on dark matter fields from N-body simulations. After training the model on a single hydrodynamical simulation, we are able to probabilistically map new dark matter only simulations to corresponding full hydrodynamical outputs. By sampling over the latent space of our VAE, we can generate posterior samples and study the variance of the mapping. We find that our reconstructed field provides an accurate representation of the target hydrodynamical fields as well as a reasonable variance estimates. This approach has promise for the rapid generation of mocks as well as for implementation in a full Bayesian inverse model of observed data.


INTRODUCTION
Understanding the large-scale structure of the universe requires simultaneous analysis of both the evolution of the underlying dark matter cosmic web and the complex hydrodynamics leading to the formation of biased tracers. Over the past thirty years, hydrodynamical simulations have become the standard tool to generate mock observable data which includes both of these effects (Evrard 1990;Cen 1992;Katz et al. 1996;Springel 2005). However, the power of these hydrodynamical simulations comes with significant computational cost, and the next generation of cosmological surveys will require unprecedented precision across a wide range of scales (e.g. Walther et al. (2021)). In this regime, computing quantities like covariance matrices (which require large numbers of simulations) becomes an increasingly daunting task, so there is a clear need for approximate methods that can ease some of the computational burden.
bhorowitz@princeton.edu Dornfest@berkeley.edu In recent years, machine learning techniques have emerged as promising surrogate models for complex hydrodynamics, as they can be used to rapidly generate hydrodynamic fields with remarkable perceptual and statistical fidelity. In Zamudio-Fernandez et al. (2019), the authors were able to generate realistic neutral hydrogen (HI) maps which reproduce the properties of hydrodynamical simulations over a range of scales. In Tröster et al. (2019), the authors used generative models to map from two dimensional dark matter maps to thermal Sunyaev Zeldovich (tSZ) maps. They were able to reproduce accurate tSZ summary statistics over a wide range of scales, given only the dark matter maps. Related work in Wadekar et al. (2020) used a more traditional feedforward architecture, HInet, to paint neutral hydrogen in all three dimensions, but this architecture does not allow exploration of posterior properties and uncertainties.
Estimation of the uncertainty of a neural network's output is critical in order to propagate errors accurately for cosmological and cosmographical analysis. However, within the astronomical community there has been relatively little work in error analysis in the context of these neural-network-based surrogate models. A promising approach in the case of low dimensional data is to make the network output be a multidimensional Gaussian distribution as opposed to a single point estimation, i.e. a Gaussian mixture model (see, for example Tsang & Schultz (2019)). However, for high dimensional outputs this approach would have difficulty capturing the full covariance in a memory efficient way.
In this work, we will instead structure the latent space of a conditional variational autoencoder (C-VAE) to learn the uncertainty in the mapping from a dark matter map to the hydrodynamical quantities. The general structure of our network is inspired by style transfer machine learning literature (Johnson et al. 2016;Esser et al. 2018) where the latent space of the C-VAE is used to capture stylistic characteristics of the mapping. A similar network has recently been used to generate a realistic distribution of galaxy images (Lanusse et al. 2020).
As we expect hydrodynamic quantities to be quasilocal, we constrain our model to maintain this property by restricting the spatial field of view of the input and use convolutional layers of different sizes in order to capture information across a range of scales. When transforming to redshift space, as is needed to match observables, this locality is broken, so rather than directly modeling observable quantities we instead focus on reconstructing the underlying (real space) baryon density, temperature, and velocity fields. Then, estimates for the target observables can trivially be computed using existing analytical tools, which also allows flexibility in modeling physical details "orthogonal" to the dark matter and hydrodynamics relationship (e.g. atomic species ratios, ionization rates, etc). For this work we will focus on Lyα flux for Lyman Alpha Forest cosmology measurements, but our model outputs are generic and could be applied to generate other target observables.
This work is a companion work to Harrington et al. (2021) which used a deterministic network to perform this mapping. While this approach successfully recovers the key Lyα observables for next generation cosmological measurements, it has difficulty in capturing stochastic processes such as shocked regions. As the approach and goal of these works is quite different we present them as two separate papers.
The paper is organized as follows: in Sect. 2 we briefly describe nyx and the simulation data set used and then describe our neural network architecture. We present our results in Sect. 3, first reviewing MAP performance and posterior accuracy. We conclude Sect. 4 with our ongoing work and future areas of interest to the community.

Hydrodynamical Simulations
We choose to obtain simulation data from nyx, a massively parallel multiphysics code, because it was developed for simulations of the intergalactic medium and has been used for many recent IGM studies (Horowitz et al. 2019;Davies et al. 2019;Walther et al. 2019), and is capable of modeling dark matter and hydrodynamic evolution in great detail. The Nyx code (Almgren et al. 2013) follows the evolution of dark matter modeled as self-gravitating Lagrangian particles, while baryons are modeled as an ideal gas on a set of rectangular Cartesian grids. The Eulerian gas dynamics equations are solved using a second-order accurate piecewise parabolic method to accurately capture shocks. Besides solving for gravity and the Euler equations, we also include the main physical processes relevant for the Lyα forest. We consider the chemistry of the gas as having a primordial composition of hydrogen and helium, include inverse Compton cooling off the microwave background and keep track of the net loss of thermal energy resulting from atomic collisional processes (Lukić et al. 2015). All cells are assumed to be optically thin to ionizing radiation, and radiative feedback is accounted for via a spatially uniform, time-varying UV background radiation given to the code as a list of photoionization and photoheating rates (Haardt & Madau 2012). We note that this type of simulation is used as a forward model in virtually any recent inference work using Lyman alpha power spectrum (Boera et al. 2019;Walther et al. 2019;Palanque-Delabrouille et al. 2020;Rogers & Peiris 2020;Walther et al. 2021). Simulations of this kind neglect the effects of inhomogeneous reionization, which produces temperature and UV background fluctuations on large scales, especially at redshifts higher than those studied in this work (z 4).
In this work we used simulations of a standard ΛCDM cosmological model, consistent with the latest cosmological constraints from the CMB (Planck Collaboration et al. 2020): Ω m =0.31, Ω Λ =0.69, Ω b =0.0487, h=0.675, σ 8 =0.82 and n s =0.965. For the hydrogen and helium mass abundances we adopted values consistent with the CMB observations and Big Bang nucleosynthesis (Coc et al. 2013): X p =0.76 and Y p =0.24. Box size of simulations is 20 h −1 Mpc, with N = 1024 3 particles and grid cells, resulting in a resolution of ∼ 20h −1 kpc, fulfilling convergence criteria for a percent-level accurate Lyα quantities (Lukić et al. 2015).
In addition to hydrodynamical simulations, we have also produced N-body simulations starting with the same initial conditions. These neglect all other forces but gravity, and all matter is considered collisionless, although baryonic effects are imprinted in the initial power  Figure 1. Schematic diagram showing the flow of our u-net conditional variational autoencoder architecture. Our model can be viewed as four interconnected parts, a encoder for the dark matter fields, an encoder for the hydrodynamical fields, a variational block, and a decoder. The network is constructed such that after training the hydrodynamical block can be removed and the latent space sampled from a unit Gaussian to generate the corresponding hydrodynamical field. Plotted is only the dark matter density and the baryon temperature, but the model also is fed as an input the dark matter velocity and outputs the baryon density and line of sight velocity. spectrum for the total matter. Baryons only minorly affect dark matter clustering in the regime relevant for the Lyα forest, but we have nevertheless produced these "companion" simulations to maintain maximum reality in our modeling when we train to infer hydrodynamical quantities from a N-body run and avoid any possible back-reaction on the dark matter field. Throughout this work we use one set of hydro and matching N-body simulation for training and a different set for validation purpose. The two differ only in the random phases of the Gaussian random field in the initial conditions.

Conditional Variational Auto-encoder Model
As neural networks are becoming an increasingly wellknown tool in astrophysics and cosmology, we just briefly summarize our model 1 and highlight how it differs from others in the literature.
We use a conditional variational auto-encoder architecture to study the posterior of hydrodynamical quantities, τ , given the dark matter field, δ. The core concept of a variational auto-encoder network is that the neural 1 Our model is implemented in tensorflow, with sample data and a saved model available at https://github.com/bhorowitz/ HyPhy-public.
network learns the probability space described by the training sample by marginalizing over a (usually bottlenecked) set of latent space parameters, Z. In this case we are interested in constraining the outputs of our model to those corresponding to a given dark matter realization which we enforce by using the dark matter realization as a "condition"; a field given to the network both during the "encoding" step (where the fields are mapped to the latent space) and the "decoding" step (where the latent space is mapped to a 3D field). For finding a given hydrodynamical realization, we are interested in calculating the following quantity, where p(τ |Z, δ) is the generator network. In order to train the network, we also need to define a connected, overlapping, encoder network, q(Z|τ, δ). We can generalize the standard evidence lower bound (ELBO) straightforwardly to include this conditioning variable where E q is the expectation value marginalized over q.
In standard style-transfer implementations, p(z|δ) plays Figure 2. Diagram showing the broadcasting of the latent space samples into the same dimensionality as the downsampled dark matter field. This stack will comprise the filters passed to the upsampling convolutional layers in the hydrodynamics decoder. This allows changing input dark matter map sizes during generation while keeping the same dimensionality of the latent space and avoiding the need for dense layers. Note that in our model each filter is three dimensional, not two dimensional as shown here.
an important role as a prior over the latent space parameters (i.e. Esser et al. (2018)). In our case, we will not be imposing any direct interpretation to our latent space parameters, and we will find that this prior distribution, p(Z|δ), will be pulled to a unit-normal Gaussian due to the loss term. We again follow the standard derivation for CVAE networks and model the probability distributions, p(τ |δ, Z) and q(Z|δ, τ ), as neural networks with associated free parameters, G θ and F φ , respectively. In Figure 1, G θ corresponds to the combination of the dark matter encoder and hydrodynamics decoder, while F φ corresponds to the combination of dark matter encoder and hydrodynamics encoder. These overlapping neural networks weights, θ and φ, are trained jointly using the standard ADAM optimizer under the associated loss where KL is the Kullback-Leibler divergence to compare distributions. Assuming we treat the generator network, G θ , as deterministic in δ and Z, we can simplify the second term with our chosen reconstruction loss. There are many possible choices for this loss function, including L1, L2, perceptual-loss, or an adversarial loss. In our case we choose the standard L1 loss term, i.e.
whereτ is the simulated true field corresponding to δ. The KL-divergence term can also be further simplified by expressing our latent space image, q φ (Z|δ, τ ), as a multidimensional Gaussian with diagonal covariance using the renormalization trick and our target distribution as a unit normal Gaussian, where µ and σ are outputs of our encoding network F φ , and correspond to the means and standard deviations of our latent space distributions from which we sample.
The variations of the hydrodynamical model are captured by properties of the dark matter field at a variety of scales. We are therefore motivated by examples of image segmentation from machine learning literature to use an altered U-Net structure for our network. We summarize this structure in Figure 1. The most notable in the structure of this network is the skip connections across the bottleneck mapping the dark matter field from the encoder to decoder side. This is designed to maximize the possible information extracted from the dark matter field, in a computationally expedient fashion, as well as to minimize the dark matter dependent structure in the latent space. In particular since we draw samples from the latent space that are fully representative of uncertainty in the mapping, we want the conditional probability of p(Z|δ) to be as close to a uniform Gaussian, N (0, 1) as possible. Structure in the latent space will result in a biased posterior sample distribution. In abstract we could implement various techniques (such as an autoregressive flow, as in Lanusse et al. (2020)), to ensure that our latent space is well mixed, but in practice we find that there is little dependency on matter properties in the latent space.

Fully Convolutional Architecture
A critical feature of our architecture is for it to be fully convolutional, thereby allowing inputs of any size during inference while being trained on smaller volumes. Traditionally, VAEs utilize dense layers to upsample the drawn latent space variable, which are then possibly followed by transposed convolutions (or other upsampling convolutions) to reach the desired output size. This design restricts the input image size to always be of the same dimension as the training set. To avoid this constraint, we broadcast the latent space parameters into a feature map whose dimension matches that of the lowest level in the U-Net, then upsample from there as is standard (see Figure 2). Since the upsampling convolutional layers are inherently local (as determined by their kernel size) we maintain the desired locality of our network, while still allowing every element of the output to "see" the full latent space.
While generated volumes can be of any dimension, the training volumes are fixed in size due to the need for dense layers in the encoder to predict the latent space distributional parameters. However, training data can easily be segmented to the proper size, and the choice of the crop size is a problem-dependent hyperparameter to be tuned. Training on boxes of limited volume means long distance correlations not captured in the dark matter distribution (e.g. a spatially fluctuating UV background) would not be well reproduced with this architecture and would require additional considerations.
For the purpose of this work, we use training boxes of approximately 4 h −1 M pc.
Since neither each individual training box or standard convolutional layer implementations are periodic, we must minimize edge effects by restricting our loss function to only compare the central region of our training sample. In this work, we train on boxes with 64 voxel side-length, of which we ignore 10 voxels on each side in order to avoid dealing with edge effects.
In order to better utilize our limited training data, we use the symmetries of our system to randomly augment training samples by performing reflections and rotations of our box over each of the three spatial dimensions. In addition, our training boxes are overlapping in space, providing some knowledge of translational symmetry of the underlying physics model.

Redshift Information
For application to real Lyα forest data, a key aspect of the mapping is to include the redshift dependence of the mapping as this dependence is used for cosmological constraints (i.e. Chabanier et al. (2019)). To include The characteristics of moderately large shocks are difficult for our network to learn, resulting in a less sharp feature. However, this uncertainty in the reconstruction is captured when sampling over latent space, as we see significant variance in the shocked region. If we look at the individual posterior samples this is even more visually apparent, see Figure 8. this in our model, we condition on a redshift field of the same size as our training volume (i.e. every pixel has an associated redshift). This allows us to vary the redshift over the box (i.e. to generate light-cones), as well as train the same model to work across cosmic time. We train our model using snapshots from the same simulation at z = 2.4, z = 3.0 and z = 4.0, as this range dominates the cosmological Lyα forest signal (Walther et al. 2021).

RESULTS
We apply our trained network to a separate simulation (not used for training), focusing on the central redshift of z = 3.0, which requires only changing the latent space projection dimension while maintaining the same network weights. We first show the results for the maximum a posterior point predicted by the network, both in the base hydrodynamical quantities and in terms of derived Lyα flux. We then discuss the posterior properties of a representative sample distribution in Subsection 3.2.
The maximum a posterior (MAP) should correspond to the highest maximum of the distribution P ( Z), which by construction will be at Z = 0 for a fully trained model. A useful starting point for our analysis will be to see how this point in the posterior space behaves to judge the quality of our reconstruction for test boxes with the same dimension as our training set. Note that we expect these reconstructions to be slightly smoothed versus generic samples.
We show test boxes both without and with significant shocked regions in Figures 3 and 4, respectively. In low to medium density regions, we recover excellent qualitative agreement across a range of scales for baryon density, baryon velocity, and temperature. For the shocked example, Figure 8, the hyphy recovered temperature field has a significantly less prominent shocked region with a much smoother boundary; this is a common phenomenon VAE-type networks (Khan et al. 2018).
We also show the model predicted variance by examining 1000 samples drawn from N (0, 1), qualitatively showing that regions of highest variance are where there is the strongest disagreement between hyphy output and the simulated true field. In particular, regions near the boundaries of structures (like filaments) have signif- x (h 1 Mpc)) Figure 5. hyphy run on a box significantly larger than the training volumes. The fully convolutional nature of the network allows the network to be run on any box that is a multiple of a base dimension (currently 8 pixels) and has the same spatial resolution. There are no obvious artifacts caused by the increased image size.
icant variance as do those with significant astrophysical shocks. We discuss these properties more in Section 3.2.

Large Volume Statistics
Next, we apply hyphy to the entire test N-body simulation at once. To match the resolution of our training sample, we uniformly down-sample the volume by a factor of two. We show the reconstructed resulting baryon density, velocity, and temperature in Figure 5.

Gas Phase Physics
A well studied aspect of cosmological hydrodynamical simulations is the relation between the gas density and gas temperature (Gunn & Peterson 1965;Sorini et al. 2016). We show our reconstructed relation in log space in Figure 6. Following Ursino et al. (2010); Martizzi et al. (2019); Galárraga-Espinosa et al. (2020), this plot can be viewed as a phase-space distribution between Warm Hot Intergalactic Medium (WHIM), warm Circumgalactic Medium (wCGM), diffuse IGM, halo gas, and "hot gas". Halo gas consists of relatively "cool" gas including ISM within galaxies, as well as more diffuse gas found in-between galaxies within halos. Sim-ilarly, warm CGM is found in dense enviroments, but has been significantly heated via shocks or feedback processes near galaxies. Diffuse IGM and WHIM are found in less dense environments, such as regions surrounding filaments. The WHIM component is of significant interest due to the "missing baryon" problem (Fukugita et al. 1998). The last component, consisting of gas at any density at a temperature in excess of 10 7 K and generally associated with massive shocks, is a vanishingly small percentage of the test volume (68 pixels out of 512 3 total pixels) and is grouped with wCGM or WHIM depending on density. In addition, some studies (e.g. Martizzi et al. (2019)) separate star-forming gas, with densities log(n H ) > −1.0, as a separate phase. Since nyx does not model star formation we do not include this phase in our analysis.
In Figure 6, we summarize the recovered volume fractions vs. the simulated truth. We find excellent recovery of the diffuse IGM and WHIM, with slightly worse performance of the wCGM component. Halo Gas, constituting a tiny fraction of the total volume, is very difficult to recover accurately, with hyphy finding a factor 19× more halo gas. However, it is important to note the vast majority of this excess is on the border of the diffuse IGM component in phase space and can probably be better described as over-estimated IGM temperature as opposed to misattributed halo gas.
A key property of interest to the IGM is the power law relation between density and temperature, which can be related to statistics of the Lyα forest (Hui & Gnedin 1997). Following the procedure in Sorini et al. (2016), we identify the gas around two bins centered at log(∆ b,0 ) = −1 and log(∆ b,1 ) = 0, with a width of 5% around the central value. We calculate the median temperature of the corresponding particles in each bin, and use those two points to determine the powerlaw relation T = T 0 ∆ γ b . For nyx, we find (T 0 , γ) = (10 4.08 K, 1.53) and for the hyphy reconstruction we find (T 0 , γ) = (10 4.08 K, 1.51).

Posterior Exploration
The main motivation of our architecture is to allow accurate unbiased posterior sampling of our hydrodynamical quantity through Gaussian sampling of our latent space variable. In this subsection we explore how accurate this mapping is.
This question is not straightforward as it is computationally difficult to calculate the true posterior for the target mapping, which would essentially require running millions of additional hydrodynamical simulations. With this in mind we are left to check if the samples drawn have the correct statistical properties of a posteriori sample, as well as examine its qualitative behaviors.
We draw 1000 random samples from a unit Gaussian distribution in latent space and predict the hydrodynamical quantities associated with a test dark matter field. To test whether or not our variance estimates are accurate, we plot the χ values for our test images in Figure 7. For a true posterior we would expect this distribution to be approximately Gaussian.
In particular, we can use the samples to construct a covariance matrix to use to test the statistical significance of deviations away form the MAP solution using the standard chi-squared formula, where in the last arrow we assume a diagonal covariance for implementation/memory reasons. If each sample is truly independent and represents the posterior, we expect that the corresponding chi-squared values from the ensemble should be Gaussian distributed with zero bias and standard deviation of one. This is arguably a necessary, but not sufficient, condition for the samples to represent a maximum a posteriori sample. We show this distribution in Figure 7.
An additional property desired of the posterior is for the samples to capture the qualitative uncertainty of shocked regions when sampled over latent space. We show one example of such a region in Figure 8. We expect dense regions at nodes of the cosmic web to have the largest hydrodynamical effects, while those in underdense environments away from cosmic structures should follow roughly power-law distributions without uncer- Figure 7. Distributions of the χ value for each of our three reconstructions. Also plotted is a unit normal Gaussian distribution, showing close agreement. Excess values along the tails indicate model failure, either to estimate proper variance or strong residual biases, or potentially importance of off-diagonal covariance terms not captured in this analysis. We find these value points constitute a very small volume fraction ( 0.1%). The field with the most variance that is not captured by our model is the temperature field, which is also the most numerically complex to calculate in the original hydrodynamical simulation due to shock physics.
tainty. In the samples in Figure 8, one sees high variability of the specifics of the shocked region, indicating the model is learning to account for hydrodynamical uncertainty. Again, it is difficult to formulate a rigorous method to test whether this variability is the "true" uncertainty without running a large suite of hydrodynamical simulations.

Predicted Lyman Alpha Forest
The Lyα forest arises from the scattering of photons at the characteristic rest frame Lyα frequency along their path from a background source, generally either a quasar or galaxy, to the observer. The fraction of the transmitted flux is given by F = exp(−τ ), where τ is the optical depth of the intervening gas. The optical depth in redshift space at a given velocity coordinate u along the line of sight is given by where u is the component of the Hubble flow velocity field u along the line-of-sight, over which the integral is calculated and b(u ) = 2k B T (u )/m p is the thermal broadening of the absorption feature. With the output of the hyphy model, we have all the components necessary to calculate the predicted Lyα forest statistics. We run the gimlet (Friesen et al. 2016) library to numerically calculate Eq. 7 on both the original simulation and our predicted output, using a mapping from hydrogen fraction to neutral hydrogen from Rahmati et al. (2013). Flux distributions are calculated along each of three axes and then averaged out. In Figure 9, we show the resulting flux distribution and in Figure 10 we show the error on the reconstructed power-spectra.

CONCLUSION
In this work, we have provided a flexible generic mapping from dark matter fields to hydrodynamical quantities, in particular the gas temperature, density, and velocity. We demonstrated that this mapping provides accurate reconstructions across a wide range of scales and captures a number of the statistical properties of the underlying truth fields. In addition to constructing this mapping, the underlying model provides posterior samples which are consistent with the underlying dark matter field. While it is computationally infeasible to analytically calculate the true maximum likelihood posterior in this case, we show that our posterior samples are a valid posterior through their variance properties.
It is important to note that we were able to construct these mappings despite only training on a single nyx simulated run, using data at z = 2.4, z = 3.0, and z = 4.0. This was possible due to various data augmentations used which exploit the symmetry of the hydrodynamical physics as well as a regularizing effect of the underlying CVAE architecture (Kamyab et al. 2019). Despite this relatively small training volume, we are still able to capture the qualitative effects and associated uncertainties in shocked regions. We do find some cosmic environments where there are constant residual biases which, despite significant testing, we were unable to reduce completely. While we recover the statistical properties of the WHIM, WCGM and diffuse IGM very well, recovering the hot halo gas contribution is significantly more difficult. We find that our model systematically under-estimates the density in these regions, while overestimating its overall volume fraction, resulting in an incorrect estimate for the phase space distribution. Gas in this region accounts for a very small small percent of the total volume, ∼ 4 × 10 −6 %, likely resulting in hyphy having difficulty capturing its properties.
For Lyα Forest analysis, these hot CGM regions do not have a noticeable effect on the forests' statistical properties so we do not focus on optimizing this aspect. Going beyond this work, one could perform various data augmentations to increase its importance in the loss to result in better model performance in this region. Approaches to dealing with such unbalanced data are well studied in the machine learning literature (Wang et al. 2016), including oversampling the minority data (Khoshgoftaar et al. 2007) and using Generative Adversarial Networks to create synthetic minority data (Kiyoiti dos Santos Tanaka & Aranha 2019).
One of the most compelling applications of the hyphy network is in forward modelling, where the underlying density field is reconstructed from observed data through an optimization process. For this application, the hyphy network could replace analytical or semianalytical approximations, such as the fluctuating Gunn Peterson approximation for Lyman alpha forest (e.g. in TARDIS, Horowitz et al. (2019Horowitz et al. ( , 2021). Our model is fully differentiable, allowing propagation through a model using first or second order methods. This is similar in spirit to work done in Modi et al. (2018), where a neural network was used to paint halo fields onto forward modeling dark matter density. However, this assumed a deterministic mapping from dark matter to galaxy light, while in hyphy this mapping is controlled by a latent space. In this approach, hydrodynamical uncertainties could be marginalized out via sampling of the latent space during optimization, or jointly optimized for and then marginalized out via variational methods.
In our loss function, we are implicitly assuming our latent space posterior can be well approximated by the dimensionality of our multivariate Gaussian latent space. It is possible a more accurate posterior would be achievable with an adversarial loss function which doesn't rely on this assumption, but such an approach comes with additional cost of difficulty in training. In addition, there is the added possibility of biasing the cosmological results due to the adversarial function implicitly learning the cosmology of the training set. For example, the adversarial loss function could implicitly learn to calculate a power spectra-like function, which would force the generative network to always produce maps with the the same power spectra as the training samples. With an L1  norm we are hopeful that our network would be transferable to other cosmological models, as these models would effect only the dark matter field without appreciably altering the hydrodynamical mapping. We expect the overall design of the network to be easily extended to other hydrodynamical properties. Simulations suites such as CAMELS (Villaescusa-Navarro et al. 2020), which include parametric feedback models, would allow the creation of an all-purpose mapping tool from dark matter to hydrodynamical models conditioned on underlying physics and redshift. While in nyx the main source of effective stochastic processes are gas shocks, other hydrodynamical simulation tools like AREPO (Springel 2010;Weinberger et al. 2020), used in CAMELS, allow for feedback processes through star formation, active galactic nuclei, supernova, etc. When studying these phenomena, it becomes very important to include a stochastic component to the network as done in hyphy. We plan on further exploring these properties in future works.

Constant redshift (z = 4.0)
Light cone redshift Constant redshift (z = 2.4) Figure 11. Figure demonstrating the generation of a light cone using hyphy. In the top two panels, we show the hyphy generated temperature corresponding to the same dark matter distribution as though it at z = 2.4 and z = 4.0. In the bottom panel we show the same distribution but with a varying redshift across the box, going from z = 4.0 (left) to z = 2.4 (right). redshift index, smoothly interpolates between the two. We plan to further apply this lightcone technique for mock generation and forward model reconstructions in future works.