Finite Volume Cells Research Articles

Rapid spatio-temporal flood modelling via hydraulics-based graph neural networks

Abstract. Numerical modelling is a reliable tool for flood simulations, but accurate solutions are computationally expensive. In recent years, researchers have explored data-driven methodologies based on neural networks to overcome this limitation. However, most models are only used for a specific case study and disregard the dynamic evolution of the flood wave. This limits their generalizability to topographies that the model was not trained on and in time-dependent applications. In this paper, we introduce shallow water equation–graph neural network (SWE–GNN), a hydraulics-inspired surrogate model based on GNNs that can be used for rapid spatio-temporal flood modelling. The model exploits the analogy between finite-volume methods used to solve SWEs and GNNs. For a computational mesh, we create a graph by considering finite-volume cells as nodes and adjacent cells as being connected by edges. The inputs are determined by the topographical properties of the domain and the initial hydraulic conditions. The GNN then determines how fluxes are exchanged between cells via a learned local function. We overcome the time-step constraints by stacking multiple GNN layers, which expand the considered space instead of increasing the time resolution. We also propose a multi-step-ahead loss function along with a curriculum learning strategy to improve the stability and performance. We validate this approach using a dataset of two-dimensional dike breach flood simulations in randomly generated digital elevation models generated with a high-fidelity numerical solver. The SWE–GNN model predicts the spatio-temporal evolution of the flood for unseen topographies with mean average errors in time of 0.04 m for water depths and 0.004 m2 s−1 for unit discharges. Moreover, it generalizes well to unseen breach locations, bigger domains, and longer periods of time compared to those of the training set, outperforming other deep-learning models. On top of this, SWE–GNN has a computational speed-up of up to 2 orders of magnitude faster than the numerical solver. Our framework opens the doors to a new approach to replace numerical solvers in time-sensitive applications with spatially dependent uncertainties.

MMC-LES of spray combustion: Analysis of mixing timescales and flame structure

Large-eddy simulations (LES) of turbulent, pilot-stabilised, dilute acetone spray flames are performed using the sparse-Lagrangian multiple mapping conditioning approach (MMC-LES). The primary objective of this work is to assess the adaptability of the recently proposed dynamically evaluated molecular mixing timescale model for MMC-LES (Sharma et al., 2022, A fully dynamic mixing timescale model for the sparse Lagrangian multiple mapping conditioning approach, Combustion and Flame, 238, 111872), referred to as the dyn-aISO model, to the non-reacting and reacting two-phase flows. This work is the first application to perform a dynamic evaluation of the mixing timescale for turbulent spray combustion using the transported filtered density function (FDF) approach. A non-reacting case and four reacting flames with varying fuel mass loading are investigated, where the flame transitions from a diffusion flame to a premixed flame structure with a decrease in fuel loading. As is conventional for MMC-LES, a sparse resolution of one stochastic particle per six Eulerian finite-volume cells is used in this study, offering a much cheaper computing expense than the conventional transported FDF approach. The liquid evaporation is represented based on a non-equilibrium model at the droplet surface. The reaction chemistry is based on an augmented GRI 3.0 mechanism containing 39 species and 226 reactions to represent the oxidation of acetone. Predictions with the dyn-aISO model are compared with the constant coefficient anisotropic mixing time scale model and the isotropic C-K model. Thorough validation of liquid-phase statistics and gas-phase temperature is performed with available experimental data and popular studies from the literature. Analysis of the flame structure, with different timescale models, for both non-premixed and premixed type flames is also carried out. Quantification of the parameters affecting the timescale, such as model constants of the sub-Lagrangian-grid scale (slgs) mixture fraction variance (Cf) and its scalar dissipation rate (Scsgs−1Cν) is provided in comparison to their conventionally chosen static coefficient values. Predictions with the dyn-aISO model are quite good for the diffusion-type flames, while some discrepancies are noticed while predicting a premixed flame structure for lower liquid loading case having near stoichiometric jet-exit mixture composition.

A stable interface-preserving reinitialization equation for conservative level set method

