Non-physical Oscillations Research Articles

Compact high order finite volume method on unstructured grids III: Variational reconstruction

This paper presents a variational reconstruction for the high order finite volume method in solving the two-dimensional Navier–Stokes equations on arbitrary unstructured grids. In the variational reconstruction, an interfacial jump integration is defined to measure the jumps of the reconstruction polynomial and its spatial derivatives on each cell interface. The system of linear equations to determine the reconstruction polynomials is derived by minimizing the total interfacial jump integration in the computational domain using the variational method. On each control volume, the derived equations are implicit relations between the coefficients of the reconstruction polynomials defined on a compact stencil involving only the current cell and its direct face-neighbors. The reconstruction and time integration coupled iteration method proposed in our previous paper is used to achieve high computational efficiency. A problem-independent shock detector and the WBAP limiter are used to suppress non-physical oscillations in the simulation of flow with discontinuities. The advantages of the finite volume method using the variational reconstruction over the compact least-squares finite volume method proposed in our previous papers are higher accuracy, higher computational efficiency, more flexible boundary treatment and non-singularity of the reconstruction matrix. A number of numerical test cases are solved to verify the accuracy, efficiency and shock-capturing capability of the finite volume method using the variational reconstruction.

Mixed Volume Element-Characteristic Fractional Step Difference Method for Contamination from Nuclear Waste Disposal

A nonlinear system with boundary-initial value conditions of convection–diffusion partial differential equations is presented to describe incompressible nuclear waste disposal contamination in porous media. The flow pressure is determined by an elliptic equation, the concentrations of brine and radionuclide are formulated by convection–diffusion equations, and the transport of temperature is defined by a heat equation. The pressure appears in convection–diffusion equations and heat equation in the form of Darcy velocity and controls the physical processes. The fluid pressure and velocity are solved by the conservative mixed volume element and the computation accuracy of Darcy velocity is improved one order. A combination method of the mixed volume element and the approximation of characteristics is applied to solve the brine and heat, where the diffusion is discretized by a mixed volume element method and the convection is treated by the method of characteristics. The characteristics can confirm strong computation stability at sharp fronts and it can avoid numerical dispersion and nonphysical oscillation. Larger time-steps along the characteristics are shown to result in smaller time-truncation errors than those resulting from standard methods. The mixed volume element method has the property of conservation on each element and it can obtain numerical solutions of the brine and adjoint vectors. The radionuclide is solved by a coupled method of characteristics and fractional step difference. The computational work is reduced greatly by decomposing a three-dimensional problem into three successive one-dimensional problems and using the algorithm of speedup. Using numerical analysis of priori estimates of differential equations, we demonstrate an optimal second order estimate in $$l^2$$ norm. Numerical data are appropriate with the scheme and it is shown that the method is a powerful tool to solve the well-known problems in porous media.

Domain Decomposition Modified with Characteristic Mixed Finite Element and Numerical Analysis for Three-Dimensional Slightly Compressible Oil-Water Seepage Displacement

A parallel algorithm is presented to solve three-dimensional slightly compressible seepage displacement where domain decomposition and characteristics-mixed finite element are combined. Decomposing the computational domain into several subdomains, we define a special function to approximate the derivative at interior boundary explicitly and obtain numerical solutions of the saturation implicitly on subdomains in parallel. The method of characteristics can confirm strong stability at the fronts, and can avoid numerical dispersion and nonphysical oscillation. It can adopt large-time step but can obtain small time truncation error. So a characteristic domain decomposition finite element scheme is put forward to compute the saturation. The flow equation is computed by the method of mixed finite element and numerical accuracy of Darcy velocity is improved one order. For a model problem we apply some techniques such as variation form, domain decomposition, the method of characteristics, the principle of energy, negative norm estimates, induction hypothesis, and the theory of priori estimates of differential equations to derive optimal error estimate in $l^2$ norm. Numerical example is given to testify theoretical analysis and numerical data show that this method is effective in solving actual applications. Then it can solve the well-known problem.

A Double-Passage Shape Correction Method for Predictions of Unsteady Flow and Aeroelasticity in Turbomachinery

