Suhas S. Jain

The reliability of machine learning in multiscale physical systems depends on how physical structure is embedded into the learning process. We investigate this in the context of turbulent multiphase flows, focusing on the prediction of subgrid interfacial area density, a key quantity governing interphase transport that remains unresolved in large-eddy simulations. In this work, we develop and evaluate two machine learning subgrid closure models to predict the three-dimensional subgrid interfacial area density: a purely data-driven 3D encoder-decoder network, and a physics-constrained variant regularized by a fractal geometric prior. Across a range of Weber numbers, the physics-based model improves predictive accuracy, reduces error variance, and suppresses nonphysical artifacts relative to purely data-driven approaches. We also show that these gains are regime-dependent: the embedded inductive bias enhances generalization in corrugation-dominated regimes where its underlying assumptions hold, but becomes ineffective in fragmentation-dominated regimes characterized by topology change and droplet breakup. These results reveal a broader principle for scientific machine learning: the utility of physics-informed models depends not only on the presence of inductive bias, but on its alignment with the governing physical regime. This suggests a path toward regime-aware learning frameworks for modeling of complex multiscale systems.

The reliability of machine learning in multiscale physical systems depends on how physical structure is embedded into the learning process. We investigate this in the context of turbulent multiphase flows, focusing on the prediction of subgrid interfacial area density, a key quantity governing interphase transport that remains unresolved in large-eddy simulations. In this work, we develop and evaluate two machine learning subgrid closure models to predict the three-dimensional subgrid interfacial area density: a purely data-driven 3D encoder-decoder network, and a physics-constrained variant regularized by a fractal geometric prior. Across a range of Weber numbers, the physics-based model improves predictive accuracy, reduces error variance, and suppresses nonphysical artifacts relative to purely data-driven approaches. We also show that these gains are regime-dependent: the embedded inductive bias enhances generalization in corrugation-dominated regimes where its underlying assumptions hold, but becomes ineffective in fragmentation-dominated regimes characterized by topology change and droplet breakup. These results reveal a broader principle for scientific machine learning: the utility of physics-informed models depends not only on the presence of inductive bias, but on its alignment with the governing physical regime. This suggests a path toward regime-aware learning frameworks for modeling of complex multiscale systems.

