Large-scale Simulations Research Articles

Efficiency and scalability of fully-resolved fluid-particle simulations on heterogeneous CPU-GPU architectures

Current supercomputers often have a heterogeneous architecture using both conventional Central Processing Units (CPUs) and Graphics Processing Units (GPUs). At the same time, numerical simulation tasks frequently involve multiphysics scenarios whose components run on different hardware due to multiple reasons, e.g., architectural requirements, pragmatism, etc. This leads naturally to a software design where different simulation modules are mapped to different subsystems of the heterogeneous architecture. We present a detailed performance analysis for such a hybrid four-way coupled simulation of a fully resolved particle-laden flow. The Eulerian representation of the flow utilizes GPUs, while the Lagrangian model for the particles runs on conventional CPUs. Two characteristic model situations involving dense and dilute particle systems are used as benchmark scenarios. First, a roofline model is employed to predict the node level performance and to show that the lattice-Boltzmann-based Eulerian fluid simulation reaches very good performance on a single GPU. Furthermore, the GPU-GPU communication for a large-scale Eulerian flow simulation results in only moderate slowdowns. This is due to the efficiency of the CUDA-aware MPI communication, combined with the use of communication hiding techniques. On 1024 A100 GPUs, an overall parallel efficiency of up to 71% is achieved. While the flow simulation has good performance characteristics, the integration of the stiff Lagrangian particle system requires frequent CPU-CPU communications that can become a bottleneck, especially when simulating the dense particle system. Additionally, special attention is paid to the CPU-GPU communication overhead since this is essential for coupling the particles to the flow simulation. However, thanks to our problem-aware co-partitioning, the CPU-GPU communication overhead is found to be negligible. As a lesson learned from this development, four criteria are postulated that a hybrid implementation must meet for the efficient use of heterogeneous supercomputers.

Nonlinear mechanical properties of irregular architected materials

Abstract Architected materials have received increasing attention due to their exotic mechanical properties including ultra-high stiffness to weight ratio, strength, energy absorption and toughness. Typically, their mechanical properties and deformation behavior arise from the periodically tessellated unit cells. Although periodicity in conventional architected materials promises homogeneity and predictability in mechanical behaviors, it imposes a strong restriction on the design space of architected materials. Inspired by bio-materials, aperiodic and disordered designs significantly expand the design space and has been proven effective in controlling and optimizing linear elastic properties. Making a step further, here we focus on the nonlinear properties of irregular lattice materials under large deformation, including the stress-strain curve and specific energy absorption. Such materials are generated by a nature-inspired virtual growth program that assembles predefined geometric building blocks in a stochastic yet controllable manner. The nonlinear properties are analyzed through quasi-static compression experiments and large-scale numerical simulations. Based on the well-agreed experimental and numerical results, through the lens of machine learning techniques, the nonlinear properties show strong correlation with the appearance frequency of the building blocks and their local connectivity, regardless of the non-deterministic nature of the microstructures. A practical constitutive model is proposed for future developments such as generative design and engineering application. Our research offers valuable insights and serves as an inspiration for deeper exploration into the intricate structure-property relationships within materials with aperiodic and disordered microstructures.

Turbulent Flow Through Sluice Gate and Weir Using Smoothed Particle Hydrodynamics: Evaluation of Turbulence Models, Boundary Conditions, and 3D Effects

Understanding flow dynamics around hydraulic structures is essential for optimizing water management systems and predicting flow behavior in real-world applications. In this study, we simulate a 3D flow control system featuring a sluice gate and a weir, commonly used in hydraulic engineering. The focus is on accurately incorporating modified dynamic boundary conditions (mDBCs) and viscosity treatment to improve the simulation of complex, turbulent flows. We assess the performance of the Smoothed Particle Hydrodynamics (SPH) method in handling these challenging conditions. Especially when the boundary conditions and applicability to industry are two of the SPH method’s grand challenges. Simulations were conducted on a Graphics Processing Unit (GPU) using the DualSPHysics code. The results were compared to theoretical predictions and experimental data found in the literature. Key hydraulic characteristics, including 3D flow effects, hydraulic jump formation, and turbulent behavior, are examined. The combination of mDBCs with the Laminar plus sub-particle scale turbulence model achieved the correct simulation results. The findings demonstrate agreement between simulations, theoretical predictions, and experimental results. This work provides a reliable framework for analyzing turbulent flows in hydraulic structures and can be used as reference data or a prototype for larger-scale simulations in both research and engineering design, particularly in contexts requiring robust and precise flow control and/or environmental management.

