Reacting Flow Simulations Research Articles

Particle deposition on the convergent nozzle walls of thermal-spray guns

The present study investigates particle depositions on convergent nozzle walls of (HVOF) spray guns, a process that reduces nozzle lifetime and degrades operation stability and reliability. We perform three-dimensional, two-way coupled simulations of reactive particle-laden compressible flow through a typical thermal spray gun nozzle. We determine critical radial particle velocities at which particles deposit on the wall. We characterize critical velocities concerning particle trajectories, particle diameter, and downstream location of particles deviating from the nozzle center. The critical velocity decreases with particle diameter and increases with the downstream location at which particles deviate from the nozzle center. We exclude the possibility that single particle–particle collision events create sufficient radial particle momentum for deposition. We use the gained insights to propose mitigating measures against particle depositions.

Adaptation of the Low Dissipation Low Dispersion Scheme for Reactive Multicomponent Flows on Unstructured Grids Using Density-Based Solvers

Abstract The low dissipation low dispersion (LD2) second-order accurate scheme is effective for scale-resolving simulations in finite volume computational fluid dynamics (CFD) solvers, thanks to its combination of a skew-symmetric split form for convective terms and matrix-valued artificial dissipation fluxes. However, reactive flow simulations face challenges due to steep gradients and varying gas properties. This study replaces the skew-symmetric scheme with the kinetic energy and entropy preserving (KEEP) scheme, utilizing quadratic and cubic split forms for convective terms, enhancing stability. The nonsmooth fluid interfaces in reactive flow simulations necessitate upwind fluxes for reactive scalars to limit total variation (TV), also requiring upwind fluxes for the mixture-dependent internal energy fluxes. Other convective terms use central discretizations from the KEEP scheme, leveraging LD2’s spatial reconstruction to minimize dispersive errors. Numerical assessments show this approach reduces spurious pressure oscillations in single and multicomponent flows. The absolute flux Jacobian for dissipation flux calculation is efficiently computed using an expanded Turkel’s approach for thermally perfect gas mixtures. Partial pressure derivatives are approximated when using the flamelet generated manifolds (FGM) combustion model. The proposed scheme is evaluated through scale-resolving simulations of the Cambridge burner flame SWB1 on an unstructured grid using the density-based solver TRACE, employing both finite rate chemistry (FRC) and FGM combustion models. Comparative analysis with the all-speed SLAU2 scheme shows the superior performance of the proposed scheme in handling turbulent reactive multicomponent flows.

Double flame dynamics in hotspot ignition

Hotspot ignition is a key problem in fundamental flame initiation, combustion control and prevention of abnormal combustion phenomena. Depending on the characteristics of the hotspot and the thermochemical condition of the reactive mixture, subsequent reaction front initiated from a hotspot can vary substantially. The current paper focus on the dynamics and chemistry for one of the most complicated scenarios - hotspot induced double flames, which include a cool and a hot flame segment. This work adopts a recently developed compressible reacting flow simulation code COGNAC to investigate the initiation and propagation of double flame dynamics from hotspot ignition in a one-dimensional spherical coordinate. The results have shown that under atmospheric pressure, double flame initiation is not feasible for sufficiently small hotspot. Under elevated pressures, double flame can be initiated within a much wider hotspot window, exhibiting complex flame dynamics, such as variations in cool and hot flame bifurcation behaviors and their initiation locations. Analysis on thermochemical structure and dynamics of double flame has shown that the interaction of cool and hot flames in a double flame structure is inherently intermutual, in that the leading cool flame enhances trailing hot flame speed through reform of the unburnt mixture, while the trailing hot flame initiation affects the leading cool flame via thermal expansion effects.

Accurate and rapid reactive flow simulations using dynamic load balancing and sparse analytical Jacobian approach

