Uniform Cartesian Grid Research Articles

Transport modeling of sedimenting particles in a turbulent pipe flow using Euler–Lagrange large eddy simulation

A volume-filtered Euler–Lagrange large eddy simulation methodology is used to predict the physics of turbulent liquid–solid slurry flow through a horizontal periodic pipe. A dynamic Smagorinsky model based on Lagrangian averaging is employed to account for the sub-filter scale effects in the liquid phase. A fully conservative immersed boundary method is used to account for the pipe geometry on a uniform cartesian grid. The liquid and solid phases are coupled through volume fraction and momentum exchange terms. Particle–particle and particle–wall collisions are modeled using a soft-sphere approach. Three simulations are performed by varying the superficial liquid velocity to be consistent with the experimental data by Dahl et al. (2003). Depending on the liquid flow rate, a particle bed can form and develop different patterns, which are discussed in light of regime diagrams proposed in the literature. The fluctuation in the height of the liquid-bed interface is characterized to understand the space and time evolution of these patterns. Statistics of engineering interest such as mean velocity, mean concentration, and mean streamwise pressure gradient driving the flow are extracted from the numerical simulations and presented. Sand hold-up calculated from the simulation results suggest that this computational strategy is capable of predicting critical deposition velocity.

Strongly coupled dynamics of fluids and rigid-body systems with the immersed boundary projection method

A strong coupling algorithm is presented for simulating the dynamic interactions between incompressible viscous flows and rigid-body systems in both two- and three-dimensional problems. In this work, the Navier–Stokes equations for incompressible flow are solved on a uniform Cartesian grid by the vorticity-based immersed boundary projection method of Colonius and Taira. Dynamical equations for arbitrary rigid-body systems are also developed. The proposed coupling method attempts to unify the treatment of constraints in the fluid and structure—the incompressibility of the fluid, the linkages in the rigid-body system, and the conditions at the interface—through the use of Lagrange multipliers. The resulting partitioned system of equations is solved with a simple relaxation scheme, based on an identification of virtual inertia from the fluid. The scheme achieves convergence in only 2 to 5 iterations per time step for a wide variety of mass ratios. The formulation requires that only a subset of the discrete fluid equations be solved in each iteration. Several two- and three-dimensional numerical tests are conducted to validate and demonstrate the method, including a falling cylinder, flapping of flexible wings, self-excited oscillations of a system of many linked plates in a free stream, and passive pivoting of a finite aspect ratio plate under the influence of gravity in a free stream. The results from the current method are compared with previous experimental and numerical results and good agreement is achieved.

Entropy stable discontinuous interfaces coupling for the three-dimensional compressible Navier–Stokes equations

Non-linear entropy stability and a summation-by-parts (SBP) framework are used to derive entropy stable interior interface coupling for the semi-discretized three-dimensional (3D) compressible Navier–Stokes equations. A complete semidiscrete entropy estimate for the interior domain is achieved combining a discontinuous entropy conservative operator of any order [1,2] with an entropy stable coupling condition for the inviscid terms, and a local discontinuous Galerkin (LDG) approach with an interior penalty (IP) procedure for the viscous terms. The viscous penalty contributions scale with the inverse of the Reynolds number (Re) so that for Re → ∞ their contributions vanish and only the entropy stable inviscid interface penalty term is recovered. This paper extends the interface couplings presented [1,2] and provides a simple and automatic way to compute the magnitude of the viscous IP term. The approach presented herein is compatible with any diagonal norm summation-by-parts (SBP) spatial operator, including finite element, finite volume, finite difference schemes and the class of high-order accurate methods which include the large family of discontinuous Galerkin discretizations and flux reconstruction schemes. This note relies on the formalism introduced in [1,3] and complements the new class of interior entropy stable SBP operators of any order for the 3D compressible Navier–Stokes equations on unstructured grids that was proposed in [1,2]. To keep the notation as simple as possible, a uniform Cartesian grid is considered in the derivation. However, the extension to generalized curvilinear coordinates and unstructured grids follows immediately if the transformation from computational to physical space preserves the semi-discrete geometric conservation [4]. The proposed interface coupling technique has been successfully combined with a high order entropy stable discretization for the simulation of two-dimensional (2D) and 3D viscous subsonic and supersonic flows presented in [1,5,6].

The Closest Point Method and Multigrid Solvers for Elliptic Equations on Surfaces