PARALLEL PERFORMANCE OF THE EXPLICIT FDTD METHODS APPLIED ON 3D ACOUSTIC WAVE EQUATION

The 3D acoustic wave equation is a second-order linear hyperbolic partial differential equation. It is widely studied in acoustics, fluid dynamics, electromagnetism, and seismology within both science and engineering fields. Among the various numerical methods, the finite difference method is considered to be the simpler to understand and easy to implement. Based on the discretize schemes finite difference time domain method described by the explicit scheme in which spatial second order derivatives are evaluated at the previous time step. The solution of 3D acoustic wave equation may be required on the high-resolution mesh points consequently more computational resources are required. In this study parallel algorithm for the numerical solution of 3D acoustic wave equation is proposed and designed by using data parallel approach and message passing schemes. The algorithm is implemented on shared memory parallel system using MATLAB parallel computing mode. The parallel performance of the designed algorithm is analyzed on different mesh sizes and time steps. It is revealed that the proposed algorithm may reduce computational time up to 3 times as compared to sequential solution algorithm. The proposed parallel algorithm remains more efficient on workers while for the efficiency of the algorithm drops because of the high communication time. The results of the proposed research could be utilized in the large-scale simulations of 3D acoustic wave equations and may enable to simulate the acoustic wave pressure at more refined meshes on high performance computing systems.

Redshift dependence of FRB host dispersion measures across cosmic epochs

We constrain the redshift dependence of (rest frame) host galaxy dispersion measures of localized Fast Radio Bursts (FRBs) by assuming it to vary as a simple power law (∝ (1+z) α ). We simultaneously fit α as well as the host dispersion measure to the data of FRBs with known redshifts. We find that α between 0 to 1 is preferred depending upon our modelling choices. Current data can constrain |α| ≲ 2 at a 68 percent confidence interval. Such constraints have implications for our understanding of galaxy formation and can be used to inform galaxy and large scale simulations.

Vortex characteristics of a large-scale Ward-type tornado simulator at Central South University

Vortex characteristics of a large-scale Ward-type tornado simulator at Central South University

Machine-Learning Modeling of Elemental Ferroelectric Bismuth Monolayer.

The bismuth monolayer has recently been experimentally identified as a novel platform for the investigation of two-dimensional single-element ferroelectric system. Here, we model the potential energy surface of a bismuth monolayer by employing a message-passing neural network and achieve an error smaller than 1.2meV per atom. Empowered by the high accuracy and fast prediction of the machine learning model, we have embarked on in-depth and large-scale atomistic simulations. These explorations are tailored to understand the temperature-dependent phase transitions, with an emphasis on the difference between free-standing monolayers and those constrained by a substrate. Furthermore, with the large system used in the simulations, we are also able to observe ferroelectric domains within these systems and shed light on their intrinsic lattice thermal conductivity.

Learning large-scale industrial physics simulations

In an industrial group like Safran, numerical simulations of physical phenomena are integral to most design processes. At Safran’s corporate research center, we enhance these processes by developing fast and reliable surrogate models for various physics. We focus here on two technologies developed in recent years. The first is a physical reduced-order modeling method for non-linear structural mechanics and thermal analysis, used for calculating the lifespan of high-pressure turbine blades and performing heat analysis of high-pressure compressors. The second technology involves learning physics simulations with non-parameterized geometrical variability using classical machine learning tools, such as Gaussian process regression. Finally, we present our contributions to the open-source and open-data community.

Anisotropic interactions for continuum modeling of protein-membrane systems.

In this work, a model for anisotropic interactions between proteins and cellular membranes is proposed for large-scale continuum simulations. The framework of the model is based on dynamic density functional theory, which provides a formalism to describe the lipid densities within the membrane as continuum fields while still maintaining the fidelity of the underlying molecular interactions. Within this framework, we extend recent results to include the anisotropic effects of protein-lipid interactions. As applications, we consider two membrane proteins of biological interest: a RAS-RAF complex tethered to the membrane and a membrane embedded G protein-coupled receptor. A strong qualitative and quantitative agreement is found between the numerical results and the corresponding molecular dynamics simulations. Combining the scope of continuum level simulations with the details from molecular level particle simulations enables research into protein-membrane behaviors at a more biologically relevant scale, which crucially can also be accessed via experiment.

Analytical Model for Atomic Relaxation in Twisted Moiré Materials.

