Divergence-free Condition Research Articles

A semi-implicit meshless method for incompressible flows in complex geometries

We present an exponentially convergent semi-implicit meshless algorithm for the solution of Navier-Stokes equations in complex domains. The algorithm discretizes partial derivatives at scattered points using radial basis functions (RBF) as interpolants. Higher-order polynomials are appended to the polyharmonic splines (PHS-RBF) and a collocation method is used to derive the interpolation coefficients. The interpolating kernels are then differentiated and the partial-differential equations are satisfied by collocation at the scattered points. The PHS-RBF interpolation is shown to be exponentially convergent with discretization errors decreasing as a high power of a representative distance between points. We present here a semi-implicit algorithm for time-dependent and steady state fluid flows in complex domains. At each time step, several iterations are performed to converge the momentum and continuity equations. A Poisson equation for pressure corrections is formulated by imposing divergence free condition on the iterated velocity field. At each time step, the momentum and pressure correction equations are repeatedly solved until the velocities and pressure converge to a pre-specified tolerance. We have demonstrated the convergence and discretization accuracy of the algorithm for two model problems and simulated four other complex problems. In all cases, the algorithm is stable for Courant numbers in excess of ten. The algorithm has the potential to accurately and efficiently solve many fluid flow and heat transfer problems in complex domains. An open source code Meshless Multi-Physics Software (MeMPhyS) [1] is available for interested users of the algorithm.

An improved consistent inflow turbulence generator for LES evaluation of wind effects on buildings

This paper proposes an improved consistent inflow turbulence generator for reproducing spatially correlated turbulence fields with specified wind profiles and spectra for large eddy simulation (LES). The inflow turbulence is developed based on the widely used Random Flow Generation (RFG) technique, and improvements include obtaining more consistent wind spectra, direct control of spatial integral length scales, mitigation of fictitious pressure fluctuations caused by nonsatisfaction of divergence-free conditions, and iteration-based adaption of velocity profiles within the future building area. These improvements are obtained from strictly theoretical derivation, and the correctness of the proposed algorithm is verified by comparing the velocity profiles and the wind spectra of the simulated turbulence atmospheric boundary layer (ABL) with target values. The proposed turbulence generator is adopted to generate inflow turbulence for analyzing wind-induced pressure distributions on isolated low-rise buildings immersed in the turbulent ABL, and the comparison between simulated results and wind tunnel experimental results obtained by Tokyo Polytechnic University (TPU) is presented. Good agreement between the numerical and experimental results verifies the applicability of the proposed generator technique to LES of wind effects on low-rise buildings.

An Extension of Gmunu: General-relativistic Resistive Magnetohydrodynamics Based on Staggered-meshed Constrained Transport with Elliptic Cleaning

We present the implementation of general-relativistic resistive magnetohydrodynamics solvers and three divergence-free handling approaches adopted in the General-relativistic multigrid numerical (Gmunu) code. In particular, implicit–explicit Runge–Kutta schemes are used to deal with the stiff terms in the evolution equations for small resistivity. The three divergence-free handling methods are (i) hyperbolic divergence cleaning (also known as the generalized Lagrange multiplier), (ii) staggered-meshed constrained transport schemes, and (iii) elliptic cleaning through a multigrid solver, which is applicable in both cell-centered and face-centered (stagger grid) magnetic fields. The implementation has been tested with a number of numerical benchmarks from special-relativistic to general-relativistic cases. We demonstrate that our code can robustly recover from the ideal magnetohydrodynamics limit to a highly resistive limit. We also illustrate the applications in modeling magnetized neutron stars, and compare how different divergence-free handling methods affect the evolution of the stars. Furthermore, we show that the preservation of the divergence-free condition of the magnetic field when using staggered-meshed constrained transport schemes can be significantly improved by applying elliptic cleaning.

Enhanced dissipation for the third component of 3D anisotropic Navier-Stokes equations

In this paper, we study the decay rates for the global small smooth solutions to 3D anisotropic incompressible Navier-Stokes equations. In particular, we prove that the horizontal components of the velocity field decay like the solutions of 2D classical Navier-Stokes equations. While the third component of the velocity field decays as the solutions of 3D Navier-Stokes equations. We remark that such enhanced decay rate for the third component is caused by the interplay between the divergence free condition of the velocity field and the horizontal Laplacian in the anisotropic Navier-Stokes equations.

Extraction Of Temperature Fields From PIV Data Of Turbulent Rayleigh-Bénard Convection Using DNS