Existing artificial compression based reinitialization approach for conservative level set method has a few drawbacks, like distortion of fluid–fluid interface and formation of unphysical fluid patches away from the interface. In this paper, a novel reinitialization approach has been presented which circumvents these limitations by reformulating the existing reinitialization equation. The new reinitialization equation by design ensures that the interface contour is kept invariant under the action of level set reinitialization. The unphysical patch formation away from the interface is also resolved here by avoiding the use of ill-conditioned contour normal vectors. Owing to the use of much simpler terms, a significant reduction in the numerical computations is achieved with the new reinitialization equation. The new reinitialization equation also enables one to choose a larger time step during the reinitialization iteration, leading to an overall reduction in computational efforts. In order to evaluate the performance of the present formulation, a set of test problems involving reinitialization of stationary level set functions is carried out first. Then, the efficacy of the proposed reinitialization approach is demonstrated using a set of standard two-dimensional scalar advection based test problems and incompressible two-phase flow problems. Finally, the ability to deal with complex mesh types is demonstrated by solving a test problem on an unstructured mesh consisting of finite volume cells having triangular and quadrilateral shapes.

Dynamic modelling of the steady state and load processing operation of a domestic refrigerator cooled through natural convection

A dynamic model of a domestic refrigerator is developed in this study and is used to analyse both the steady state (cyclic) and warm load cooling (load processing) operation. The compressor is modelled semi-empirically using loss parameters, while the capillary tube is modelled by incorporating both friction and momentum pressure drop. Heat exchangers (condenser and evaporator) are modelled as tubes divided into finite volume cells. The capillary tube-suction line heat exchanger is treated by connecting part of the capillary tube to the suction line via a heat resistor. The cabinet model consists of submodels for walls, air and shelves connected in a network to allow for thermal interaction. All the components are then arranged in a cycle, and cyclic operation is validated at 32 °C using a fixed speed compressor. A sensitivity analysis is then carried out by varying the compressor speed. In the second part of the study, the model is modified to include the effect of placing a warm load of water in the refrigerator and simulating the subsequent processing operation. This is validated using an experiment in which 20 L of water at 27 °C is placed inside the Fresh Food compartment. Results show that the power consumption value is captured well, with the average power consumption deviating less than 10% for all tests. The energy consumption is within ±5% for all tests. For load processing, the deviation in water temperature was a maximum of 1.2 °C for the whole cooling period.

Inverse distance interpolation for used in unstructured mesh finite volume solver

This article discusses adjusting inverse distance interpolation for use in unstructured mesh finite volume solutions. The adjustment was made on the weight function of the inverse distance interpolation using the Laplacian of the flow variable inside a Voronoi-dual of finite volume cells. We tested the accuracy of the adjusted inverse distance interpolation on two-dimensional potential flows. It was found that the adjusted and standard inverse distance interpolations have a similar degree of accuracy when used in unstructured, Delaunay based, finite volume mesh. However, the L1 norm error of the adjusted version of the inverse distance interpolation was much smaller than the L1 norm error of the standard version.

A multiresolution local-timestepping scheme for particle–laden multiphase flow simulations using a level-set and point-particle approach

We present a novel storage-and-access approach for Lagrangian particles in a multiresolution local-timestepping framework for multiphase flow simulations. The proposed method applies a block-particle mapping strategy for efficient access of all particles on a specific refinement level while traversing through the multiresolution tree. This allows for extending the local timestepping algorithm with its refinement level-dependent timestep sizes to the evolution of particles, which results in significant speed-up in comparison to standard timestepping. For multi-scale multiphase flow simulations, the particle model is combined with a level-set based multiphase model on a Cartesian grid. To maintain robustness of fluid-state interpolation for particles near the level set-based fluid–fluid interface, WENO-based interpolation is applied which includes both real- and ghost-fluid cells. This enforces the sharp interface property also for interpolated fluid states, and suppresses spurious oscillations in the event of discontinuities.We validate the particle model with a one-dimensional simulation of a single particle in an air–helium shock tube for one-way coupling, and with two-dimensional simulations of a particle injected in a quiescent domain for the feedback force. Simulations of two-dimensional aerodynamic fragmentation in shear-induced entrainment and Rayleigh–Taylor piercing regimes use Lagrangian particles as sub-grid scale representation of small droplets post-breakup. Finally, three-dimensional massively-parallel simulations of single- and triple-bubble collapse near a wall coated with free-floating particles are presented. The particle cleaning radius of the single-bubble setup agrees reasonably well with experimental reference data. These simulations consider 106 particles and an Eulerian grid with effectively 1012 finite-volume cells at compression rates of more than 90% for the particulate phase. This underlines the advantageous effect of embedding the particles in the multiresolution tree with its spatial and temporal adaptivity, which is necessary for performing such large scale simulations efficiently.