AbstractIn this paper, a double-passage shape correction (DPSC) method is presented for simulation of unsteady flows around vibrating blades and aeroelastic prediction. Based on the idea of phase-lagged boundary conditions, the shape correction method was proposed aimed at efficiently dealing with unsteady flow problems in turbomachinery. However, the original single-passage shape correction (SPSC) may show the disadvantage of slow convergence of unsteady solutions and even produce nonphysical oscillation. The reason is found to be related with the disturbances on the circumferential boundaries that can not be damped by numerical schemes. To overcome these difficulties, the DPSC method is adopted here, in which the Fourier coefficients are computed from flow variables at implicit boundaries instead of circumferential boundaries in the SPSC method. This treatment actually reduces the interaction between the calculation of Fourier coefficients and the update of flow variables. Therefore a faster convergence speed could be achieved and also the solution stability is improved. The present method is developed to be suitable for viscous and turbulent flows. And for real three-dimensional (3D) problems, the rotating effects are also considered. For validation, a 2D oscillating turbine cascade, a 3D oscillating flat plate cascade and a 3D practical transonic fan rotor are investigated. Comparisons with experimental data or other solutions and relevant discussions are presented in detail. Numerical results show that the solution accuracy of DPSC method is favorable and at least comparable to the SPSC method. However, fewer iteration cycles are needed to get a converged and stable unsteady solution, which greatly improves the computational efficiency.

A High Order Material Point Method

The classical material point method (MPM) developed in the 90s is known for drawbacks which affect the quality of results. The movement of material points from one element to another leads to non-physical oscillations known as ‘grid crossing errors’. Furthermore, the use of material points as integration points renders a numerical quadrature rule of limited quality. Different solutions have been proposed in recent years to overcome these drawbacks. In this paper the approach of combining quadratic B-spline basis functions with a reconstruction based quadrature rule is pursued to solve these numerical problems. High-order B-spline basis functions solve the problem of grid crossing completely, whereas the considered reconstruction based quadrature rule reduces the quadrature error observed with MPM. In addition, the use of quadratic B-splines leads to a more accurate piecewise linear approximation of the stress field compared to the piecewise constant one obtained with linear Lagrangian basis functions commonly used with MPM. Two 1D benchmarks are considered involving large deformations, a vibrating bar and a column under self-weight. They render excellent results when adopting this high-order MPM.

An Essentially Non-oscillatory Crank–Nicolson Procedure for the Simulation of Convection-Dominated Flows

The Crank–Nicolson (CN) time-stepping procedure incorporating the second-order central spatial scheme is unconditionally stable and strictly non-dissipative for linear convection flows; however, its numerical solution in practice can be oscillatory for nonsmooth solutions. This article studies variants of the CN method for the simulation of linear convection-dominated diffusion flows, in which the explicit convection part is approximated by an upwind scheme, to effectively suppress nonphysical oscillations. The second-order essentially non-oscillatory scheme incorporated in the CN procedure (ENO-CN) has been found effective for a non-oscillatory numerical solution of minimum numerical dissipation. A stability analysis is provided for ENO-CN, which turns out to be unconditionally stable for problems of nonzero diffusion. However, for purely convective flows, it is stable only when the CFL condition is satisfied. Numerical results are presented to demonstrate its stability and accuracy.

A Globally Coupled Pressure Method for the Discretization of the Tensor-Pressure Equation on Non-K-orthogonal Grids

Summary Complex permeability tensors together with general nonorthogonal and unstructured grids pose great challenges to reservoir simulation. The widely used two-point flux approximation (TPFA) is inadequate for a rigorous discretization of the flow equations on such challenging grids. Multipoint flux approximation (MPFA) methods have been proposed to meet the challenges and are currently being deployed in next-generation simulators. In this work, we propose an alternative flux-continuous cell-centered finite-volume method called the globally coupled pressure (GCP) method to discretize the pressure equation on general grids with full permeability tensors. To accurately construct fluxes through control-volume interfaces, pressure at the centroid of all interfaces is introduced as auxiliary unknowns. Flux continuity across each interface gives one equation. Assembling all the flux-continuity equations together gives a system of linear equations that can be solved simultaneously for all the auxiliary unknowns. Flux across control-volume interfaces can then be approximated with the pressure values at control-volume centers only. The fundamental difference between the GCP method and MPFA methods is that, in the latter, auxiliary unknowns are locally coupled within an interaction region and then eliminated in a local stencil by imposing flux-continuity conditions, whereas in the former, all the auxiliary unknowns are globally coupled and can only be eliminated in a global stencil. Consequently, control volumes in our GCP method are directly associated with the edges of the original grid and not by means of a dual grid overlaid and allied with the centers of the grid. Two variants of the GCP method are presented here, and extensive numerical experiments are conducted to test the performance of the GCP methods. The results show that both variants of our GCP method are in good agreement with the classical MPFA-O method on non-K-orthogonal grids for less-challenging problems. Convergence studies reveal that the first variant of our GCP method has slower convergence rates than the MPFA-O method for some problems. However, the second variant of GCP has comparable, and in some cases, better convergence properties compared with the MPFA-O method. With numerical experiments, we further investigate monotonicity properties of our GCP method on highly anisotropic media. For Dirichlet boundary conditions, our GCP methods also suffer from nonphysical oscillations, with some degrees of improvement over the MFPA-O method. When no-flow boundary conditions are used, our GCP method is much more robust and does not produce spurious boundary extrema as MPFA methods do. Finally, we extend our GCP method to fully unstructured grids, and the results show that it is also more robust than the MPFA-O method on unstructured grids.