To gain insights into the dynamical processes that govern the problem of turbulent Rayleigh-Bénard convection, information about the temperature field---next to the velocity field---is required. In the present study, we exploit the governing flow equations to extract the temperature field from a velocity field measured via tomographic particle image velocimetry. The methods used within this study are derived from algorithms originally used in the framework of direct numerical simulations. First, we enforce the divergence-free condition on the measured velocity field using a fractional-step method based on a Poisson equation. Second, we extract the temperature field from the velocity field via the momentum equation. Third, we extract the temperature field from the velocity field via the solution of the stationary convection-diffusion equation. The extracted velocity field---together with the corresponding divergence-free velocity field---can serve as an initial condition for a direct numerical simulation of the problem. Validating the extraction method based on the stationary convection-diffusion equation with data from direct numerical simulations shows a good correlation coefficient between the extracted and the simulated temperature field (0.86). Although the third method incorporates the ambitious assumption of stationarity, temperature fields are successfully extracted from the measured velocity fields. However, with a correlation coefficient of 0.57, the agreement between the average extracted temperature field and the average temperature field measured in a subvolume of the flow domain is somewhat less than the one found in the simulated data. The latter is presumably attributed not only to the extraction method, but also to the underlying measurement data. Thus, an increase in measurement accuracy in future applications is likely to improve the accuracy of the extraction methods as well.

A Three-Order, Divergence-Free Scheme for the Simulation of Solar Wind

In this paper, we present a three-order, divergence-free finite volume scheme to simulate the steady state solar wind ambient. The divergence-free condition of the magnetic field is preserved by the constrained transport (CT) method. The CT method can keep the magnetic fields divergence free if the magnetic fields is divergence free initially. Thus, a least-squares reconstruction of magnetic field with the divergence free constraints is used to make the magnetic fields global solenoidality initially. High order spatial accuracy is obtained through a non-oscillatory hierarchical reconstruction, while a high order time discretization is achieved using a three-order Runge–Kutta scheme. This new model of three order in space and time is validated by numerical results for Carrington rotation 2207. The numerical results show its capability for producing stable reliable results for structured solar wind. The high-order, divergence-free properties of this method make it an ideal tool for the simulations of coronal mass ejection in future.

An Electric-field-driven Global Coronal Magnetohydrodynamics Simulation Model Using Helioseismic and Magnetic Imager Vector-magnetic-field Synoptic Map Data

We present the simulation methodology and results of our new data-driven global coronal magnetohydrodynamics (MHD) simulation model. In this model, the solar-surface electric field is first calculated such that the curl will satisfy both the induction equation and the given temporal variations of the solar-surface magnetic field. We use the synoptic maps of the Helioseismic and Magnetic Imager three-component vector-magnetic-field data to specify the solar-surface magnetic-field vector for a period from Carrington Rotations (CRs) 2106 to 2110. A set of whole-Sun three-component electric-field maps are obtained for each CR transition interval of about 27.3 days. Using the inverted electric field as the driving variable, our new global coronal MHD model, with the angular resolution of π/64, can trace the evolution of the three-dimensional coronal magnetic field that matches the specified time-dependent solar-surface magnetic-field maps and simultaneously satisfies the divergence-free condition. A set of additional boundary treatments are introduced to control the contribution of the horizontal components of the magnetic field at the weak-field regions. The strength of the solar-surface magnetic field is limited to 20 Gauss for the sake of computational stability in this study. With these numerical treatments, the nonpotential coronal features, such as twisted loop structures, and their eruptive outward motions are obtained. This present model, capable of introducing three-component solar-surface magnetic-field observation data to coronal MHD simulations, is our first step toward a better model framework for the solar corona and hence solar wind.

High-order curvilinear finite element magneto-hydrodynamics I: A conservative Lagrangian scheme

Magneto-hydrodynamics is one of the foremost models in plasma physics with applications in inertial confinement fusion, astrophysics and elsewhere. Advanced numerical methods are needed to get an insight into the complex physical phenomena. The classical Lagrangian methods are typically limited to the low orders of convergence and suffer from violation of the divergence-free condition for magnetic field or conservation of the invariants. This paper is the first part of a new series about high-order non-ideal magneto-hydrodynamics, where a multi-dimensional conservative Lagrangian method based on curvilinear finite elements is presented. The condition on zero divergence of magnetic field and conservation of mass, momentum, magnetic flux and the total energy are satisfied exactly. The curvilinear elements prevent entangling of the computational mesh and its imprinting into the solution. A high-order conservative time integration is applied, where an arbitrary order of convergence is attained for problems of ideal magneto-hydrodynamics. The resistive magnetic field diffusion is solved by an implicit scheme. Description of the method is given and multiple test problems demonstrating properties of the scheme are performed. The construction of the method and possible future directions of development are discussed.