The present study introduces a rapid and accurate customized solver on the OpenFOAM platform for large-scale industrial computations. Specifically, a sparse analytical Jacobian approach utilizing the SpeedCHEM library was implemented to enhance the efficiency of the ordinary differential equation solver. The dynamic load balancing code was used to distribute computational workloads uniformly across multiple processes. Optimization continued with open multi-processing to improve parallel computing efficiency and the local time stepping scheme to maximize individual cell time steps. The effectiveness and robustness of the customized solver were first validated using Sandia flames D–F as benchmarks. The results showed that the customized solver exhibited better strong scaling characteristics and led to a speed increase of up to 30 times for two-dimensional Sandia flame D calculations. The numerical predictions for temperature and species distribution closely matched the experimental trends, confirming the accuracy of the solver. Subsequently, a three-dimensional numerical study on a 10 kW ammonia co-combustion furnace was conducted, exploring the performance of the solver in large-scale reactive simulations. Results analysis indicated that the acceleration capability was reduced due to increased communication overhead between processors, achieving up to 7.06 times speed-up. However, as the size of the reaction mechanism increases, better acceleration capabilities can be demonstrated. The numerical predictions could closely replicate experimental trends, effectively predicting NO emission trends within the combustion furnace. This study offers one viable solution for rapid and accurate calculations in the OpenFOAM platform, which could be applied in the subsequent ammonia industrial combustion processes.

A chemistry load balancing model for OpenFOAM

Efficient simulation tools are crucial for studying complex systems such as reacting flows where computational costs of computing chemical reaction rates can vastly exceed the costs for the integration of the convective and diffusive transport terms. Load imbalance in parallel computing poses a significant challenge for massively parallel reacting flow simulations. In response, a novel load balancing library has been developed to enhance OpenFOAM's solver performance in parallel environments. This library seamlessly integrates with OpenFOAM, offering ease of use and applicability to any OpenFOAM reacting solver incorporating finite-rate chemistry. In addition, it supports the standard and the dynamic adaptive chemistry model (TDAC) of OpenFOAM.The newly developed load-balanced standard and TDAC models address significant load imbalances by exchanging information between processes via MPI calls and tracking ODE solution times on a cell level. The TDAC model introduces dual tables on each core and enables immediate addition of computed solutions, enhancing computational efficiency. Validation on various test cases, including simulations on the HLRS Hawk supercomputer with up to 8000 cores, confirms identical results compared to the original unbalanced models, with notable speed-up factors of up to 6 for the standard and 5 for the TDAC model. Despite non-linear scaling at lower cell count per processor, load-balanced models consistently outperform unbalanced counterparts, making them the preferred choice for reacting flow simulations in OpenFOAM.

Modeling of turbulence, chemistry, and heat transfer of a rocket nozzle

Calculating the composition of the combustion products serves as the starting point for the high-area ratio (HAR) rocket nozzle. The HAR rocket nozzle finds application in next-generation hypersonic vehicles and missiles. The gas dynamics, gas radiation, and turbulence-chemistry interactions (TCI) play a crucial role in simulations of hypersonic turbulent reactive flows through a rocket nozzle. The interaction of turbulence and chemistry leads to random species mixing and enhanced heat transfer. The k-omega shear stress transport (SST) model is employed to calculate the average flow fields. The eddy dissipation concept (EDC) with reduced chemistry is applied to capture the TCI in the nozzle. The spectral properties of H2O are calculated using the full spectrum k-distribution (FSK) model in conjunction with the first-order spherical harmonics method. The mass fractions of combustion products obtained from Particle Swarm Optimization (PSO) and Monte Carlo (MC) methods agree well with the Chemical Equilibrium with Applications (CEA) code. The flow fields obtained with the current solver are validated with the published experimental data. A significant rise and fall of 7–8% have been observed in the centerline flow properties when not including the TCI models in the flow solver. The centerline flow fields obtained with different discretization schemes are in good agreement. The radiative heat flux from the H2O species is more dominant than the convective.

Investigation On A Novel Injector Concept For Spinning Combustion Technology In High-Pressure Conditions By Advanced Laser-Based Diagnostics

Safran Helicopter Engines has recently patented the spinning combustion technology in which the burnt gases from one injector travel tangentially along the combustor annulus towards the neighboring injectors. Compared to conventional designs, the new kerosene injection systems are dedicated to improve air/fuel mixture ignition but also to further reduce NOx and soot particle emissions. Experimental studies are performed on these fuel injectors in a high-pressure/high-temperature combustion facility designed by the CORIA research laboratory. This test bench is able to reproduce the same operating conditions encountered in a helicopter combustor over the entire range of nominal operating conditions and has large optical accesses for the implementation of laser-based diagnostics. In the current paper, we present results concerning flame structure and NO formation in the primary zone under pressure conditions of up to 14 bar, using simultaneous OH-PLIF, NO-PLIF and kerosene-PLIF laser diagnostics. These experimental studies were supplemented by high-speed PIV measurements. A good spatial correlation between the distribution of liquid and vapour kerosene and the location of the flame front was observed. Depending on the operating conditions in terms of fuel/air ratio, mass flow rates and pressure, different flame structures resulting from the modification of the interaction between fuel injection and aerodynamics are observed. Furthermore, it was found that the Zeldovich pathway mainly controls the formation of NO in the vicinity of the flame front. In addition, the effects of FAR and pressure also have a significant impact on NO production. All these results are now intended to serve as a comprehensive validation database for the development and testing of high-fidelity LES tools dedicated to the simulation of reactive flows in aero-engine combustion chambers.