Sheath formation in an oblique magnetic field - some comments on length scales and source terms

A one-dimensional fluid model is presented and used for the analysis of the potential formation in front of a negative, large, planar electrode immersed in magnetized plasma. The magnetic field lines can form an arbitrary angle with the electrode surface, with the exception that they must not be completely parallel to the surface. The model includes elastic ion collisions and 3 source terms that describe ionization mechanisms. If the exponential source term is used, the pre-sheath is longer than in the case of the constant source term. The zero source term results in qualitatively different solutions that may even exhibit nonphysical oscillations. The problem involves 4 characteristic length scales. Any of them can be used for the normalization of the space coordinate and this has no effect to the solutions.

Numerical modeling of subsidence in saturated porous media: A mass conservative method

In this paper, a second order accurate cell-centered finite volume method (FVM) is coupled with a finite element method (FEM) to solve the deformation of a saturated porous layer based on Biot’s consolidation model. The proposed numerical technique is applied to the fully unstructured triangular grids to simulate actual geological formations. To reconstruct the pressure gradient at control volume faces, the diamond scheme is implemented as a multipoint flux approximation method. Also the least square algorithm is used to interpolate pressure at the vertices from the cell-center values. The stability of this numerical model is studied in comparison to the different FEMs through various examples. It is shown that, although the Taylor-Hood FEM has been introduced as a remedy for violation of the inf-sup condition, it does not entirely remove the non-physical oscillations. Contrary to the linear and Taylor-Hood FEMs, the proposed discretization model provides monotonic solution without imposing any restriction on the mesh or time step size. Compared to the mixed FEM, the method achieves local mass balance with fewer degrees of freedom. To couple the flow and mechanical sub-problems, the fixed-stress operator split is implemented as an iterative sequential method, due to its unconditional stability, accuracy and high rate of convergence. The accuracy of the proposed model is verified via a range of examples including analytical and numerical solutions. The performance of this methodology is assessed through modeling of subsidence in an aquifer-interbed system. This problem illustrates the capability of the model in providing stable solution in heterogeneous domains with complicated shapes.

Simulation of nonlinear propagation of biomedical ultrasound using pzflex and the Khokhlov-Zabolotskaya-Kuznetsov Texas code.

Nonlinear acoustics plays an important role in both diagnostic and therapeutic applications of biomedical ultrasound and a number of research and commercial software packages are available. In this manuscript, predictions of two solvers available in a commercial software package, pzflex, one using the finite-element-method (FEM) and the other a pseudo-spectral method, spectralflex, are compared with measurements and the Khokhlov-Zabolotskaya-Kuznetsov (KZK) Texas code (a finite-difference time-domain algorithm). The pzflex methods solve the continuity equation, momentum equation and equation of state where they account for nonlinearity to second order whereas the KZK code solves a nonlinear wave equation with a paraxial approximation for diffraction. Measurements of the field from a single element 3.3 MHz focused transducer were compared with the simulations and there was good agreement for the fundamental frequency and the harmonics; however the FEM pzflex solver incurred a high computational cost to achieve equivalent accuracy. In addition, pzflex results exhibited non-physical oscillations in the spatial distribution of harmonics when the amplitudes were relatively low. It was found that spectralflex was able to accurately capture the nonlinear fields at reasonable computational cost. These results emphasize the need to benchmark nonlinear simulations before using codes as predictive tools.

Compact grid-characteristic schemes of higher orders of accuracy for a 3D linear transport equation

A numerical solution of a 3D linear transport equation on parallelepipedic computational grids is considered. By the technique of splitting in coordinates, compact grid-characteristic schemes of higher orders of accuracy are generalized to the 3D case. The influence of particular steps of the computational algorithm on the accuracy of the resulting scheme is investigated. The approach for retaining the order of convergence of a scheme on a smooth solution and minimizing nonphysical oscillations on a discontinuous solution in the 3D case is proposed.

Comparison and modification: TVD schemes for scalar transport on an unstructured grid