A physics-inspired alternative to spatial filtering for large-eddy simulations of turbulent flows

Large-eddy simulations (LES) are widely used for computing high Reynolds number turbulent flows. Spatial filtering theory for LES is not without its shortcomings, including how to define filtering for wall-bounded flows, commutation errors for non-uniform filters and extensibility to flows with additional complexity, such as multiphase flows. In this paper, the theory for LES is reimagined using a coarsening procedure that imitates nature. This physics-inspired coarsening approach is equivalent to Gaussian filtering for single-phase wall-free flows but opens up new insights for both physical understanding and modelling even in that simple case. For example, an alternative to the Germano identity is introduced and used to define a dynamic procedure without the need for a test filter. Non-uniform resolution can be represented in this framework without commutation errors, and the divergence-free condition is retained for incompressible flows. Potential extensions of the theory to more complex physics such as multiphase flows are briefly discussed.

A Semi-Implicit Meshless Method for Incompressible Flows in Complex Geometries

We present an exponentially convergent semi-implicit meshless algorithm for the solution of Navier-Stokes equations in complex domains. The algorithm discretizes partial derivatives at scattered points using radial basis functions as interpolants. Higher-order polynomials are appended to polyharmonic splines (PHS-RBF) and a collocation method is used to derive the interpolation coefficients. The interpolating kernels are then differentiated and the partial-differential equations are satisfied by collocation at the scattered points. The PHS-RBF interpolation is shown to be exponentially convergent with discretization errors decreasing as a high power of a representative distance between points. We present here a semi-implicit algorithm for time-dependent and steady state fluid flows in complex domains. At each time step, several iterations are performed to converge the momentum and continuity equations. A Poisson equation for pressure corrections is formulated by imposing divergence free condition on the iterated velocity field. At each time step, the momentum and pressure correction equations are repeatedly solved until the velocities and pressure converge to a pre-specified tolerance. We have demonstrated the convergence and discretization accuracy of the algorithm for two model problems and simulated three other complex problems. In all cases, the algorithm is stable for Courant numbers in excess of ten. The algorithm has the potential to accurately and efficiently solve many fluid flow and heat transfer problems in complex domains. An open source code Meshless Multi-Physics Software (MeMPhyS) is available for interested users of the algorithm.

Estimates of the Distance to the Solution of an Evolutionary Problem Obtained by Linearization of the Navier–Stokes Equation

The paper is concerned with linearization of the Navier–Stokes equation in the space-time cylinder QT . The main goal is to deduce computable estimates of the distance between the exact solution and a function in the energy admissible class of vector-valued functions. First, the estimates are derived for the case, where this class contains only divergence free (solenoidal) functions. In the next section, estimates of the distance to sets of divergence free functions depending on the space and time variables are considered. These results are used to extend earlier derived estimates to non-solenoidal approximations. The corresponding estimates contain an additional term which can be viewed as a penalty for the violation of the divergence free condition.

Development of explicit moving particle simulation method with applications

Moving particle simulation (MPS) is a Lagrange meshless method capable of simulating complex fluid problems. To improve computational efficiency with a simple algorithm, explicit algorithms have been proposed such as the Weakly Compressible MPS (WCMPS) method. Although accuracy can be fulfilled in WCMPS, instability arises due to pressure fluctuations of the moving particles. In this study, an explicit Poisson Pressure Equation (PPE) was introduced to address particle instability by improving the pressure calculation, referred to as the Weakly Compressible Moving Particle Explicit (WCMPE) method. By predicting pressure through the equation of state for WCMPS, divergence-free condition is conserved without the implementation of small time-steps. Several numerical cases were investigated by the WCMPE method, including 2D jet imping on a solid plate case, oscillating drop, lid driven cavity, rolling tank experiment, sloshing tank simulation and jet impinging from the bottom of a tank. In general, smooth pressure distribution and improved particle stability were obtained with the incorporation of WCMPE.

On the impact of the numerical method on magnetic reconnection and particle acceleration – I. The MHD case