Multiresolution analysis as a criterion for effective dynamic mesh adaptation – A case study for Euler equations in the SAMR framework AMROC

Dynamic mesh adaptation methods require suitable refinement indicators. In the absence of a comprehensive error estimation theory, adaptive mesh refinement (AMR) for nonlinear hyperbolic conservation laws, e.g. compressible Euler equations, in practice utilizes mainly heuristic smoothness indicators like combinations of scaled gradient criteria. As an alternative, we describe in detail an easy to implement and computationally inexpensive criterion built on a two-level wavelet transform that applies projection and prediction operators from multiresolution analysis. The core idea is the use of the amplitude of the wavelet coefficients as smoothness indicator, as it can be related to the local regularity of the solution. Implemented within the fully parallelized and structured adaptive mesh refinement (SAMR) software system AMROC (Adaptive Mesh Refinement in Object-oriented C++), the proposed criterion is tested and comprehensively compared to results obtained by applying the scaled gradient approach. A rigorous quantification technique in terms of numerical adaptation error versus used finite volume cells is developed and applied to study typical two- and three-dimensional problems from compressible gas dynamics. It is found that the proposed multiresolution approach is considerably more efficient and also identifies – besides discontinuous shock and contact waves – in particular smooth rarefaction waves and their interaction as well as small-scale disturbances much more reliably. Aside from pathological cases consisting solely of planar shock waves, the majority of realistic cases show reductions in the number of used finite volume cells between 20 to 40%, while the numerical error remains basically unaltered.

Adaptive moving mesh upwind scheme for the two-species chemotaxis model

Chemotaxis systems are used to model the propagation, aggregation and pattern formation of bacteria/cells in response to an external stimulus, usually a chemical one. A common property of all chemotaxis systems is their ability to model a concentration phenomenon—rapid growth of the cell density in small neighborhoods of concentration points/curves. More precisely, the solution may develop singular, spiky structures, or even blow up in finite time. Therefore, the development of accurate and computationally efficient numerical methods for the chemotaxis models is a challenging task.We study the two-species Patlak–Keller–Segel type chemotaxis system, in which the two species do not compete, but have different chemotactic sensitivities, which may lead to a significantly difference in cell density growth rates. This phenomenon was numerically investigated in Kurganov and Lukáčová-Medviďová (2014) and Chertock et al. (2018), where second- and higher-order methods on uniform Cartesian grids were developed. However, in order to achieve high resolution of the density spikes developed by the species with a lower chemotactic sensitivity, a very fine mesh had to be utilized and thus the efficiency of the numerical method was affected.In this work, we consider an alternative approach relying on mesh adaptation, which helps to improve the approximation of the singular structures evolved by chemotaxis models. We develop, in particular, an adaptive moving mesh (AMM) finite-volume semi-discrete upwind method for the two-species chemotaxis system. The proposed AMM technique allows one to increase the density of mesh nodes at the blowup regions. This helps to substantially improve the resolution while using a relatively small number of finite-volume cells.

A conservative interface-interaction method for compressible multi-material flows

In this paper we develop a conservative interface-interaction method dedicated to simulating multiple compressible fluids with sharp interfaces. Numerical models for finite-volume cells cut by more than two material-interface are proposed. First, we simplify the interface interaction inside such a cell to avoid the need for explicit interface reconstruction and very complex flux calculation. Second, conservation is strictly preserved by an efficient conservation correction procedure for the cut cell. To improve robustness, a multi-material scale separation model is developed to remove consistently non-resolved interface scales. In addition, a multi-resolution method and a local time-stepping scheme are incorporated into the proposed multi-material method to speed up high-resolution simulations. Various numerical test cases, including the multi-material shock tube problem, inertial confinement fusion implosion, triple-point shock interaction and shock interaction with multi-material bubbles, show that the method is suitable for a wide range of complex compressible multi-material flows.

Assessment of circulation and inter-basin transport in the Salish Sea including Johnstone Strait and Discovery Islands pathways