This study experimentally investigates bubble size evolution and void fraction redistribution in an unexplored, coalescence-dominated regime of a decaying turbulent bubbly flow. The flow is generated downstream of a regenerative pump in a duct, with bulk Reynolds number ( Re ) $\sim \mathcal{O}(10^5)$ , Taylor-scale Reynolds number ( Re $_\lambda$ ) $\sim \mathcal{O}(10^3)$ and void fraction ( $\phi$ ) $\sim \mathcal{O}(1\,\%)$ , where the inlet turbulence is extremely intense (turbulence intensity $\gt 30\,\%$ ) but decays rapidly along the duct. Shadowgraph imaging and particle shadow velocimetry are used for measurements. The experimentally obtained turbulent dissipation in the duct flow decays as $\varepsilon \sim \mathcal{L}^{-2}$ , where $\mathcal{L}$ is the axial position, in close agreement with the homogeneous isotropic turbulence prediction of $\varepsilon \sim \mathcal{L}^{-2.2}$ . High-speed imaging and statistical analysis reveal that bubble coalescence dominates over breakup across most of the domain, leading to monotonic growth in the Sauter mean diameter ( $d_{32}$ ) and progressive broadening of the bubble size distribution. The normalised extreme-to-mean diameter ratio ( $\mathcal{D}$ ) increases axially and asymptotically from ${\sim} 1.9$ (breakup regime) and saturates at ${\sim} 2.2$ (coalescence regime), indicating the emergence of a quasi-self-similar bubble size distribution. The probability density function of the bubble diameter exhibits a dual power-law tail with exponents $-10/3$ and $-3/2$ near the duct inlet. However, after a few hydraulic diameters, a single $-3/2$ power-law scaling emerges, indicating a regime of pure coalescence in which all bubbles are smaller than the Hinze scale. The cumulative distribution plotted against $d/d_{32}$ shows that the slope decreases and the distribution width increases with both axial position and void fraction $(\phi )$ . Although classical Hinze scaling gives $d_{\textit{H}} \propto \mathcal{L}^{0.9}$ , our theory for $d_{32

We present direct numerical simulations (DNS) of laminar to turbulent transition in Herschel-Bulkley (HB) yield-stress fluids flowing through pipes and rectangular channels. The simulations employ a Herschel-Bulkley formulation that captures the yield-stress-driven plug, its breakdown, and the emergence of near-wall turbulent structures, enabling direct resolution of the transition mechanisms. The DNS cover a broad range of generalized Reynolds numbers, Re_G = 378 to 5300, allowing us to resolve plug formation, transition onset, and fully turbulent regimes. In pipe flow, the simulations reproduce the characteristic transition sequence, which includes a strong plug and negligible turbulence at low Re_G, a sharp rise in turbulence intensity and u'rms within a narrow transitional window (Re_G ~ 2000 to 3000), and wall-dominated turbulence with a weakened core at higher Re_G. Transition occurs only when local Reynolds stresses exceed the yield stress. The resulting regime boundaries (Re_G < 1735 laminar, 1735 < Re_G < 2920 transitional, and Re_G > 2920 turbulent) align with trends reported for Carbopol fluids. This work provides the first DNS resolving the complete laminar to turbulent transition in HB fluids for both pipe and channel configurations, offering unified insight into plug breakdown, turbulence localization, and the role of yield stress in transition mechanisms. Experimental validation using a 3.6 m acrylic channel with particle image velocimetry (PIV) is planned to further assess the DNS predictions and quantify geometry-dependent transition thresholds.

We present direct numerical simulations (DNS) of laminar to turbulent transition in Herschel-Bulkley (HB) yield-stress fluids flowing through pipes and rectangular channels. The simulations employ a Herschel-Bulkley formulation that captures the yield-stress-driven plug, its breakdown, and the emergence of near-wall turbulent structures, enabling direct resolution of the transition mechanisms. The DNS cover a broad range of generalized Reynolds numbers, Re_G = 378 to 5300, allowing us to resolve plug formation, transition onset, and fully turbulent regimes. In pipe flow, the simulations reproduce the characteristic transition sequence, which includes a strong plug and negligible turbulence at low Re_G, a sharp rise in turbulence intensity and u'rms within a narrow transitional window (Re_G ~ 2000 to 3000), and wall-dominated turbulence with a weakened core at higher Re_G. Transition occurs only when local Reynolds stresses exceed the yield stress. The resulting regime boundaries (Re_G < 1735 laminar, 1735 < Re_G < 2920 transitional, and Re_G > 2920 turbulent) align with trends reported for Carbopol fluids. This work provides the first DNS resolving the complete laminar to turbulent transition in HB fluids for both pipe and channel configurations, offering unified insight into plug breakdown, turbulence localization, and the role of yield stress in transition mechanisms. Experimental validation using a 3.6 m acrylic channel with particle image velocimetry (PIV) is planned to further assess the DNS predictions and quantify geometry-dependent transition thresholds.

Classical diffuse-interface methods for non-equilibrium formulations of two-phase flows are not capable of preserving interface thickness. Furthermore, when coupled with shock-capturing schemes, the interface may diffuse in its entirety. Thus, it is of utmost importance for accurate simulation of two-phase flows to develop interface-preserving formulations. However, standard phase field models are observed to remove any initial boundedness of small amounts of other phase added for robustness of diffuse interface methods. Therefore, the previously proposed conservative interface-regularization formulation based on the phase-field method is modified to guarantee that the initial volume fraction bound is retained throughout the simulations. The formulation is shown to retain the interface thickness and the small amount of the other phase added even when discretized with low-cost and low-dissipation kinetic energy–and entropy preserving schemes. Surface tension is coupled with the compressible two-phase formulation and validated against a standard test case. Lastly, robustness of the formulation is also assessed for highly compressible cases, which uses a flux blending strategy with a Lax-Friedrichs flux to achieve a positivity-preserving and entropy stable formulation. In highly compressible scenarios, such as shock-water column interaction, the framework is demonstrated to be robust and qualitative agreement is observed.

This work presents a high-fidelity, multi-physics framework for investigating the heat transfer of aircraft icing under realistic conditions. The methodology integrates a compressible Wall-Modeled Large-Eddy Simulation (WMLES) solver for aerodynamics, a Lagrangian particle tracking method for droplet impingement, and an unsteady Conjugate Heat Transfer (CHT) solver for the solid ice domain. This framework is applied to a NACA 23012 airfoil featuring laser-scanned rime ice roughness to analyze the interplay between surface topology, droplet impingement, and thermal conduction driven by the solidification of supercooled droplets. Results reveal that surface roughness induces high-frequency spatial variations in the local collection efficiency, leading to a heterogeneous heat source distribution from droplet solidification. However, lateral conduction within the ice layer diffuses this localized heat generation, resulting in a surface temperature field that is notably smoother than the impingement profile. Analysis of the interfacial energy budget at the stagnation region indicates that approximately 30\% of the latent heat released by droplet freezing is conducted into the ice layer. Temperature profiles within the solid reveal that near the inner airfoil skin, the conduction is mostly one-dimensional. However, near the highly irregular ice surface, the three-dimensional nature of the roughness elements drives multi-dimensional conduction, redirecting heat laterally from warm windward faces to cool leeward regions of individual roughness elements. These findings suggest that traditional icing models relying on isothermal wall assumptions or simplified one-dimensional conduction may not adequately capture the local thermodynamics of rough ice accretions.

Two-phase turbulent flows are ubiquitous with applications in aerospace, energy, and process industries. The interfacial area is an important physical parameter in such flows, since it is used to calculate mass, momentum, and energy transfer between the phases. High-fidelity simulation approaches, like direct numerical simulations (DNS), provide an accurate estimate of the interfacial area, albeit at very high and often unfeasible computational costs. Lower-and mid-fidelity approaches estimate the closure term, namely the subgrid interfacial area, using highly simplified empirical models. In this work, we present a physics-aided machine learning model that predicts the subgrid interfacial area, which can be used as a closure model for the simulation of two-phase turbulent flows.

This work presents high-fidelity numerical simulations of incompressible two-phase mixing layers with a focus on interfacial dynamics and turbulence. Two sets of simulations are conducted with two-dimensional (2D) and three-dimensional (3D) configurations using a conservative phase-field method. Grid refinement studies are carried out from N = 256 to N = 2048 grid points in both configurations to study convergence. Growth of momentum and volume-fraction-based thicknesses is examined, revealing initial linear behavior across all resolutions. However, interface breakup is increasingly delayed at higher resolutions in two-dimensional setup with the finest grid displaying unexpectedly similar behavior to coarser grids. This suggests that physical breakup mechanisms cannot be fully captured in 2D due to the absence of 3D instabilities. Visualizations further reveal significant resolution-dependent ligament breakup. A grid convergence study in 3D (with a significantly longer domain in the third dimension) was performed, obtaining a converged growth rate with increasing resolution.

In this work, we perform numerical simulations of forced two-phase isotropic turbulence to study the stationary states of a two-phase mixture. We first formulate three different approaches to force a two-phase turbulent flow that maintains a constant (a) mixture kinetic energy, (b) mixture kinetic energy + surface energy, and (c) individual phase kinetic energy, and study their effect on the final stationary state of the system. We show that forcing the two-phase system eliminates the arbitrariness associated with the initialization of the second phase in two-phase turbulence simulations. Next, we study the global statistics of the two-phase mixture in the stationary state for varying density ratios ( 1 0 − 3 to 1 0 3 ) and viscosity ratios ( 1 0 − 3 to 1 0 3 ), void fraction (dilute to dense regimes), and Weber numbers (breakup and non-breakup regimes). We utilize a spatial filtering approach, that is applicable for the general case of incompressible/compressible two-phase flows, to compute the energy spectra of each of the phases and study the turbulent kinetic energy distribution across scales. We find that the interface behaves as a soft wall, inhibiting the accommodation of eddies inside the dispersed phase. Finally, we show that bubbles and drops exhibit different turbulent characteristics for the same global conditions of the carrier phase, with bubbles exhibiting higher turbulent intensities inside them and droplets behaving ballistically, similar to solid particles. This study acts as a foundational work on stationary two-phase turbulence, which can be used for further analysis and the development of subgrid-scale models for larger-scale simulations of two-phase turbulent flows.

In this paper, we address the challenge of accurately calculating droplet/bubble properties (e.g., volume, number) in diffuse-interface two-phase flow simulations. Currently, flood-fill algorithms can truncate a significant portion of the volume of droplets/bubbles contained within the diffuse interface region or artificially merge multiple droplets/bubbles. This error is also dependent on the volume fraction cutoff value, which is typically chosen to be 0.5 arbitrarily, in the flood-fill algorithms. We propose a simple volume-correction approach that incorporates an analytical approximation of the truncated volume to correct for the missing droplet/bubble volumes. This proposed method results in accurately recovering the dispersed phase volumes with minimal volume error over a wide range of volume fraction cutoff values, and hence, can also accurately recover the number of droplets/bubbles. This can be a valuable tool for accurate calculation of drop/bubble size distributions for analysis and for Eulerian-to-Lagrangian conversion of the dispersed phase in multi-scale modeling approaches.

Large-eddy simulations (LES) are used to investigate the aerodynamic performance of a semispan swept-wing model based on the NASA Common Research Model. Two ice shapes are considered: 1) a real ice shape from laser-scanning accreted ice in the NASA Icing Research Tunnel and 2) an artificial ice shape derived from the maximum combined cross-section of the 3D ice accretion. Simulations are performed in free air at Reynolds numbers, [Formula: see text], of [Formula: see text], [Formula: see text], and [Formula: see text], with a Mach number of 0.18. Integrated loads, pressure distributions, and wall-shear stress visualizations are compared to experiments over various angles of attack. Grid resolution studies show good agreement for the real ice geometry, with lift coefficient differences within three counts ([Formula: see text]). The artificial ice geometry, however, shows greater sensitivity to grid resolution. Despite good agreement in lift for the artificial ice shape, discrepancies in drag and moments are observed. Pressure profiles suggest error cancellation, explaining the lift agreement. The smooth surfaces in the artificial ice shape introduce laminar boundary layers, laminar separation, and turbulent transition, unlike real ice, where roughness triggers immediate transition. These findings highlight the importance of incorporating appropriate roughness scales into artificial or predicted ice shapes for practical LES in iced aircraft configurations.