ABSTRACT We present 2D magnetohydrodynamics numerical simulations of tearing-unstable current sheets coupled to a population of non-thermal test particles, in order to address the problem of numerical convergence with respect to grid resolution, numerical method, and physical resistivity. Numerical simulations are performed with the pluto code for astrophysical fluid dynamics through different combinations of Riemann solvers, reconstruction methods, and grid resolutions at various Lundquist numbers. The constrained transport method is employed to control the divergence-free condition of magnetic field. Our results indicate that the reconnection rate of the background tearing-unstable plasma converges only for finite values of the Lundquist number and for sufficiently large grid resolutions. In general, it is found that (for a second-order scheme) the minimum threshold for numerical convergence during the linear phases requires the number of computational zones covering the initial current sheet width to scale roughly as $\sim \sqrt{\bar{S}}$, where $\bar{S}$ is the Lundquist number defined on the current sheet width. On the other hand, the process of particle acceleration is found to be nearly independent of the underlying numerical details inasmuch as the system becomes tearing-unstable and enters in its non-linear stages. In the limit of large $\bar{S}$, the ensuing power-law index quickly converge to p ≈ 1.7, consistently with the fast reconnection regime.

Computational geometric methods for preferential clustering of particle suspensions

A geometric numerical method for simulating suspensions of spherical and non-spherical particles with Stokes drag is proposed. The method combines divergence-free matrix-valued radial basis function interpolation of the fluid velocity field with a splitting method integrator that preserves the sum of the Lyapunov spectrum while mimicking the centrifuge effect of the exact solution. We discuss how breaking the divergence-free condition in the interpolation step can erroneously affect how the volume of the particulate phase evolves under numerical methods. The methods are tested on suspensions of 104 particles evolving in a discrete cellular flow field. The results are that the proposed geometric methods generate more accurate and cost-effective particle distributions compared to conventional methods.

Towards performance portability in the Spark astrophysical magnetohydrodynamics solver in the Flash-X simulation framework

Simulations of core-collapse supernovae, and other astrophysical phenomena, are quintessential extreme-scale computing challenges. For core-collapse supernova simulations to be carried out by the ExaStar project under the Exascale Computing Project umbrella, a robust, efficient, and state-of-the-art magnetohydrodynamics solver is a critical requirement. In Flash-X, the primary software instrument for ExaStar, a new magnetohydrodynamics solver has been designed and implemented from the ground up to achieve accuracy and efficiency for simulations of complex astrophysical flows. This new solver, dubbed Spark, uses high-order spatial reconstruction, Runge–Kutta time integration, and an efficient cell-centered approach to satisfying the divergence-free condition for the magnetic fields. Spark was written to be optimized for data locality in cache hierarchy of CPUs. Since data locality optimizations for cache hierarchy are not directly compatible with those of accelerators, we have taken the approach of using program synthesis to avoid massive amounts of code replication that would be necessary if we were to maintain two different versions of the solver. Our program synthesis relies on a simple key-dictionary approach, implemented in python, that enables us to assemble the version of the solver suitable for the target hardware from code fragments identified by specific keys. In this paper, we describe the data locality optimizations of the solver for CPUs and accelerators and the program synthesis tools that enable this portability. We also detail the parallel performance of Spark for both CPUs and accelerators.

A simple augmented IIM for 3D incompressible two-phase Stokes flows with interfaces and singular forces

In this paper, a Stokes solver on a MAC staggered grid for three-dimensional incompressible two-phase flows with interfaces and singular forces is presented. The main purpose is to develop the simplified IIM to ensure the accuracy near the interface represented by a level-set function. The idea of the simplified IIM is to apply Taylor expansion only along the normal direction and incorporate the jump conditions up to the second normal derivatives into the finite difference scheme. A Newton iteration method is applied to approximate the orthogonal projections for points near the interface. The augmented variables are introduced to decouple the coupled jump conditions of velocity and pressure. The determination of the augmented variables is finally reduced to a linear system, which is solved by the GMRES method. A CG-Uzawa type method is used to solve the discetized Stokes equation with some correction terms. Numerical results show that the proposed algorithm is efficient and can reach global second-order accuracy. The present algorithm can also preserve the area conservation and the discrete divergence free condition very well.

Solution of Volume Integral Equation Using the SWG-Edge Hybrid Basis Functions for Inhomogeneous Dielectric Objects With Multiboundary