The Salish Sea consisting of Puget Sound and Georgia Basin in U.S and Canadian waters has been the subject of several independent data collection and modeling studies. However, these interconnected basins and their hydrodynamic interactions have not received attention as a contiguous unit. The Strait of Juan de Fuca is the primary pathway through which Pacific Ocean water enters the Salish Sea but the role played by Johnstone Strait and the complex channels northeast of Vancouver Island, connecting the Salish Sea and the Pacific Ocean, on overall Salish Sea circulation has not been characterized. In this paper we present a modeling-based assessment of the two-layer circulation and transport through the multiple interconnected sub-basins within the Salish Sea including the effect of exchange via Johnstone Strait and Discovery Islands. The Salish Sea Model previously developed using the finite volume community ocean model (FVCOM) was expanded over the continental shelf for this assessment encircling Vancouver Island, including Discovery Islands, Johnstone Strait, Broughton Archipelago and the associated waterways. A computational technique was developed to allow summation of volume fluxes across arbitrary transects through unstructured finite volume cells. Tidally averaged volume fluxes were computed at multiple transects. The results were used to validate the classic model of Circulation in Embracing Sills for Puget Sound and to provide quantitative estimates of the lateral distribution of tidally averaged transport through the system. Sensitivity tests with and without exchanges through Johnstone Strait demonstrate that it is a pathway for Georgia Basin runoff and Fraser River water to exit the Salish Sea and for Pacific Ocean inflow. However the relative impact of this exchange on circulation and flushing in Puget Sound Basin is small.

A High-Order Piecewise Polynomial Reconstruction for Finite-Volume Methods Solving Convection and Diffusion Equations

In finite-volume methods for fluid flows, the average of field variables over local mesh cells are the unknowns that are integrated in time based on the integral conservation laws. In order to compute the cell-face fluxes, nodal variables and derivative values at cell faces are needed, which ultimately determine the accuracy of the finite-volume method. In this work, a piecewise fourth-order polynomial reconstruction model based on volume averages is developed for smoothly varying flow fields over finite-volume cells. The obtained polynomial is used to extrapolate the cell-face variables and derivative values with high-order accuracy. The extrapolated cell-face quantities are used directly to compute the convection or diffusion integrals to construct high-order finite-volume methods. Some one-dimensional examples are shown to demonstrate the fifth-order accuracy of the proposed approach when solving the linear convection equation, and the fourth-order solution accuracy when solving the linear diffusion equation.

An operator-splitting scheme for the stream function–vorticity formulation of the unsteady Navier–Stokes equations

The stream function–vorticity formulation of the (Navier–)Stokes equations yields a coupled system of a parabolic equation for the vorticity and an elliptic equation for the stream function. The essential coupling between them occurs through the boundary conditions which in case of a Dirichlet boundary involve only the stream function. Therefore, the boundary condition for the vorticity must be derived from them and thus the vorticity equation must be coupled to the stream function equation via its boundary condition. In this paper we propose an unconditionally stable splitting scheme for the unsteady Stokes equations in a stream function–vorticity formulation, that decouples the vorticity and stream function computations at each time step. The spatial discretization is based on a finite volume discretization on (generally) unstructured Delaunay grids and corresponding Voronoi finite volume cells. A generalization of the well-known Thom vorticity boundary condition is derived for such grids and the corresponding discrete problem is decoupled by a two-step splitting scheme which results in a decoupled discrete parabolic problem for the vorticity and an elliptic problem for the stream function. Furthermore, the scheme is extended to the unsteady Navier–Stokes equations. Finally, the stability and accuracy of the resulting schemes are demonstrated on numerical examples.

Unconstrained Voronoi Grids for Densely Spaced Complex Wells in Full-Field Reservoir Simulation