By virtue of being atomically thin, the electronic properties of heterostructures built from two-dimensional materials are strongly influenced by atomic relaxation. The atomic layers behave as flexible membranes rather than rigid crystals. Here we develop an analytical theory of lattice relaxation in twisted moiré materials. We obtain analytical results for the lattice displacements and corresponding pseudo gauge fields, as a function of twist angle. We benchmark our results for twisted bilayer graphene and twisted WSe_{2} bilayers using large-scale molecular dynamics simulations. Our single-parameter theory is valid in graphene bilayers for twist angles θ≳0.7°, and in twisted WSe_{2} for θ≳1.6°. We also investigate how relaxation alters the electronic structure in twisted bilayer graphene, providing a simple extension to the continuum model to account for lattice relaxation.

Distinct deformation mechanisms of silicate glasses under nanoindentation: The critical role of structure

We use large-scale molecular dynamics simulations to investigate the indentation response of three silica-based glasses with varying compositional complexities. Our primary goal is to clarify the roles of the typical network-modifying species, namely, sodium, and the secondary network-forming species, namely, boron, in influencing the mechanical behavior of the glasses under localized stress. The distinct mechanical responses of the glasses are linked to structural features such as bond strength, network connectivity, and atomic packing density. The enhanced nanoscale ductility of sodium silicate and sodium borosilicate glasses, compared to silica, is attributed to the structural flexibility induced by Na atoms, which depolymerize the network, and by B species in mixed coordination. We also find that shear flow, driven by network flexibility, is the dominant deformation mechanism in the sodium silicate and sodium borosilicate glasses, while densification dominates in silica due to its low packing density. The evolution of short-to-intermediate-range structures is responsible for the distinct deformation behaviors of the glasses. These results highlight the critical role of structure in determining the deformation mechanisms of silicate glasses under sharp contact loads, providing insights for improving the mechanical performance of these materials.

Enhanced Entanglement in the Measurement-Altered Quantum Ising Chain

Understanding the influence of measurements on the properties of many-body systems is a fundamental problem in quantum mechanics and for quantum technologies. This paper explores how a finite density of stochastic local measurement modifies a given state’s entanglement structure. Considering various measurement protocols, we explore the typical quantum correlations of their associated projected ensembles arising from the ground state of the quantum Ising model. Using large-scale numerical simulations, we demonstrate substantial differences among inequivalent measurement protocols. Surprisingly, we observe that forced on-site measurements can enhance both bipartite and multipartite entanglement. We present a phenomenological toy model and perturbative calculations to analytically support these results. Furthermore, we extend these considerations to the non-Hermitian Ising model, naturally arising in optically monitored systems, and we show that its qualitative entanglement features are not altered by a finite density of projective measurements. Overall, these results reveal a complex phenomenology where local quantum measurements do not simply disentangle degrees of freedom, but may actually strengthen the entanglement in the system.

Role of lattice distortion on spallation of CoCrCuFeNi high-entropy alloy

High-entropy alloys (HEAs), known for their high strength and enhanced ductility, have promising applications across various fields. Lattice distortion is a key factor in their strengthening, yet its role in dynamic fracture strength or spall strength remains unclear. This study employs large-scale nonequilibrium molecular dynamics simulations to investigate the dynamic responses of equiatomic CoCrCuFeNi HEA under shock velocities ranging from 0.6 to 1.45 km/s. By comparing the systems described using an average-atom interatomic potential, we uncover the role of lattice distortion. Our results reveal that spall strength exhibits complex behavior depending on the HEA's shock response. As shock velocity increases, the deformation mechanism transitions from elastic to dislocation and stacking fault (SF) dominated, eventually leading to a face-centered cubic to body-centered cubic phase transition. Lattice distortion significantly alters the active slip planes of dislocations and SFs, resulting in more SF intersections, while its effect on compression-induced phase transition is minor. During shock-induced spallation, residual defects after decompression significantly affect spall strength. Lattice distortion introduces additional stress and strain concentration sites, facilitating void formation and reducing spall strength. The temperature at the spall region is identified as a major factor governing spall strength variation under different shock velocities. Although lattice distortion can mitigate the softening effect of elevated temperature, it ultimately reduces spall strength, challenging the traditional views on its strengthening role. Moreover, the effects of lattice distortions on spall strength are quantified in terms of lattice misfit under varying loading strain rates and initial temperatures.

Hydrostatic mass bias for galaxy groups and clusters in the FLAMINGO simulations