Jacobian-scaled K-means clustering for physics-informed segmentation of reacting flows

This work introduces Jacobian-scaled K-means (JSK-means) clustering, which is a physics-informed clustering strategy centered on the K-means framework. The method allows for the injection of underlying physical knowledge into the clustering procedure through a distance function modification: instead of leveraging conventional Euclidean distance vectors, the JSK-means procedure operates on distance vectors scaled by matrices obtained from dynamical system Jacobians evaluated at the cluster centroids. The goal of this work is to show how the JSK-means algorithm – without modifying the input dataset – produces clusters that capture regions of dynamical similarity, in that the clusters are redistributed towards high-sensitivity regions in phase space and are described by similarity in the source terms of samples instead of the samples themselves. The algorithm is demonstrated on a complex reacting flow simulation dataset (a channel detonation configuration), where the dynamics in the thermochemical composition space are known through the highly nonlinear and stiff Arrhenius-based chemical source terms. Interpretations of cluster partitions in both physical space and composition space reveal how JSK-means shifts clusters produced by standard K-means towards regions of high chemical sensitivity (e.g., towards regions of peak heat release rate near the detonation reaction zone). The findings presented here illustrate the benefits of utilizing Jacobian-scaled distances in clustering techniques, and the JSK-means method in particular displays promising potential for improving former partition-based modeling strategies in reacting flow (and other multi-physics) applications.

Acceleration of the complex reacting flow simulation with a generalizable neural network based on meta-learning

A chemical decoupled method is proposed in this paper. To accelerate the complex reacting flow simulation, the chemistry and flow computation is decoupled by a Neural Network (NN). A chemical kinetics database including 2,106,496 flow states and corresponding reacting source terms has been generated. Nevertheless, the data presents imbalance and spans a wide range of magnitudes, which challenges the NN training. The normalization method, NN’s construction, and performance evaluation indicators are optimized to address this challenge. Another problem is that different processes, such as transport or reaction, dominate different flow regions. This phenomenon requires the NN sufficient generalization to accurately map the complex reactions and flow states, even when the dominant process is different. The generalization of the NN is improved through the meta-learning strategy, the residual connection, and the inter-layer sine normalization. The effect of these enhancements is tested by noisy data and simulations with different turbulence model and boundary conditions. The method has been demonstrated in a supersonic combustion simulation with 16 species and 137 reactions. The maximum mean square relative error between the simulation with NN and the traditional simulation is 13.59 %. The calculation time is reduced from 6415.989 min (about 4.456 days) of the traditional simulation to 1983.576 min (about 1.377 days) of the simulation with NN. This method has application potential in accelerating the complex reacting flow simulation.

Numerical simulation of chemical reactive flow of Boger fluid over a sheet with heat source and local thermal non-equilibrium conditions

In this article, the effect of nonlinear heat source on chemical reactive flow of magnetized Boger fluid over a sheet with Marangoni convection and porous medium is explored. In the modelling, local thermal non-equilibrium (LTNE) conditions and heat sources are taken into account. The unique heat transport for both liquid and solid phases is provided by the LTNE model, which is based on thermal equations. As a result, various temperature profiles are used in this work for both the solid and fluid phases. It is essential for streamlining industrial procedures, improving material processing methods, comprehending biological systems like blood flow for medicine delivery, determining the environmental impact of pollutant dispersion, and increasing the effectiveness of oil and gas transportation. This model supports crucial decision-making processes in numerous disciplines by offering insights into the intricate interplay between reactive fluids, heat transmission, and chemical reactions. The appropriate similarity variables are used to compress the model equation system, and the shooting method and the conventional bvp4c solver are then used to solve the system numerically. Optimal homotopy approach is used to achieve residual errors. There are suggested best estimates for auxiliary variables. Results show that for larger values of the solvent fraction parameter, Boger fluid exhibits an enhanced velocity profile and a declining thermal profile of the liquid phase. The concentration and thermal profiles of both the solid and liquid phases are declined and enhanced the velocity of the Boger fluid for higher values of the Marangoni convection parameter.