Summary Accurate representation of near-well flow is an important subject matter in reservoir simulation. In today's field-scale reservoir simulation, cell-centered structured grids remain the predominant practice. Typically, well-inflow performance of the perforated cells is connected to the finite-volume solution by means of well indices that may not be well-defined when the wellbore intersects the finite-volume cells in a complex trajectory. Fine gridding is also required to resolve the flow dynamics in the near-well regions. Strong grid-orientation sensitivities can also contribute to the numerical errors and may require significant local grid refinement to alleviate. There are ongoing resesarch-and-development (R&D) efforts on applying unstructured grids to better represent the near-well flow in reservoir simulation, but their applications are mainly in single-well study or sector modeling with a few wells. Some of the reasons cited for this include (1) the lack of an effective, easy-to-use full-field complex well-gridding tool; (2) the lack of supporting unstructured workflow for full-cycle reservoir simulation; (3) the cost of unstructured-grid simulation; and (4) the availability of post-analysis and visualization tools for unstructured-grid simulation. The paper describes a novel method to automatically generate unstructured grids that conform to complex well paths in field-scale simulation. The method uses a multilevel approach to place cells optimally within the solution domain on the basis of the “regions of interests.” The wellbore geometry is honored by means of the construction of a near-well grid that is complemented with multilevel quad-tree (Fig. 1) refinements to achieve the desired resolution in grid transition zones. The method includes an algorithm to remove small cells and pinching cells on the basis of local grid quality measures and cell prioritization to honor well paths. The gridding process forms a component of a production-level reservoir-simulation workflow. The use of unstructured grid results in computational savings by placing cells where the resolution is needed. An in-house massively parallel simulator is used to run the unstructured-grid models. Simulation examples for full-field applications with hundreds of complex wells by use of both structured grids and unstructured grids will be used to compare results, accuracy, and performance of the gridding method for reservoir simulation.

A novel application of barycentric coordinates to Lagrange two-way coupled spray calculations

Modern CFD codes make possible the application of two-way coupled Lagrange particle tracking to complex unstructured meshes. State of the art tracking algorithms require the handling of arbitrary-shaped finite volume cells. By developing a new universal module for particle tracking, a novel technique has been introduced: association of the instantaneous particle location to the nearest vertex within a cell. The nearest-node method is then used to store the source terms arising from the spray calculation into vertex-based fields. In this novel method, barycentric coordinates are used for defining the vicinity of a vertex within a cell. The major benefits of this approach are its high computational efficiency and its ability to suppress the influence of mesh skewness on the quality of the results.

Evaluation of diffusive flux across faces of arbitrary shaped finite volume cells

A new procedure based on principles of directional derivatives of vector calculus has been presented for the computation of normal gradient across faces of arbitrarily skewed finite volume cells. The equations developed are compact, computationally less expensive and enable easy development of set of algebraic equations for implicit methods. The uniqueness and boundedness of the derivatives have been examined and confirmed with the help of typical test functions for a second order accurate implementation of the scheme. The accuracy of the gradient evaluation on arbitrarily skewed finite volume cells has been demonstrated using the order of convergence of solution to Poisson’s equation and viscous flow in a driven square cavity. A straight forward extension of the scheme for three dimensional problems has been illustrated by computing three dimensional potential flow over a sphere in Cartesian coordinates and viscous flow in a cubic cavity. A new interpretation to Coirier diamond path integration presented here shows that the diamond path can be integrated in a compact manner identical to that of a typical implementation of the present approach.

Transient Method for Convective Heat Transfer Measurement With Lateral Conduction—Part II: Application to an Isolated Spherical Roughness Element

The effect of lateral conduction on convective heat transfer measurements using a transient infrared technique over an isolated spherical roughness element (bump) is evaluated. Comparisons are made between a full 3D finite-volume analysis and a simpler 1D transient conduction model. The surface temperature history was measured with a high resolution infrared camera during an impulsively started hot-gas flow at a flow Reynolds number of 860,000. The boundary layer was turbulent with the bump heights equivalent to 0.75, 1.5, and 3 times the boundary layer momentum thickness. When considering transient conduction effects only in the bump wake, the 1D approximate method underestimates the actual Stanton number estimated with the 3D model. This discrepancy is only 10% for a 75% change in St number occurring over a surface distance of 10 mm (the half-width of the wake). When the actual bump topology is accounted for in estimating the Stanton number on the bump itself with the 3D analysis technique, the increased surface area of the finite-volume cells on the protruding bump actually decreases the predicted value of St locally. The net result is that the two effects can cancel each other, and in some cases the 1D approximate technique can provide a reasonably accurate estimate of the surface heat transfer without the added complexity of the 3D finite-volume method. For the case of the largest bump tested, with maximum surface angularity exceeding 60 deg, the correction for 3D topology yields a 1D St estimate that is within 20–30% of the 3D estimate over much of the bump surface. These observed effects are valid for transient measurement techniques while the opposite is true for steady-state measurement techniques.