Elliptic partial differential equations are important from both application and analysis points of view. In this paper we apply the closest point method to solve elliptic equations on general curved surfaces. Based on the closest point representation of the underlying surface, we formulate an embedding equation for the surface elliptic problem, then discretize it using standard finite differences and interpolation schemes on banded but uniform Cartesian grids. We prove the convergence of the difference scheme for the Poisson's equation on a smooth closed curve. In order to solve the resulting large sparse linear systems, we propose a specific geometric multigrid method in the setting of the closest point method. Convergence studies in both the accuracy of the difference scheme and the speed of the multigrid algorithm show that our approaches are effective.

Simulation of forced deformable bodies interacting with two-dimensional incompressible flows: Application to fish-like swimming

We present an efficient algorithm for simulation of deformable bodies interacting with two-dimensional incompressible fluid flows. The temporal and spatial discretizations of the Navier–Stokes equations in vorticity stream-function formulation are based on classical fourth-order Runge–Kutta scheme and compact finite differences, respectively. Using a uniform Cartesian grid we benefit from the advantage of a new fourth-order direct solver for the Poisson equation to ensure the incompressibility constraint down to machine zero over an optimal grid. For introducing a deformable body in fluid flow, the volume penalization method is used. A Lagrangian structured grid with prescribed motion covers the deformable body which is interacting with the surrounding fluid due to the hydrodynamic forces and the torque calculated on the Eulerian reference grid. An efficient law for controlling the curvature of an anguilliform fish, swimming toward a prescribed goal, is proposed which is based on the geometrically exact theory of nonlinear beams and quaternions. Validation of the developed method shows the efficiency and expected accuracy of the algorithm for fish-like swimming and also for a variety of fluid/solid interaction problems.

SPGrid

We introduce a new method for fluid simulation on high-resolution adaptive grids which rivals the throughput and parallelism potential of methods based on uniform grids. Our enabling contribution is SPGrid , a new data structure for compact storage and efficient stream processing of sparsely populated uniform Cartesian grids. SPGrid leverages the extensive hardware acceleration mechanisms inherent in the x86 Virtual Memory Management system to deliver sequential and stencil access bandwidth comparable to dense uniform grids. Second, we eschew tree-based adaptive data structures in favor of storing simulation variables in a pyramid of sparsely populated uniform grids, thus avoiding the cost of indirect memory access associated with pointer-based representations. We show how the costliest algorithmic kernels of fluid simulation can be implemented as a composition of two kernel types: (a) stencil operations on a single sparse uniform grid, and (b) structured data transfers between adjacent levels of resolution, even when modeling non-graded octrees. Finally, we demonstrate an adaptive multigrid-preconditioned Conjugate Gradient solver that achieves resolution-independent convergence rates while admitting a lightweight implementation with a modest memory footprint. Our method is complemented by a new interpolation scheme that reduces dissipative effects and simplifies dynamic grid adaptation. We demonstrate the efficacy of our method in end-to-end simulations of smoke flow.

A high-order Cartesian-grid finite-volume method for aeroacoustics simulations

A moving-least-square based finite-volume method is developed to simulate acoustic wave propagation and scattering from complicated solid geometries. This hybrid method solves the linearized perturbed compressible equations as the governing equations of the acoustic field. The solid boundaries are embedded in a uniform Cartesian grid and represented using level set fields. Thus, the current approach avoids unstructured grid generation for the irregular geometries. The desired boundary conditions are imposed sharply on the immersed boundaries using a ghost fluid method. The scope of the implementation of the moving moving-least-square approach in the current solver is threefold: reconstruction of the field variables on cell faces for high-order flux construction, population of the ghost cells based on the desired boundary condition, and filtering the high wave number modes near the immersed boundaries. The computational stencils away from the boundaries are identical; hence, only one moving-least-square shape-function is computed and stored with its underlying grid pattern for all the interior cells. This feature significantly reduces the memory requirement of the acoustic solver compared to similar finite-volume method on irregular unstructured mesh. The acoustic solver is validated against several benchmark problems.

A well-balanced and asymptotic-preserving scheme for the one-dimensional linear Dirac equation