In this paper, we address the challenge of accurately calculating droplet/bubble properties (e.g., volume, number) in diffuse-interface two-phase flow simulations. Currently, flood-fill algorithms can truncate a significant portion of the volume of droplets/bubbles contained within the diffuse interface region or artificially merge multiple droplets/bubbles. This error is also dependent on the volume fraction cutoff value, which is typically chosen to be 0.5 arbitrarily, in the flood-fill algorithms. We propose a simple volume-correction approach that incorporates an analytical approximation of the truncated volume to correct for the missing droplet/bubble volumes. This proposed method results in accurately recovering the dispersed phase volumes with minimal volume error over a wide range of volume fraction cutoff values; and hence, can also accurately recover the number of droplets/bubbles. This can be a valuable tool for accurate calculation of drop/bubble size distributions for analysis and for Eulerian-to-Lagrangian conversion of the dispersed phase in multi-scale modeling approaches.

In this work, we perform numerical simulations of forced two-phase isotropic turbulence to study the stationary states of a two-phase mixture. We first formulate three different approaches to force a two-phase turbulent flow that maintains a constant (a) mixture kinetic energy, (b) mixture kinetic energy + surface energy, and (c) individual phase kinetic energy, and study their effect on the final stationary state of the system. We show that forcing the two-phase system eliminates the arbitrariness associated with the initialization of the second phase in two-phase turbulence simulations. Next, we study the global statistics of the two-phase mixture in the stationary state for varying density ratios ($10^{-3}$ to $10^3$) and viscosity ratios ($10^{-3}$ to $10^3$), void fraction (dilute to dense regimes), and Weber numbers (breakup and non-breakup regimes). We utilize a spatial filtering approach, that is applicable for the general case of incompressible/compressible two-phase flows, to compute the energy spectra of each of the phases and study the turbulent kinetic energy distribution across scales. We find that the interface behaves as a soft wall, inhibiting the accommodation of eddies inside the dispersed phase. Finally, we show that bubbles and drops exhibit different turbulent characteristics for the same global conditions of the carrier phase, with bubbles exhibiting higher turbulent intensities inside them and droplets behaving ballistically, similar to solid particles. This study acts as a foundational work on stationary two-phase turbulence, which can be used for further analysis and the development of subgrid-scale models for larger-scale simulations of two-phase turbulent flows.