In most TVD schemes, the r-factors were proposed according to the cell-centered (CC) finite volume method (FVM) framework for the numerical approximation to the convective term. However, it is questionable whether those r-factors would be appropriate and effective for the vertex-centered (VC) FVM. In the paper, we collected five kinds of r-factor formulae and found out that only three of those, respectively by Bruner (1996), Darwish and Moukalled (2003) and Cassuli and Zanolli (2005) can be formally extended to a context of the VC FVM. Numerical tests indicate that the TVD schemes and r-factors, after being extended and introduced to a context of the VC FVM, maintained their similar characteristics as in a context of the CC FVM. However, when the gradient-based r-factors and the SUPERBEE scheme were applied simultaneously, non-physical oscillations near the sharp step would appear. In the transient case, the oscillations were weaker in a context of the VC FVM than those in a context of the CC FVM, while the effect was reversed in the steady case. To eliminate disadvantages in the gradient-based r-factor formula, a new modification method by limiting values on the virtual node, namely ϕU in the paper, was validated by the tests to effectively dissipate spurious oscillations.

Numerical Modeling of 2-D and 3-D Flows using Artificial Compressibility Method and Collocated Mesh

In this paper, applications of a numerical model on simulation of two and three-dimensional flows are presented. This model solves Navier-Stokes equations using finite volume method and large eddy simulation (LES) in a collocated mesh. Artificial compressibility method with dual t ime stepping is used to solve the time dependent equations. Also a modified m omentum i nterpolation method (MIM) based on the unsteady flows i s deployed t o overcome t he non-physical pressure oscillation. Capability of the presented numerical code for flow s imulation, i s a ssessed by a pplication f or twodimensional square and three-dimensional lid-driven cavity flows. Numerical r esults of cavity flow presents very good agreement with the numerical and experimental data of other existent researches.

Study of high-order-accurate limiters for time-dependent contact discontinuity and shock capturing

ABSTRACTIn most applications of computational fluid dynamics (CFD), it is very important to correctly project the interaction phenomena when we compute the typical discontinuous structure. Numerical methods for the conservation law are usually divided into two kinds: One is the “single stepping method”; another is the “semidiscrete scheme”. When the semidiscrete scheme handles the issues of high accuracy and operational process without oscillations, it is a very successful method. In this article, we consider the one-dimensional Euler equations for the Riemann problem of a typical shock tube without external active force to improve the nonphysical numerical oscillation problem near the contact discontinuity interface and shock wave. We try to improve the computational accuracy of numerical resolution of the Riemann problem and reduce nonphysical numerical oscillations through the AUSMDV numerical flux scheme with implementing distinct flux limiters and time-stepping methods. Using the conditions of single gas/dual gases, we carry out a series of comparisons for the accuracy of shock capturing as well as the capability of improving numerical oscillations.

Numerical solution of unsteady advection dispersion equation arising in contaminant transport through porous media using neural networks

A soft computing approach based on artificial neural network (ANN) and optimization is presented for the numerical solution of the unsteady one-dimensional advection–dispersion equation (ADE) arising in contaminant transport through porous media. A length factor ANN method, based on automatic satisfaction of arbitrary boundary conditions (BCs) was chosen for the numerical solution of ADE. The strength of ANN is exploited to construct a trial approximate solution (TAS) for ADE in a way that it satisfies the initial or BCs exactly. An unsupervised error is constructed in approximating the solution of ADE which is minimized by training ANN using gradient descent algorithm (GDA). Two challenging test problems of ADE are considered in this paper, in which, the first problem has steep boundary layers near x=1 and many numerical methods create non-physical oscillation near steep boundaries. Also for the second problem many numerical schemes suffer from computational noise and instability issues. The proposed method is advantageous as it does not require temporal discretization for the solution of the ADEs as well as it does not suffer from numerical instability. The reliability and effectiveness of the presented algorithm is investigated by sufficient large number of independent runs and comparison of results with other existing numerical methods. The results show that the present method removes the difficulties arising in the solution of the ADEs and provides solution with good accuracy.

A Robust Immersed Boundary-Lattice Boltzmann Method for Simulation of Fluid-Structure Interaction Problems

AbstractA robust immersed boundary-lattice Boltzmann method (IB-LBM) is proposed to simulate fluid-structure interaction (FSI) problems in this work. Compared with the conventional IB-LBM, the current method employs the fractional step technique to solve the lattice Boltzmann equation (LBE) with a forcing term. Consequently, the non-physical oscillation of body force calculation, which is frequently encountered in the traditional IB-LBM, is suppressed greatly. It is of importance for the simulation of FSI problems. In the meanwhile, the no-slip boundary condition is strictly satisfied by using the velocity correction scheme. Moreover, based on the relationship between the velocity correction and forcing term, the boundary force can be calculated accurately and easily. A few test cases are first performed to validate the current method. Subsequently, a series of FSI problems, including the vortex-induced vibration of a circular cylinder, an elastic filament flapping in the wake of a fixed cylinder and sedimentation of particles, are simulated. Based on the good agreement between the current results and those in the literature, it is demonstrated that the proposed IB-LBM has the capability to handle various FSI problems effectively.