FINITE VOLUME METHODS FOR MULTI-SYMPLECTIC PDES ∗

We investigate the application of a cell-vertex finite volume discretization to multi-symplectic PDEs. The investigated discretization reduces to the Preissman box scheme when used on a rectangular grid. Concerning arbitrary quadrilateral grids, we show that only methods with parallelogram-like finite volume cells lead to a multi-symplectic discretization; i.e., to a method that preserves a discrete conservation law of symplecticity. One of the advantages of finite volume methods is that they can be easily adjusted to variable meshes. But, although the implementation of moving mesh finite volume methods for multi-symplectic PDEs is rather straightforward, the restriction to parallelogram-like cells implies that only meshes moving with a constant speed are multi-symplectic. To overcome this restriction, we suggest the implementation of reversible moving mesh methods based on a semi-Lagrangian approach. Numerical experiments are presented for a one dimensional dispersive shallow-water system.

A Finite Volume Extension of the Lax-Friedrichs and Nessyahu-Tadmor Schemes for Conservation Laws on Unstructured Grids

The non-oscillatory central difference scheme of Nessyahu and Tadmor, in which the resolution or Riemann problems at the cell interfaces is by-passed thanks to the use of the staggered Lax-Friedrichs scheme, is extended here to a two-step, two-dimensional non-oscillatory centered scheme in finite volume formulation. The construction of the scheme rests on a finite volume extension of the Lax-Friedrichs scheme, in which the finite volume cells are the barycentric cells constructed around the nodes of an FEM triangulation, for odd time steps, and some quadrilateral cells associated with this triangulation, for even time steps. Piecewise linear cell interpolants using least-squares gradients combined with a van Leer-type slope limiting allow for an oscillation-free second-order resolution. Some preliminary numerical experiments suggest that two-dimensional problems can be handled very efficiently by the method presented here.

An implicit mixed finite-volume-finite-element method for solving 3D turbulent compressible flows

The development of new aeronautic projects require accurate and efficient simulations of compressible flows in complex geometries. It is well known that most flows of interest are at least locally turbulent and that the modelling of this turbulence is critical for the reliability of the computations. A turbulence closure model which is both cheap and reasonably accurate is an essential part of a compressible code. An implicit algorithm to solve the 2D and 3D compressible Navier-Stokes equations on unstructured triangular/tetrahedral grids has been extended to turbulent flows. This numerical scheme is based on second-order finite element-finite volume discretization: the diffusive and source terms of the Navier-Stokes equations are computed using a finite element method, while the other terms are computed with a finite volume method. Finite volume cells are built around each node by means of the medians. The convective fluxes are evaluated with the approximate Riemann solver of Roe coupled with the van Albada limiter. The standard κ-e model has been introduced to take into account turbulence. Implicit integration schemes with efficient numerical methods (CGS, GMRES and various preconditioning techniques) have also been implemented

Three-Dimensional Coupled Fluid-Structure Simulation of Pericardial Bioprosthetic Aortic Valve Function

A computational, three-dimensional coupled fluid-structure dynamics model was developed for a generic pericardial aortic valve in a rigid aortic root graft with physiologic sinuses. Valve geometry was based on that of the natural valve. Blood flow was modeled as pulsatile, laminar, Newtonian, incompressible flow. The structural model accounted for material and geometric nonlinearities and also simulated leaflet coaptation. A body fitted grid was used to subdivide the flow domain into computational finite volume cells. Shell finite elements were used to discretize the leaflet volume. A finite volume computational fluid dynamics code and finite element structure dynamics code were used to solve the flow and structure equations, respectively. The fluid flow and structural equations were coupled using an implicit "influence coefficient" technique. Physiologic ventricular and aortic pressure waveforms were prescribed as the flow boundary conditions. The aortic flow field, valve structural configuration, and leaflet stresses were computed at 2 msec intervals. Model predictions on aortic flow and transient variation in valve orifice area were in close agreement with corresponding experimental in vitro data. These findings suggest that the computer model has potential for being a powerful design tool for bioprosthetic aortic valves.