Classical diffuse interface methods for non-equilibrium formulations of two-phase flows are not capable of preserving interface-thickness. Furthermore, when coupled with shock-capturing schemes, the interface may diffuse in its entirety. Thus, it is of utmost importance for accurate simulation of two-phase flows to describe interface-preserving formulations. Herein, a conservative interface-regularization formulation based on the phase-field method is proposed, and it analytically guarantees boundedness of volume fraction, absence of spurious oscillations of phasic kinetic energy and mixture entropy preservation in the absence of relaxation. The formulation is shown to stabilize dispersion errors originated due to the interface jumps and to preserve interface thickness when discretized with low lost and low dissipation second-order schemes. A possible issue with this formulation and discretization is the possibility of failing to preserve positivity of quantities of interests, hence, a blending strategy with a Lax-Friedrichs flux is also proposed and tested. In highly compressible scenarios such as shock-drop interaction, the interface regularization when applied to the positive-preserving flux is demonstrated to improve qualitative results on coarser grids.

Accurate modeling of ice accretion is important for the safe and efficient design of aircraft and propulsion systems. The first step in calculating the volumes of ice accumulated on aircraft wings is to estimate the impingement rate of droplets on the lifting surfaces. When these droplets exceed diameters of 40 ��m, the splashing effects become relevant in the calculation of the collection efficiency. In this work, we present large-eddy simulations (LES) of supercooled large droplets (SLD) impingement rates using a Lagrangian particle approach. This technique is used to represent droplet clouds with and without the presence of SLDs. Best practices to accelerate the convergence rate of the collection efficiency using LES and lagrangian particles are presented. It is shown that simulating droplet clouds with significant variations in droplet sizes leads to a very slow time convergence of the collection efficiency distribution. Calculating the collection efficiency using a combination of monodisperse cloud simulations leads to a 7-fold reduction in computation time as opposed to using a single polydisperse cloud simulation. Splashing models are evaluated to assess their accuracy within the LES framework. The results presented in this article demonstrate that LES coupled with Lagrangian particle tracking is highly effective in accurately modeling the impingement rate of droplets on aerodynamic surfaces.