A hybrid discretization scheme for solution of volume integral equation (VIE) by method of moments (MoM) for electromagnetic scattering from dielectric objects is proposed in this article. The Schaubert-Wilton-Glisson and edge (SWG-Edge) hybrid basis functions are used in this discretization scheme. According to the divergence-free condition of electric displacement vector, a kind of edge basis functions defined in elements including boundary faces which separate a dielectric object from the background is derived. As a result, we get a SWG-Edge hybrid basis set. Details for the calculation of the corresponding matrix elements for the edge basis and testing functions are presented. Numerical results show the validity and accuracy of the hybrid discretization scheme. Finally, the proposed method is used for efficient solution of VIE for inhomogeneous dielectric objects with multiboundary. It is shown that for multiboundary problems, the number of unknowns of the hybrid basis is only about 71% of the traditional SWG basis. This means that the memory for solution of VIE by the traditional SWG basis functions can be reduced by half. Therefore, the SWG-edge hybrid basis is much more efficient than the traditional SWG basis for solution of VIE with multiboundary problems.

One-stage simplified lattice Boltzmann method for two- and three-dimensional magnetohydrodynamic flows

In this paper, we propose a new simplified lattice Boltzmann method (SLBM) for magnetohydrodynamic flows that outperforms the classical one in terms of accuracy, while preserving its advantages. A very recent paper [De Rosis et al., “Double-D2Q9 lattice Boltzmann models with extended equilibrium for two-dimensional magnetohydrodynamic flows,” Phys. Fluids 33, 035143 (2021)] demonstrated that the SLBM enforces the divergence-free condition of the magnetic field in an excellent manner and involves the lowest amount of virtual memory. However, the SLBM is characterized by the poorest accuracy. Here, the two-stage algorithm that is typical of the SLBM is replaced by a one-stage procedure following the approach devised for non-conductive fluids in a very recent effort [Delgado-Gutierrez et al., “A single-step and simplified graphics processing unit lattice Boltzmann method for high turbulent flows,” Int. J. Numer. Methods Fluids 93, 2339 (2021)]. The Chapman–Enskog expansion formally demonstrates the consistency of the present scheme. The resultant algorithm is very compact and easily implemented. Given all these features, we believe that the proposed approach is an excellent candidate to perform numerical simulations of two- and three-dimensional magnetohydrodynamic flows.

The divergence-free condition in axisymmetric magnetohydrodynamic models

ABSTRACT Axisymmetric magnetohydrodynamic (MHD) models are useful in studies of magnetized winds and non-linear Alfvén waves in solar and stellar atmospheres. We demonstrate that a condition often used in these models for the determination of a nearly vertical magnetic field is applicable to a radial field instead. A general divergence-free condition in curvilinear coordinates is self-consistently derived and used to obtain the correct condition for the variation of a nearly vertical magnetic field. The obtained general divergence-free condition along with the transfield equation completes the set of MHD equations in curvilinear coordinates for axisymmetric motions and could be useful in studies of magnetized stellar winds and non-linear Alfvén waves.

Provably physical-constraint-preserving discontinuous Galerkin methods for multidimensional relativistic MHD equations

We propose and analyze a class of robust, uniformly high-order accurate discontinuous Galerkin (DG) schemes for multidimensional relativistic magnetohydrodynamics (RMHD) on general meshes. A distinct feature of the schemes is their physical-constraint-preserving (PCP) property, i.e., they are proven to preserve the subluminal constraint on the fluid velocity and the positivity of density, pressure, and internal energy. This is the first time that provably PCP high-order schemes are achieved for multidimensional RMHD. Developing PCP high-order schemes for RMHD is highly desirable but remains a challenging task, especially in the multidimensional cases, due to the inherent strong nonlinearity in the constraints and the effect of the magnetic divergence-free condition. Inspired by some crucial observations at the PDE level, we construct the provably PCP schemes by using the locally divergence-free DG schemes of the recently proposed symmetrizable RMHD equations as the base schemes, a limiting technique to enforce the PCP property of the DG solutions, and the strong-stability-preserving methods for time discretization. We rigorously prove the PCP property by using a novel “quasi-linearization” approach to handle the highly nonlinear physical constraints, technical splitting to offset the influence of divergence error, and sophisticated estimates to analyze the beneficial effect of the additional source term in the symmetrizable RMHD system. Several two-dimensional numerical examples are provided to further confirm the PCP property and to demonstrate the accuracy, effectiveness and robustness of the proposed PCP schemes.