Graphics processing unit/artificial neural network-accelerated large-eddy simulation of swirling premixed flames

Within the scope of reacting flow simulations, the real-time direct integration (DI) of stiff ordinary differential equations for the computation of chemical kinetics stands as the primary demand on computational resources. Meanwhile, as the number of transport equations that need to be solved increases, the computational cost grows more substantially, particularly for those combustion models involving direct coupling of chemistry and flow such as the transported probability density function model. In the current study, an integrated graphics processing unit-artificial neural network (GPU-ANN) framework is introduced to comply with heavy computational costs while maintaining high fidelity. Within this framework, a GPU-based solver is employed to solve partial differential equations and compute thermal and transport properties, and an ANN is utilized to replace the calculation of reaction rates. Large eddy simulations of two swirling flames provide a robust validation, affirming and extending the GPU-ANN approach's applicability to challenging scenarios. The simulation results demonstrate a strong correlation in the macro flame structure and statistical characteristics between the GPU-ANN approach and the traditional central processing unit (CPU)-based solver with DI. This comparison indicates that the GPU-ANN approach is capable of attaining the same degree of precision as the conventional CPU-DI solver, even in more complex scenarios. In addition, the overall speed-up factor for the GPU-ANN approach is over two orders of magnitude. This study establishes the potential groundwork for widespread application of the proposed GPU-ANN approach in combustion simulations, addressing various and complex scenarios based on detailed chemistry, while significantly reducing computational costs.

Multiphase reactive flow and heat transfer simulations in the unlocking components of an ultra-low shock pyrotechnic separation nut

This study presents a multiphase reactive flow and conjugate heat transfer model for the unlocking components in an ultra-low shock pyrotechnic separation nut. The process of the shape memory alloy (SMA) tube heated by the heat released from the combustion of the micro gas pyrotechnic charges was simulated. The established model was validated by comparing the simulation with experimental results for the temperature of the outer wall of the SMA tube. Then, the detailed ignition and combustion of micro gas pyrotechnic charges, transient flow of combustion products, and heat transfer in the SMA tube were analyzed. The effects of the amount and burning rate of micro gas pyrotechnic charges on the unlocking time were investigated. It is found that the unlocking time decreases with an increase in the amount of micro gas pyrotechnic charges. Moreover, the reduction range of the unlocking time decreases with the increase in the amount of micro gas pyrotechnic charges. Reducing the differences in the burning rate contributes to enhancing the consistency of unlocking time.

Delayed Detached-Eddy Simulations of Rough-Wall Turbulent Reactive Flows in a Supersonic Combustor

Reactive delayed detached-eddy simulations (DDESs) of hydrogen injection into a confined transverse supersonic flow of vitiated air are conducted. The corresponding conditions were studied in the LAPCAT-II combustor, which—due to thermal coating—exhibits non-negligible wall roughness. Its effects are taken into account with the equivalent sand-grain roughness model, and, to the best of the authors’ knowledge, it is the first time that this modeling framework is considered within the DDES framework to simulate turbulent reactive flows. Two operating conditions are considered, which differ in the value of the momentum ratio between the hydrogen and vitiated air inlet streams, thus leading to two distinct values of the operating equivalence ratio (ER). For its smallest value, smooth combustion is subject to a preliminary thermal runaway period, while for its largest value, combustion is more strongly intertwined with shock wave dynamics and boundary-layer separation. Since it features large subsonic regions, which are characteristic of situations close to thermal choking, the latter is referred to as the sudden or partially choked combustion mode. Computational results are assessed through comparisons with available experimental data for both operating conditions. The main mechanism through which wall roughness alters the combustion development lies in the reduction of the effective cross-sectional area that is induced by boundary-layer thickening. For the largest ER value, DDES conducted with the sand-grain roughness model does recover the partially choked combustion mode, while smooth wall (i.e., standard) DDES does not.