For a turbulent rough-wall flow, it is known that the outer layer of a turbulent boundary layer is independent of wall roughness effects except in the roughness's role in setting both the friction velocity and boundary layer thickness (Jimenez, 2004; Kadivar et al., 2021). To accurately model these effects, we leverage a recently developed parameterized velocity transformation for rough-wall boundary layers (Bornhoft et al., 2023). This transformation is utilized to develop a new wall model for marginally resolved rough surfaces, modifying the equilibrium wall model (EQWM). The new wall model is validated in a turbulent channel flow and compared against available DNS data, as well as on an early-time rime ice geometry against both experimental and wall-resolved LES data. In both cases, significant improvements are observed in local quantities such as velocity profiles, boundary layer growth, and local friction coefficients, as well as in integrated quantities such as lift and drag.

A localized artificial-diffusivity method is developed for capturing discontinuities, such as shocks and contacts, in compressible flows. A new sensor for contact discontinuity makes the method more localized, and a discretely consistent formulation eliminates the need for filtering the solution or filtering the sensors to obtain robust solutions. Improved predictions are observed in canonical shock-tube problems and large-eddy simulations of homogeneous isotropic turbulence.

In this study, we employ wall-modeled large-eddy simulations (WMLES) to analyze the aerodynamic performance of an 8.9%-scaled semi-span swept wing model based on the NASA Common Research Model (CRM). The simulations consider two ice shapes: (1) a real-ice shape obtained from laser-scanning accreted ice in the NASA Icing Research Tunnel (IRT), and (2) a processed artificial-ice shape derived by constructing the maximum combined cross-section of the 3D ice accretion. Both ice shapes are simulated in a free-air configuration at a mean aerodynamic chord-based Reynolds number, Re_MAC, of 1.6x10^6 and Mach number of 0.18. Integrated loads (lift, drag, and moment coefficients), pressure distributions, and wall-shear stress visualizations are compared against experimental measurements for several angles of attack in both pre- and post-stall regions. The simulations are conducted at three different grid resolutions: 24 million control volumes (CV), 85 million CV, and 320 million CV. We observe good agreement across all cases for the real-ice geometry for both integrated loads and pressure distributions. The lift coefficient results closely match the experimental data, with differences within three lift counts (∆C_L=0.03) across all grid resolutions. However, the artificial-ice geometry displays greater sensitivity to grid resolution. Specifically, while the lift coefficient aligns well with the artificial-ice shape, disparities with the experimental data are observed in drag and moments. A closer examination through pressure profiles reveals that the apparent agreement in lift results from error cancellation. Introducing smooth surfaces in the artificial ice shape introduces physics not present in real ice flows, such as laminar boundary layers, laminar separation, and turbulent transition. These issues do not exist for the real ice geometry, as the roughness on the leading edge of the ice triggers bypass transition. For the artificial ice geometry, we run an additional grid resolution where we model the smooth leading edge ice shape with a no-slip condition, resulting in a grid with 1.5 billion CVs. This modification yields good agreement when comparing integrated loads and pressure profiles. This work underscores the challenges of modeling artificial ice shapes using state-of-the-art simulation approaches. We emphasize the importance of incorporating appropriate roughness scales into artificial or predicted ice shapes to ensure practical computational costs for the application of WMLES in iced aircraft configurations.