The numerical approximation of one-dimensional relativistic Dirac wave equations is considered within the recent framework consisting in deriving local scattering matrices at each interface of the uniform Cartesian computational grid. For a Courant number equal to unity, it is rigorously shown that such a discretization preserves exactly the $$L^2$$ norm despite being explicit in time. This construction is well-suited for particles for which the reference velocity is of the order of $$c$$ , the speed of light. Moreover, when $$c$$ diverges, that is to say, for slow particles (the characteristic scale of the motion is non-relativistic), Dirac equations are naturally written so as to let a “diffusive limit” emerge numerically, like for discrete 2-velocity kinetic models. It is shown that an asymptotic-preserving scheme can be deduced from the aforementioned well-balanced one, with the following properties: it yields unconditionally a classical Schrodinger equation for free particles, but it handles the more intricate case with an external potential only conditionally (the grid should be such that $$c \Delta x\rightarrow 0$$ ). Such a stringent restriction on the computational grid can be circumvented easily in order to derive a seemingly original Schrodinger scheme still containing tiny relativistic features. Numerical tests (on both linear and nonlinear equations) are displayed.

A Fourth Order Finite Difference Method for Numerical Solution of the Goursat Problem

In this article, we consider the finite difference method for numerically solving the Goursat Problem, using uniform Cartesian grids on the square region. Numerical examples are considered to ensure accuracy of the developed method on both linear and nonlinear Goursat Problems of second order partial differential equations. The results obtained for these numerical examples validate the efficiency and expected fourth order accuracy of the method.

A Study of Three-Dimensional Acoustic Scattering by Arbitrary Distribution Multibodies Using Extended Immersed Boundary Method

A three-dimensional computational model for acoustic scattering with complex geometries is presented, which employs the immersed boundary technique to deal with the effect of both hard and soft wall boundary conditions on the acoustic fields. In this numerical model, the acoustic field is solved on uniform Cartesian grids, together with simple triangle meshes to partition the immersed body surface. A direct force at the Lagrangian points is calculated from an influence matrix system, and then spreads to the neighboring Cartesian grid points to make the acoustic field satisfy the required boundary condition. This method applies a uniform stencil on the whole domain except at the computational boundary, which has the benefit of low dispersion and dissipation errors of the used scheme. The method has been used to simulate two benchmark problems to validate its effectiveness and good agreements with the analytical solutions are achieved. No matter how complex the geometries are, single body or multibodies, complex geometries do not pose any difficulty in this model. Furthermore, a simple implementation of time-domain impedance boundary condition is reported and demonstrates the versatility of the computational model.

Numerical study of single bubble motion in liquid metal exposed to a longitudinal magnetic field

The paper presents numerical simulations modeling the ascent of an argon bubble in liquid metal with and without an external magnetic field. The governing equations for the fluid and the electric potential are discretized on a uniform Cartesian grid and the bubble is represented with a highly efficient immersed boundary method. The simulations performed were conducted matching experiments under the same conditions so that sound validation is possible. The three-dimensional trajectory of the bubble is analyzed quantitatively and related to the flow structures in the wake. Indeed, the substantial impact of the magnetic field in the bubble trajectory results from its influence on the wake. Quantitative data describing the selective damping of vortex structures are provided and discussed. As a result of applying a longitudinal field, the time-averaged bubble rise velocity increases for large bubbles, it reaches a maximum and then decreases when further increasing the magnetic interaction parameter. For small bubbles, the time-averaged bubble rise velocity decreases when increasing the magnetic field. The bubble Strouhal number as a dimensionless frequency is reduced with the application of a magnetic field for all bubbles considered and the zig–zag trajectory of the bubble becomes more rectilinear.

A Finite Difference Method for Numerical Solution of Goursat Problem of Partial Differential Equation

In this article, we report the finite difference method for numerically solving the Goursat Problem, using uniform Cartesian grids on the square region. We have considered both linear and nonlinear Goursat problems of partial differential equations for the numerical solution, to ensure the accuracy of the developed method. The results obtained for these numerical examples validate the efficiency, expected order and accuracy of the method.

Fast and Accurate Evaluation of Nonlocal Coulomb and Dipole-Dipole Interactions via the Nonuniform FFT