ABSTRACT The masses of galaxy clusters are commonly measured from X-ray observations under the assumption of hydrostatic equilibrium (HSE). This technique is known to underestimate the true mass systematically. The fiducial FLAMINGO (Full-hydro Large-scale structure simulations with All-sky Mapping for the Interpretation of Next Generation Observations) cosmological hydrodynamical simulation predicts the median hydrostatic mass bias to increase from $b_\text{HSE} \equiv (M_\text{HSE,500c}-M_\text{500c})/M_\text{500c} \approx -0.1$ to −0.2 when the true mass increases from group to cluster mass scales. However, the bias is nearly independent of the hydrostatic mass. The scatter at fixed true mass is minimum for $M_\text{500c}\sim 10^{14}~\text{M}_\odot$, where $\sigma (b_\text{HSE})\approx 0.1$, but increases rapidly towards lower and higher masses. At a fixed true mass, the hydrostatic masses increase (decrease) with redshift on group (cluster) scales, and the scatter increases. The bias is insensitive to the choice of analytic functions assumed to represent the density and temperature profiles, but it is sensitive to the goodness of fit, with poorer fits corresponding to a stronger median bias and a larger scatter. The bias is also sensitive to the strength of stellar and active galactic nucleus feedback. Models predicting lower gas fractions yield more (less) biased masses for groups (clusters). The scatter in the bias at fixed true mass is due to differences in the pressure gradients rather than in the temperature at $R_\text{500c}$. The total kinetic energies within $r_\text{500c}$ in low- and high-mass clusters are sub- and supervirial, respectively, though all become subvirial when external pressure is accounted for. Analyses of the terms in the virial and Euler equations suggest that non-thermal motions, including rotation, account for most of the hydrostatic mass bias. However, we find that the mass bias estimated from X-ray luminosity weighted profiles strongly overestimates the deviations from HSE.

Bending Moiré in Twisted Bilayer Graphene.

Moiré potentials caused by lattice mismatches significantly alter electrons in two-dimensional materials, inspiring the discovery of numerous unique physical properties. While strain modulates the structure and symmetry of the moiré potential, serving as a tuning mechanism for interactions, the impact of out-of-plane deformation, e.g., bending, on the moiré superlattice remains unknown. Here, we performed large-scale molecular dynamics simulations to study the evolution of the moiré superlattice of twisted bilayer graphene under out-of-plane bending deformation. Our findings indicated that curvature-dependent bending caused both global and local lattice structure modifications in the moiré superlattice. We revealed a linear relationship between lattice displacement and bending curvature across varying initial twist angles along with precise regulation of local interlayer rotation. Additionally, the atomic potential energy landscape revealed that the localized atomic stacks underwent a whirlpool-like transformation, becoming a relaxed superlattice. This work opens up new opportunities for tailoring moiré superlattices by using out-of-plane bending engineering.

A framework for validating efficiency models of thermal separation columns with tray internals

Abstract Tray columns are globally used for distillation processes, and their performance must be optimized to reduce their cost and energy expenditures as well as carbon emissions. A prerequisite to achieve these targets demands a realistic account of the tray performance based on tray efficiency model application. The proof of concept of a Refined Residence Time Distribution (RRTD) model, recently proposed by the authors, is presented. The multifaceted challenges in acquiring hydrodynamic and performance data suitable for model validation are discussed in this work. In that regard, an unprecedented experimental campaign was performed in a representative large-scale air-water tray column simulator based on the application of multiplex flow profiler, new chemical system, and novel data processing schemes. Using these data and additional case studies, the capabilities of the RRTD model were demonstrated. This model successfully accounted for the impacts of varying local liquid flow characteristics and vapor flow maldistributions on the tray efficiency thereby advancing the state of the art of the tray efficiency models. This work particularly aimed at proposing a constructive framework to evaluate tray efficiency models in the future campaigns, and those campaigns will benefit from the learnings of the present work.

Large-scale dual AGN in large-scale cosmological hydrodynamical simulations