Cluster and void formations are key processes in the dynamics of particle-laden turbulence. In this work, we assess the performance of various neural network models for synthesizing preferential concentration fields of particles in turbulence. A database of direct numerical simulations of homogeneous isotropic two-dimensional turbulence with one-way coupled inertial point particles, is used to train the models using vorticity as the input to predict the particle number density fields. We compare autoencoder, U--Net, generative adversarial network (GAN), and diffusion model approaches, and assess the statistical properties of the generated particle number density fields. We find that the GANs are superior in predicting clusters and voids, and therefore result in the best performance. Additionally, we explore a concept of ``supersampling", where neural networks can be used to predict full particle data using only the information of few particles, which yields promising perspectives for reducing the computational cost of expensive DNS computations by avoiding the tracking of millions of particles. We also explore the inverse problem of synthesizing the enstrophy fields using the particle number density distribution as the input at different Stokes numbers. Hence, our study also indicates the potential use of neural networks to predict turbulent flow statistics using experimental measurements of inertial particles.

In this work, we propose a novel scalar-transport model for the simulation of scalar quantities that are confined to the interface in two-phase flows. In a two-phase flow, the scalar quantities, such as salts and surfactants, can reside at the interface and can modify the properties of the interface, in the time scales of interest. This confinement of the scalars leads to the formation of sharp gradients of the scalar concentration values at the interface, presenting a serious challenge for its numerical simulations. To overcome this challenge, we propose a computational model for the transport of scalars that maintains the confinement condition for these quantities. The model is discretized using a central-difference scheme, which leads to a non-dissipative implementation that is crucial for the simulation of turbulent flows. The model is used with the ACDI diffuse-interface method (Jain, J. Comput. Phys., 2022), but can also be used with other algebraic-based interface-capturing methods. Furthermore, the provable strengths of the proposed model are: (a) the model maintains the positivity property of the scalar concentration field, a physical realizability requirement for the simulation of scalars, when the proposed criterion is satisfied, (b) the proposed model is such that the transport of the scalar concentration field is consistent with the transport of the volume fraction field, which results in effective discrete confinement of the scalar at the interface; and therefore, prevents the artificial numerical diffusion of the scalar into the bulk region of the two phases. Finally, we present numerical simulations using the proposed model for both one-dimensional and multidimensional cases and assess: the accuracy and robustness of the model, the validity of the positivity property of the scalar concentration field, and the confinement of the scalar at the interface.

A phase-field method for unstructured grids that is accurate, conservative, and robust is proposed in this work. The proposed method also results in bounded transport of volume fraction, and the interface thickness adapts automatically to local grid size. In addition to this, we present a novel formulation for two-phase flows on collocated grids that is provably energy stable, which is a critical feature for robust simulations of two-phase turbulent flows. The proposed method is first evaluated on canonical test cases, including scenarios like drop advection and drop in a shear flow. The accuracy of the proposed method in the presence of grid transitions, which is important for simulations in complex geometries using unstructured grids, is also evaluated. We assess the robustness of our method by performing simulations of a drop in homogeneous isotropic turbulence at infinite Reynolds number with varying density ratios. Furthermore, as a step toward verification and validation, simulations of test cases in complex geometries and conditions, such as, a damped surface wave, infinitesimal interface distortion on a liquid jet, and the atomization of a liquid jet from the engine combustion network's Spray A nozzle are presented. These simulations illustrate the accuracy, robustness, and applicability of the proposed method in various kinds of complex two-phase flow applications of engineering interest.

In this study, five ice shapes generated at NASA Glenn's Icing Research Tunnel (IRT) are simulated at multiple angles of attack (Broeren et al., J. of Aircraft, 2018). These geometries target different icing environments, both early-time and longer-duration glaze and rime ice exposure events, including a geometry that results from using a thermal ice-protection system. Using the laser-scanned geometries, detailed representations of the three-dimensional ice geometries are resolved on the grid and simulated using wall-modeled LES. Integrated loads (lift, drag, and moment coefficients) and pressure distributions are compared against experimental measurements in both clean and iced conditions for several angles of attack in both pre-and post-stall regions. The relevant comparisons to the experimental results show that qualitative and acceptable quantitative agreement with the data is observed across all geometries. Glaze ice formations exhibit larger and highly nonuniform ice features, such as `horns', in contrast to rime ice formations characterized by smaller, uniformly distributed roughness elements. In wall-modeled LES, it was observed that larger roughness scales in the glaze ice that trigger transition can be accurately resolved. Therefore, it is possible for WMLES to accurately capture the aerodynamics of glaze ice shapes without the need for additional modeling. In contrast, rime ice geometries required additional resolution to accurately represent the aerodynamic loads. This study demonstrates the effectiveness of the wall-modeled LES technique in simulating the complex aerodynamic effects of iced airfoils, providing valuable insights for aircraft design in icing environments and highlighting the importance of accurately representing ice geometries and roughness scales in simulations.

