Rayleigh-Taylor mixing is a classical hydrodynamic instability that occurs when a light fluid pushes against a heavy fluid. The two main sources of nonideal behavior in Rayleigh-Taylor (RT) mixing are regularizations (physical and numerical), which produce deviations from a pure Euler equation, scale invariant formulation, and nonideal (i.e., experimental) initial conditions. The Kolmogorov theory of turbulence predicts stirring at all length scales for the Euler fluid equations without regularization. We interpret mathematical theories of existence and nonuniqueness in this context, and we provide numerical evidence for dependence of the RT mixing rate on nonideal regularizations; in other words, indeterminacy when modeled by Euler equations. Operationally, indeterminacy shows up as nonunique solutions for RT mixing, parametrized by Schmidt and Prandtl numbers, in the large Reynolds number (Euler equation) limit. Verification and validation evidence is presented for the large eddy simulation algorithm used here. Mesh convergence depends on breaking the nonuniqueness with explicit use of the laminar Schmidt and Prandtl numbers and their turbulent counterparts, defined in terms of subgrid scale models. The dependence of the mixing rate on the Schmidt and Prandtl numbers and other physical parameters will be illustrated. We demonstrate numerically the influence of initial conditions on the mixing rate. Both the dominant short wavelength initial conditions and long wavelength perturbations are observed to play a role. By examination of two classes of experiments, we observe the absence of a single universal explanation, with long and short wavelength initial conditions, and the various physical and numerical regularizations contributing in different proportions in these two different contexts.large eddy simulations | subgrid scale models | turbulence R ayleigh-Taylor (RT) instability is a classical hydrodynamic instability (1, 2) that occurs when a light fluid pushes against a heavy fluid. The density contrast is measured dimensionlessly by the Atwood number A ¼ ðρ 2 − ρ 1 Þ∕ðρ 2 þ ρ 1 Þ, with ρ i the density of fluid i. The mixing rate can be characterized by the dimensionless parameter α, defined through the equationwhere h is the penetration distance of the light fluid into the heavy and g is the acceleration. The mixing rapidly becomes turbulent, with Reynolds numbers Re ≈ 50;000 observed in typical experiments. For this reason, simulations based on compressible codes, with hyperbolic time step (CFL) restrictions, are unable to resolve the viscous length scales, and such codes are run in an under-resolved manner. Such simulations are called large eddy simulations (LES). LES require subgrid scale (SGS) models to describe the effect of the small (subgrid) scales on the large (grid-resolved) ones. Models for the Reynolds stress based on eddy viscosity are the most familiar of this type.Much of the work on RT mixing has been influenced by ideas of universality and self similarity. If α is to be a universal physical...