We present a fast and accurate algorithm for the evaluation of nonlocal (long-range) Coulomb and dipole-dipole interactions in free space. The governing potential is simply the convolution of an interaction kernel $U(\mathbf{x})$ and a density function $\rho(\mathbf{x})=|\psi(\mathbf{x})|^2$ for some complex-valued wave function $\psi(\mathbf{x})$, permitting the formal use of Fourier methods. These are hampered by the fact that the Fourier transform of the interaction kernel $\widehat{U}(\mathbf{k})$ has a singularity and/or $\widehat{\rho}(\mathbf{k})\ne0$ at the origin $\mathbf{k}={\bf 0}$ in Fourier (phase) space. Thus, accuracy is lost when using a uniform Cartesian grid in $\mathbf{k}$ which would otherwise permit the use of the FFT for evaluating the convolution. Here, we make use of a high-order discretization of the Fourier integral, accelerated by the nonuniform fast Fourier transform (NUFFT). By adopting spherical and polar phase-space discretizations in three and two dimensions, respectively, the singularity in $\hat{U}(\mathbf{k})$ at the origin is canceled so that only a modest number of degrees of freedom are required to evaluate the Fourier integral, assuming that the density function $(\rho\mathbf{x})$ is smooth and decays sufficiently quickly as $|\mathbf{x}| \rightarrow \infty$. More precisely, the calculation requires $O(N\log N)$ operations, where $N$ is the total number of discretization points in the computational domain. Numerical examples are presented to demonstrate the performance of the algorithm.

A Two-Dimensional Version of the Godunov Scheme for Scalar Balance Laws

A Godunov scheme is derived for two-dimensional scalar conservation laws without or with source terms following ideas originally proposed by Boukadida and LeRoux [Math. Comput., 63 (1994), pp. 541--553] in the context of a staggered Lax--Friedrichs scheme. In both situations, the numerical fluxes are obtained at each interface of a uniform Cartesian computational grid just by means of the “external waves” involved in the entropy solution of the elementary two-dimensional (2D) Riemann problems; in particular, all the wave-interaction phenomena are overlooked. This restriction of the wave pattern suffices for deriving the exact numerical fluxes of the staggered Lax--Friedrichs scheme, but it furnishes only an approximation for the Godunov scheme: we show that under convenient assumptions, these flux functions are smooth and the resulting discretization process is stable under nearly optimal CFL restriction. A well-balanced extension is presented, relying on the Curl-free component of the Helmholtz decomposi...

VIE-FG-FFT for Analyzing EM Scattering from Inhomogeneous Nonmagnetic Dielectric Objects

A new realization of the volume integral equation (VIE) in combination with the fast Fourier transform (FFT) is established by fitting Green’s function (FG) onto the nodes of a uniform Cartesian grid for analyzing EM scattering from inhomogeneous nonmagnetic dielectric objects. The accuracy of the proposed method is the same as that of the P-FFT and higher than that of the AIM and the IE-FFT especially when increasing the grid spacing size. Besides, the preprocessing time of the proposed method is obviously less than that of the P-FFT for inhomogeneous nonmagnetic dielectric objects. Numerical examples are provided to demonstrate the accuracy and efficiency of the proposed method.

Time asynchronous relative dimension in space method for multi-scale problems in fluid dynamics

A novel computational method is presented for solving fluid dynamics equations in the multi-scale framework when the system size is an important parameter of the governing equations. The method (TARDIS) is based on a concurrent transformation of the governing equations in space and time and solving the transformed equations on a uniform Cartesian grid with the corresponding causality conditions at the grid interfaces. For implementation in the framework of TARDIS, the second-order CABARET scheme of Karabasov and Goloviznin [1] is selected for it provides a good combination of numerical accuracy, computational efficiency and simplicity of realisation. Numerical examples are first provided for several isothermal gas dynamics test problems and then for modelling of molecular fluctuations inside a microscopic flow channel and ultrasound wave propagation through a nano-scale region of molecular fluctuations.

Combining triangle Gaussian integration and modified NUFFT for evaluating two-dimensional Fourier transform integrals

The regular fast Fourier transform (FFT) requires a uniform Cartesian orthogonal grid which has considerable stair-casing errors when dealing with the function having an arbitrary shape boundary. The recently proposed two-dimensional discontinuous fast Fourier transform (2D-DFFT) can overcome this problem by using triangle mesh discretization and Gaussian numerical integration. However, the interpolation is used for the function data in the original 2D-DFFT, which reduces the accuracy performance especially for the case of oscillating functions. This work presents a useful modification of the original 2D-DFFT by removing the requirement of function interpolation to obtain significant accuracy improvement. In addition, the modified 2D nonuniform fast Fourier transform (NUFFT) with real-valued least-square interpolation coefficients are developed to speed up the computation of numerical Fourier transform over the triangle mesh. Numerical experiments are conducted to demonstrate the effectiveness and advantages of the proposed algorithms.

A Highly Accurate FGG-FG-FFT for the Combined Field Integral Equation