International audience

In this work, a localized artificial-viscosity/diffusivity method is proposed for accurately capturing discontinuities in compressible flows. There have been numerous efforts to improve the artificial diffusivity formulation in the last two decades, through appropriate localization of the artificial bulk viscosity for capturing shocks. However, for capturing contact discontinuities, either a density or internal energy variable is used as a detector. An issue with this sensor is that it not only detects contact discontinuities, but also falsely detects the regions of shocks and vortical motions. Using this detector to add artificial mass/thermal diffusivity for capturing contact discontinuities is hence unnecessarily dissipative. To overcome this issue, we propose a sensor similar to the Ducros sensor (for shocks) to detect contact discontinuities, and further localize artificial mass/thermal diffusivity for capturing contact discontinuities. The proposed method contains coefficients that are less sensitive to the choice of the flow problem. This is achieved by improved localization of the artificial diffusivity in the present method. A discretely consistent dissipative flux formulation is presented and is coupled with a robust low-dissipative scheme, which eliminates the need for filtering the solution variables. The proposed method also does not require filtering for the discontinuity detector/sensor functions, which is typically done to smear out the artificial fluid properties and obtain stable solutions. Hence, the challenges associated with extending the filtering procedure for unstructured grids is eliminated, thereby, making the proposed method easily applicable for unstructured grids. Finally, a straightforward extension of the proposed method to two-phase flows is also presented.

&lt;div class="section abstract"&gt;&lt;div class="htmlview paragraph"&gt;In this work, ice accretion is investigated on a fundamental level using a novel Eulerian phase field approach that captures the phase interface. This method, unlike the Allen-Cahn method, does not lead to spurious phase change (artificial mass loss). This method is also straightforward to implement and avoids normal vector reconstructions along the interface or ghost cells. Additionally, it has well-defined and novel stiffness constraints for accuracy and stability that define parameters in the model such as the kinetic coefficient &lt;i&gt;μ&lt;/i&gt; and the interface regularization coefficient &lt;i&gt;γ&lt;/i&gt;. An incompressible solver is constructed and used to verify the new method using an analytical Stefan problem solution in both 1D and 2D domains.&lt;/div&gt;&lt;/div&gt;

&lt;div class="section abstract"&gt;&lt;div class="htmlview paragraph"&gt;Modeling of icing is important for the design of aircraft lifting surfaces and for the design of efficient propulsion systems. The computational modeling of ice accretion prediction is important to replace the expensive experimental techniques for calculating the ice shapes in Icing tunnels, and the first step toward modeling ice accretion is to accurately compute the droplet collection efficiency which acts as the input to the accretion model. In this work, we perform large-eddy simulations of supercooled droplet transport and impingement onto complex aircraft geometries using a Lagrangian particle approach. We assess the improvement in modeling droplet impingement by computing the droplet collection efficiency and by comparing with the existing experimental data.&lt;/div&gt;&lt;/div&gt;

&lt;div class="section abstract"&gt;&lt;div class="htmlview paragraph"&gt;Predicting the aerodynamic performance of an aircraft in icing conditions is critical as failures in an aircraft’s ice protection system can compromise flight safety. Aerodynamic effects of icing have typically relied on RANS modeling, which usually struggles to predict stall behavior, including those induced by surface roughness. Encouraged by recent studies using LES that demonstrate the ability to predict stall characteristics on full aircraft with smooth wings at an affordable cost [&lt;span class="xref"&gt;1&lt;/span&gt;], this study seeks to apply this methodology to icing conditions. Measurements of lift, drag, and pitching moments of a NACA23012 airfoil under clean and iced conditions are collected at Re = 1.8M. Using laser scanned, detailed representations of the icing geometries, LES calculations are conducted to compare integrated loads against experimental measurements in both clean and iced conditions at various angles of attack through the onset of stall [&lt;span class="xref"&gt;2&lt;/span&gt;]. This study will explore several critical ice shapes to validate our approach. These include early-time rime, early-time glaze, and horn ice shapes.&lt;/div&gt;&lt;/div&gt;