Abstract Detecting dual active galactic nuclei (DAGN) in observations and understanding theoretically which massive black holes (MBHs) compose them and in which galactic and large-scale environment they reside are becoming increasingly important questions as we enter the multi-messenger era of MBH astronomy. This paper presents the abundance and properties of DAGN produced in nine large-scale cosmological hydrodynamical simulations. We focus on DAGN powered by AGN with $L_{\rm bol}\geqslant 10^{43}\, \rm erg\, s^{-1}$ and belonging to distinct galaxies, i.e. pairs that can be characterised with current and near-future electromagnetic observations. We find that the number density of DAGN separated by a few to 30 proper kpc varies from 10−8 (or none) to $10^{-3} \, \rm comoving\, Mpc^{3}$ in the redshift range z = 0–7. At a given redshift, the densities of the DAGN numbers vary by up to two orders of magnitude from one simulation to another. However, for all simulations, the DAGN peak is in the range z = 1–3, right before the peak of cosmic star formation or cosmic AGN activity. The corresponding fractions of DAGN (with respect to the total number of AGN) range from 0 to 6 %. We find that simulations could produce too few DAGN at z = 0 (or merge pairs too quickly) compared to current observational constraints while being consistent with preliminary constraints at high redshift (z ∼ 3). Next-generation observatories (e.g., AXIS) will be of paramount importance to detect DAGN across cosmic times. We predict the detectability of DAGN with future X-ray telescopes and discuss DAGN as progenitors for future LISA gravitational wave detections.

From surface roughness to crater formation in a 2D multi-scale simulation of ultrashort pulse laser ablation

Surface roughness plays a critical role in ultrashort pulse laser ablation, particularly for industrial applications using burst mode operations, multi-pulse laser processing, and the generation of laser-induced periodic surface structures. Hence, we address the impact of surface roughness on the resulting laser ablation topography, comparing predictions from a simulation model to experimental results. We present a comprehensive multi-scale simulation framework that first employs finite-difference-time-domain simulations for calculating the surface fluence distribution on a rough surface measured by atomic-force-microscopy followed by the two-temperature model coupled with hydrodynamic/solid mechanics simulation for the initial material heating. Lastly, a computational fluid dynamics model for material relaxation and fluid flow is developed and employed. Final state results of aluminum and AISI 304 stainless steel simulations demonstrated alignment with established ablation models and crater dimension prediction. Notably, Al exhibited significant optical scattering effects due to initial surface roughness of 15 nm—being 70 times below the laser wavelength -leading to localized, selective ablation processes and substantially altered crater topography compared to idealized conditions. Contrary, AISI 304 with Rq surface roughness of 2 nm showed no difference. Hence, we highlight the necessity of incorporating realistic, material-specific surface roughness values into large-scale ablation simulations. Furthermore, the induced local fluence variations demonstrated the inadequacy of neglecting lateral heat transport effects in this context.

Faster Sampling in Molecular Dynamics Simulations with TIP3P-F Water.

The need for short time steps currently limits routine atomistic molecular dynamics (MD) simulations to the microsecond time scale. For long time steps, the numerical integration of the equations of motion becomes unstable, resulting in catastrophic crashes. Here, we combine mass repartitioning and rescaling to construct a water model that increases the sampling efficiency in biomolecular simulations without compromising integration stability and with preserved structural and thermodynamic properties. The resulting "fast water" is then used with a time step as before in combination with standard force fields. The reduced water viscosity and faster diffusion result in proportionally faster sampling of the larger-scale motions in the conformation space of both solute and solvent. We illustrate this approach by developing TIP3P-F based on the popular TIP3P model of water. A roughly 2-fold boost in the sampling efficiency at minimal cost in accuracy is substantial and helps lower the energy impact of large-scale MD simulations. The approach is general and can readily be applied to other water models and different types of solvents.

Uncertainty quantification in coupled wildfire-atmosphere simulations at scale.

Uncertainties in wildfire simulations pose a major challenge for making decisions about fire management, mitigation, and evacuations. However, ensemble calculations to quantify uncertainties are prohibitively expensive with high-fidelity models that are needed to capture today's ever-more intense and severe wildfires. This work shows that surrogate models trained on related data enable scaling multifidelity uncertainty quantification to high-fidelity wildfire simulations of unprecedented scale with billions of degrees of freedom. The key insight is that correlation is all that matters while bias is irrelevant for speeding up uncertainty quantification when surrogate models are combined with high-fidelity models in multifidelity approaches. This allows the surrogate models to be trained on abundantly available or cheaply generated related data samples that can be strongly biased as long as they are correlated to predictions of high-fidelity simulations. Numerical results with scenarios of the Tubbs 2017 wildfire demonstrate that surrogate models trained on related data make multifidelity uncertainty quantification in large-scale wildfire simulations practical by reducing the training time by several orders of magnitude from 3 months to under 3 h and predicting the burned area at least twice as accurately compared with using high-fidelity simulations alone for a fixed computational budget. More generally, the results suggest that leveraging related data can greatly extend the scope of surrogate modeling, potentially benefiting other fields that require uncertainty quantification in computationally expensive high-fidelity simulations.