In this paper, a novel realization of the Integral Equation in combination with the fast Fourier transform for the CFIE is established by Fitting both the Green's function and its Gradient onto the nodes of a uniform Cartesian grid. The new method has been compared with several existing popular FFT-based methods, including the AIM, the IE-FFT, and the p-FFT. The accuracy of the proposed method is significantly higher than other FFT-based methods, and the method is not sensitive to both the grid spacing and the expansion order. The outstanding merit of the proposed method is that the fitting procedure is independent of the basis functions. Therefore, when the higher order basis functions would be adopted in the method of moments, only one fitting procedure for the Green's function and its gradient on a basis function support is needed to meet all of basis functions defined on this support. Some numerical examples are provided in this paper to demonstrate the accuracy and efficiency of the proposed method.

A particle-tracking approach for accurate material derivative measurements with tomographic PIV

The evaluation of the instantaneous 3D pressure field from tomographic PIV data relies on the accurate estimate of the fluid velocity material derivative, i.e., the velocity time rate of change following a given fluid element. To date, techniques that reconstruct the fluid parcel trajectory from a time sequence of 3D velocity fields obtained with Tomo-PIV have already been introduced. However, an accurate evaluation of the fluid element acceleration requires trajectory reconstruction over a relatively long observation time, which reduces random errors. On the other hand, simple integration and finite difference techniques suffer from increasing truncation errors when complex trajectories need to be reconstructed over a long time interval. In principle, particle-tracking velocimetry techniques (3D-PTV) enable the accurate reconstruction of single particle trajectories over a long observation time. Nevertheless, PTV can be reliably performed only at limited particle image number density due to errors caused by overlapping particles. The particle image density can be substantially increased by use of tomographic PIV. In the present study, a technique to combine the higher information density of tomographic PIV and the accurate trajectory reconstruction of PTV is proposed (Tomo-3D-PTV). The particle-tracking algorithm is applied to the tracers detected in the 3D domain obtained by tomographic reconstruction. The 3D particle information is highly sparse and intersection of trajectories is virtually impossible. As a result, ambiguities in the particle path identification over subsequent recordings are easily avoided. Polynomial fitting functions are introduced that describe the particle position in time with sequences based on several recordings, leading to the reduction in truncation errors for complex trajectories. Moreover, the polynomial regression approach provides a reduction in the random errors due to the particle position measurement. Finally, the acceleration can be evaluated analytically, which greatly reduces the truncation errors due to finite differences. The approach is first assessed using computer-generated data of an advecting vortex ring. Precision errors in the material derivative can be reduced with a factor 2–3. This is achieved when a long sequence is considered (e.g. 15–20 recordings). Similarly, truncation errors typically associated with direct integration and finite differences from the PIV-based technique are almost eliminated. It is shown that the material derivative information obtained at the scattered locations in the 3D domain can be reduced to a uniform Cartesian grid by means of a second-order spatial regression with no significant artefact. The technique is applied to a benchmark Tomo-PIV experiment of a transitional jet in water. The results confirm the conclusions obtained with the numerical study. Moreover, it is shown that the evaluation of the instantaneous 3D pressure field can be retrieved with significant reduction in artefacts associated with random and truncation errors.

A two-grid adaptive volume of fluid approach for dendritic solidification

A two-grid approach to simulating dendritic growth in the solidification of pure metals is presented. The approach is based on using a uniform Cartesian grid to solve the energy equation and an adaptive Cartesian grid of higher-resolution to solve the interface evolution, providing a more accurate representation of the interface. The energy equation is solved using a diffuse-interface method, and the advection equation of a discretized solid fraction function is solved using an unsplit VOF (volume of fluid) method, with PLIC (piecewise linear interface calculation) reconstruction of the interface. The thermal gradients at both sides of the interface, which are needed to obtain the front velocity, are calculated with the aid of the distance function to the reconstructed interface, which is obtained using an efficient method described in detail in this work. The influence of the grid resolution used to solve the advection equation on the accuracy of the method is analyzed. For the dendritic growth cases considered in this work, particularly, it was found that using a grid resolution for the advection equation two times higher than that used for the energy equation considerably improves the results, while keeping the CPU time consumed at an acceptable level. To underline the importance of an accurate estimation of the interface curvature, the results obtained using three techniques commonly used in VOF methods of different degrees of accuracy (a convolved VOF technique, a distance function technique and an improved height function technique, which provides second-order accuracy) are compared. The proposed methodology is assessed by comparing the numerical results with analytical solutions and with results obtained by different authors for the formation of complex dendritic structures in two and three dimensions.