Compact high order finite volume method on unstructured grids I: Basic formulations and one-dimensional schemes

The large reconstruction stencil has been the major bottleneck problem in developing high order finite volume schemes on unstructured grids. This paper presents a compact reconstruction procedure for arbitrarily high order finite volume method on unstructured grids to overcome this shortcoming. In this procedure, a set of constitutive relations are constructed by requiring the reconstruction polynomial and its derivatives on the control volume of interest to conserve their averages on face-neighboring cells. These relations result in an over-determined linear equation system, which, in the sense of least-squares, can be reduced to a block-tridiagonal system in the one-dimensional case. The one-dimensional formulations of the reconstruction are discussed in detail and a Fourier analysis is presented to study the dispersion/dissipation and stability properties. The WBAP limiter based on the secondary reconstruction is used to suppress the non-physical oscillations near discontinuities while achieve high order accuracy in smooth regions of the solution. Numerical results demonstrate the method's high order accuracy, robustness and shock capturing capability. A new reconstruction procedure for high order finite volume scheme is proposed.This procedure can achieve arbitrary order of accuracy using a compact stencil consisting of only face neighboring cells.The spectral properties of the proposed finite volume schemes are analyzed.Some test cases are solved numerically to show the superior performances of the proposed schemes.

A new fifth order finite difference WENO scheme for solving hyperbolic conservation laws

In this paper a new simple fifth order weighted essentially non-oscillatory (WENO) scheme is presented in the finite difference framework for solving the hyperbolic conservation laws. The new WENO scheme is a convex combination of a fourth degree polynomial with two linear polynomials in a traditional WENO fashion. This new fifth order WENO scheme uses the same five-point information as the classical fifth order WENO scheme 14,20, could get less absolute truncation errors in L 1 and L ∞ norms, and obtain the same accuracy order in smooth region containing complicated numerical solution structures simultaneously escaping nonphysical oscillations adjacent strong shocks or contact discontinuities. The associated linear weights are artificially set to be any random positive numbers with the only requirement that their sum equals one. New nonlinear weights are proposed for the purpose of sustaining the optimal fifth order accuracy. The new WENO scheme has advantages over the classical WENO scheme 14,20 in its simplicity and easy extension to higher dimensions. Some benchmark numerical tests are performed to illustrate the capability of this new fifth order WENO scheme.

Influence of momentum interpolation methods on the accuracy and convergence of pressure–velocity coupling algorithms in OpenFOAM®

Checkerboard pressure is one of the most significant problems arising from the use of collocated grids for fluid dynamic simulations using the finite-volume method. The original Rhie–Chow momentum interpolation technique, termed Original Momentum Interpolation Method (OMIM) was proposed to eliminate the non-physical saw-tooth pressure oscillations. However, it was soon proved that the steady-state solutions obtained with this technique were under-relaxation factor dependent. Nevertheless, standard OMIM is still commonly used in several CFD codes like OpenFOAM®, the software used for this work. In this paper the OMIM and a possible correction for under-relaxation dependency, which has been implemented in OpenFOAM®, are discussed in detail. The proposed methodology is compared and contrasted with OMIM in terms of accuracy of the solution and speed of convergence for several classical pressure–velocity segregated algorithms for steady state solvers; namely SIMPLE, SIMPLE-C, SIMPLE-R and PISO. A classical laminar 2D cavity is used as the base test-case. The study is then extended to a more complex 2D airfoil profile (NACA0012). In the cases considered mesh uniformity and orthogonality are progressively reduced and turbulence starts playing an important role, limiting therefore the convergence of the cases and the performance of the correction.

Phase structure of cold magnetized color superconducting quark matter

The influence of intense magnetic fields on the behavior of color superconducting cold quark matter is investigated using an SU(2)f NJL-type model for which a novel regulation scheme is introduced. In such a scheme the contributions which are explicitly dependent on the magnetic field turn out to be finite and, thus, do not require to be regularized. As a result of this, non-physical oscillations that arise from regularizing magnetic field dependent terms are naturally removed, and oscillations that are actually physical can be better appreciated. The phase diagrams in the ẽB – μ plane are presented for different values of the diquark coupling.