Goal‐Oriented Two‐Layered Kernel Models as Automated Surrogates for Surface Kinetics in Reactor Simulations

AbstractMulti‐scale modeling allows the description of real reactive systems under industrially relevant conditions. However, its application to rational catalyst and reactor design is hindered by the prohibitively high computational cost associated with the chemical kinetics on the catalyst scale. Here, the computational cost is drastically reduced by introducing goal‐oriented kernel models that serve as surrogates for the chemical kinetics. This special model type allows for automated training set design and reliable results, even outside the training region. Therefore, it can be envisioned as a plug‐and‐play solution for accelerating reactive flow simulations with guaranteed accuracy.

Increasing Computational Efficiency of CFD Simulations of Reactive Flows at Catalyst Surfaces through Dynamic Load Balancing.

We propose a numerical strategy based on dynamic load balancing (DLB) aimed at enhancing the computational efficiency of multiscale CFD simulation of reactive flows at catalyst surfaces. Our approach employs DLB combined with a hybrid parallelization technique, integrating both MPI and OpenMP protocols. This results in an optimized distribution of the computational load associated with the chemistry solution across processors, thereby minimizing computational overheads. Through assessments conducted on fixed and fluidized bed reactor simulations, we demonstrated a remarkable improvement of the parallel efficiency from 19 to 87% and from 19 to 91% for the fixed and fluidized bed, respectively. Owing to this improved parallel efficiency, we observe a significant computational speed-up of 1.9 and 2.1 in the fixed and fluidized bed reactor simulations, respectively, compared to simulations without DLB. All in all, the proposed approach is able to improve the computational efficiency of multiscale CFD simulations paving the way for a more efficient exploitation of high-performance computing resources and expanding the current boundaries of feasible simulations.

Turbulence-radiation interaction in large eddy simulations of turbulent reactive flows: a review

Turbulence-radiation interaction (TRI) arises from the highly non-linear coupling between fluctuations of radiative intensity and fluctuations of temperature and chemical composition of the medium. It is well recognized that neglecting turbulent fluctuations of the latter scalars when computing the mean radiation field, as is standard in Reynolds-Averaged Navier-Stokes (RANS) simulations, leads to mispredictions of the radiative heat fluxes and the radiant fraction in as much as one order of magnitude. Conversely, the importance of TRI in large eddy simulations (LES)—i.e., the effect of subfilter-scale (SFS) fluctuations on the filtered radiation field—has only started to be investigated in the last decade. The present text reviews the knowledge gathered so far on the subject. The instantaneous and the spatially filtered forms of the radiative transfer equation (RTE) are presented, and the unresolved correlations resulting from the filtering operation are discussed. Models for accounting for SFS-TRI are surveyed, with the most popular ones (namely, probability density function-based methods and empirical approximations) reviewed in detail. An appraisal of the literature on the importance and characteristics of SFS-TRI is performed, with special attention given to turbulent jet diffusion flames and pool fires, which are the two configurations for which SFS-TRI has been studied the most. The main findings of these studies are summarized, and recommendations are given on directions of future works on the subject.

A Computational Study of the Influence of Drag Models and Heat Transfer Correlations on the Simulations of Reactive Polydisperse Flows in Bubbling Fluidized Beds

In this work, the influence of gas–solid drag and heat transfer coefficient models on the prediction capacity of the Euler–Euler approach to simulate reactive bubbling fluidized bed flows is studied. Three different cases are considered, a non-reactive bidisperse bubbling fluidized bed flow (Case 1), and two reactive polydisperse flows in bubbling fluidized beds, one for biomass gasification (Case 2), and the other for biomass pyrolysis (Case 3). The Gidaspow, Syamlal–O’Brien, and BVK gas–solid drag models and the Gunn, Ranz–Marshall, and Li–Mason gas–solid heat transfer correlations are investigated. A Eulerian multiphase approach in a two-dimensional Cartesian domain is employed for the simulations. Computational results for the three cases are compared with experimental data from the literature. The results obtained here contribute to a better understanding of the impacts of such closure models on the prediction ability of the Euler–Euler approach to simulate reactive flows. The results indicate that, for the simulation of reactive flows in bubbling fluidized bed reactors, the kinetic modeling of the reactions has a global effect, which superposes with the influence of the drag and heat transfer coefficient models. Nevertheless, local parameters can be noticeably affected by the choice of the interface closure models. Finally, this work also identifies the models that lead to the best results for the cases analyzed here, and thus proposes the use of such selected models for gasification and pyrolysis processes occurring in bubbling fluidized bed reactors.