The predictive capability of two way--coupled point-particle Euler-Lagrange model in accurately capturing particle-flow interactions under grid refinement, wherein the particle size can be comparable to the grid size, is systematically evaluated. Two situations are considered, (i) uniform flow over a stationary particle, and (ii) decaying isotropic turbulence laden with Kolmogorov-scale particles. Particle-fluid interactions are modeled using only the standard drag law, typical of large density-ratio systems. A zonal, advection-diffusion-reaction (Zonal-ADR) model is used to obtain the undisturbed fluid velocity needed in the drag closure. Two main types of interpolation kernels, grid-based and particle size--based, are employed. The effect of interpolation kernels on capturing the particle-fluid interactions, kinetic energy, dissipation rate, and particle acceleration statistics are evaluated in detail. It is shown that the interpolation kernels whose width scales with the particle size perform significantly better under grid refinement than kernels whose width scales with the grid size. Convergence with respect to spatial resolution is obtained with the particle size--based kernels with and without correcting for the self-disturbance effect. While the use of particle size--based interpolation kernels provide spatial convergence and perform better than kernels that scale based on grid size, small differences can still be seen in the converged results with and without correcting for the particle self-disturbance. Such differences indicate the need for self-disturbance correction to obtain the best results, especially when the particles are larger than the grid size.

Interfacial Rayleigh–Taylor mixing is crucial to describing important natural and engineering processes, such as exploding supernovae, laser micromachining, hot spots in inertial confinement fusion, and optical telecommunications. These require the characterization of the time dependence of the driving acceleration. We compare our theoretical formulation based on group theory foundations with interface-capturing numerical simulations for linear and nonlinear two-dimensional Rayleigh–Taylor instabilities in a finite-sized domain with time-varying acceleration over broad ranges of Atwood numbers and acceleration exponents. Detailed corroboration between theory and simulations is provided for this foundational case. Both demonstrate the strong interfacial nature of Rayleigh–Taylor instabilities, which suggests that practical flow fields can be reconstructed from the derived fluid potential using the proposed theory. A robust agreement is also obtained for the early and late-time evolution of the amplitudes of the bubble and spike, which demonstrate that the Rayleigh–Taylor flow can transition to the mixing regime even for a single-mode initial perturbation. Corroboration with experiments of high energy density plasmas motivated by studies of supernovae is also achieved. In addition, a long-standing puzzle in Rayleigh–Taylor dynamics on the interplay between the acceleration, the shear, and the interface morphology in the theory and simulations is resolved by accounting for finite viscosity of the fluids. The characterization of Rayleigh–Taylor instabilities as a highly interfacial phenomenon provides valuable insight into its multiscale nature, which enhances the design and understanding of numerous processes of practical interest.

Multiphase flows have a wide range of applications in natural and engineering processes. In this talk, I’ll present a novel diffuse-interface model and robust numerical methods for simulations of compressible multiphase flows. I’ll first present the accurate conservative diffuse-interface/phase-field (ACDI) model for the simulation of multiphase flows. This method conserves the mass of each of the phases, and results in bounded transport of the volume fraction while maintaining the interface thickness on the order of only one to two grid points. I’ll present results from the canonical test cases, showing the improvement in the accuracy of interface shape and surface tension forces over the commonly used second-order conservative phase-field method. The capability of the ACDI model to maintain such sharp interfaces without the need for any special geometric treatment, unlike the sharp-interface methods, makes it a highly attractive interface-capturing method for accurate simulation of multiphase flows at an affordable cost. Next, for the simulation of compressible multiphase flows, a five-equation model that consists of transport equations for the volume fraction, the mass of each phase, momentum, and total energy is used. Starting from this baseline five-equation model, I’ll present modifications to the model in such a way that the resulting system of equations can be discretized using a non-dissipative central scheme that is suitable for the simulation of turbulent flows. The resulting model is conservative, accurate, scalable, and maintains a constant interface thickness throughout the simulation. I will present simulations of canonical and complex turbulent compressible multiphase flows. For stable and accurate numerical simulations of compressible flows, particularly at higher Reynolds numbers (Re), it is known that a discrete entropy condition needs to be satisfied in addition to the discrete conservation of kinetic energy. I’ll present a numerical flux formulation for the five-equation model that satisfies this condition (a KEEP scheme) and show that this formulation results in stable numerical simulations of compressible turbulent multiphase flows at high Re. Finally, I’ll briefly highlight some of the related research efforts on the use of these methods for the simulation of shock-interface interaction problems, phase change, and multi-material systems.