Compressible Large Eddy Simulation of Thermoacoustic Instabilities in the PRECCINSTA Combustor Using Flamelet Generated Manifolds With Dynamic Thickened Flame Model

Abstract The fully compressible, density-based CFD-solver TRACE has been extended for simulations of turbulent reacting flows in aero engine gas turbine combustors. The flamelet generated manifolds combustion model is utilized to account for detailed chemical kinetics and combined with the dynamically thickened flame model to resolve the flame front on the large eddy simulation (LES) mesh. The chemistry tabulation is coupled with the LES solver by inversion of the transported energy equation using tabulated mixture averaged NASA polynomial coefficients. LES of the PRECCINSTA test case, a lean, partially premixed swirl combustor are performed and the two distinctive regimes are correctly predicted: a stable regime with a “quite” stable flame and an unstable regime with an oscillating flame driven by self-excited thermoacoustic instabilities. Statistics collected from the simulations, mean, and root-mean-square values are in good agreement with the experimental reference data for both operating conditions. The dominant frequency of the unstable flame deviates from the measurement by about 100 Hz and requires further investigation. The results demonstrate the general suitability of the simulation framework for reacting flow simulations in gas turbine combustion systems and the prediction of self-excited thermoacoustic oscillations.

Transient Thermomechanical Stress Analysis of Hot Surface Ignition Device Using Sequentially Coupled Computational Fluid Dynamics–Finite Element Analysis Approach

Abstract Energy addition using a hot surface ignition device is required for reliable ignition of aircraft compression ignition engines running on fuel variations and at altitude conditions. Thus, durability of the hot surface ignition device is crucial for application in these engines. Thermomechanical stress is one of the key parameters that determine durability, which requires an accurate prediction of the transient temperature field based on well-defined boundary conditions representing the dynamic and complex fluid flow inside engines. To meet this requirement, the present study focuses on transient thermomechanical stress analysis using a sequentially coupled computational fluid dynamics (CFD)–finite element analysis (FEA) approach to understand transient thermomechanical responses of the hot surface ignition device. A three-dimensional transient reacting flow simulation was conducted first using converge software, the results of which were exported to map thermal and pressure boundary conditions onto a structural finite element mesh. Transient thermomechanical stress analysis was performed sequentially using abaqus software utilizing the mapped boundary conditions. The results such as transient temperature history, resultant thermomechanical stress, displacement, potential failure modes, etc., were critically reviewed, which can provide helpful information for further design improvement.

Adaptive physics-informed neural operator for coarse-grained non-equilibrium flows

This work proposes a new machine learning (ML)-based paradigm aiming to enhance the computational efficiency of non-equilibrium reacting flow simulations while ensuring compliance with the underlying physics. The framework combines dimensionality reduction and neural operators through a hierarchical and adaptive deep learning strategy to learn the solution of multi-scale coarse-grained governing equations for chemical kinetics. The proposed surrogate’s architecture is structured as a tree, with leaf nodes representing separate neural operator blocks where physics is embedded in the form of multiple soft and hard constraints. The hierarchical attribute has two advantages: (i) It allows the simplification of the training phase via transfer learning, starting from the slowest temporal scales; (ii) It accelerates the prediction step by enabling adaptivity as the surrogate’s evaluation is limited to the necessary leaf nodes based on the local degree of non-equilibrium of the gas. The model is applied to the study of chemical kinetics relevant for application to hypersonic flight, and it is tested here on pure oxygen gas mixtures. In 0-textrm{D} scenarios, the proposed ML framework can adaptively predict the dynamics of almost thirty species with a maximum relative error of 4.5% for a wide range of initial conditions. Furthermore, when employed in 1-textrm{D} shock simulations, the approach shows accuracy ranging from 1% to 4.5% and a speedup of one order of magnitude compared to conventional implicit schemes employed in an operator-splitting integration framework. Given the results presented in the paper, this work lays the foundation for constructing an efficient ML-based surrogate coupled with reactive Navier-Stokes solvers for accurately characterizing non-equilibrium phenomena in multi-dimensional computational fluid dynamics simulations.