Simulation results and applications of an advection bounded scheme to practical flows

Simulation results and applications of an advection bounded scheme to practical flows

Valdemir G. Ferreira^I; Rafael A.B. de Queiroz^I; Miguel A. Caro Candezano^I; Giseli A.B. Lima^I; Laís Corrêa^I; Cassio M. Oishi^II; Fernando L.P. Santos^III

^IDepartamento de Matemática Aplicada e Estatística Instituto de Ciências Matemáticas e de Computação - USP, São Carlos, SP, Brazil. E-mails: [e-mail] / [e-mail] / [e-mail] / [e-mail] / [e-mail]

^IIDepartamento de Matemática Computação Universidade Estadual Júlio de Mesquita Filho - UNESP, Presidente Prudente, SP, Brazil. Email: [e-mail]

^IIIDepartamento de Bioestatística, Instituto de Biociências de Botucatu - UNESP, Botucatu, SP, Brazil. Email: [e-mail]

ABSTRACT

This paper reports experiments on the use of a recently introduced advection bounded upwinding scheme, namely TOPUS (Computers & Fluids 57 (2012) 208-224), for flows of practical interest. The numerical results are compared against analytical, numerical and experimental data and show good agreement with them. It is concluded that the TOPUS scheme is a competent, powerful and generic scheme for complex flow phenomena.

Mathematical subject classification: Primary: 06B10; Secondary: 06D05.

Key words: convection term discretization, convective upwinding scheme, Navier-Stokes equations, Euler equations, time-dependent fluid flow, advection modelling.

1 Introduction

The need to solve advection-dominated PDEs (Partial Differential Equations) is ubiquitous throughout computational fluid dynamics applications. In order to achieve physically relevant numerical solutions for these equations, one inevitably requires to equip the computational algorithm with some high resolution upwind (bias) scheme for the convective fluxes.

High resolution upwind schemes, an extension of the monotonicity preserving first-order upwind scheme by the use of non-linear limiters, have successful been employed for the simulation of a variety of PDEs; since advection of scalars to non-linear conservation laws (see, for instance, [1]). However, their adaptation to PDEs for predicting flow field in the presence of shocks or steep gradients is not so common in the literature and, in particular, their application to incompressible free surface flows at high Reynolds numbers is hindered by the moving boundary. For achieving this goal, we have described in [2] a high degree polynomial upwind-based TVD (Total Variation Diminishing) finite difference scheme, called TOPUS (Third-Order Polynomial Upwind Scheme).

TOPUSis based on the application of theTVD/CBC(Convection Boundedness Criterion) stability criteria combined with the conditions of Leonard [3]. This scheme has been presented (see [2]) in both the normalized variables and also as a flux limiting technique, and has been shown to possess three important features: simplicity, robustness and generality of application. The main point of that paper was to demonstrate that the TOPUS scheme can be employed to solve a wide range of linear and non-linear PDE, preserving total variation as time integration evolves. The authors have also made a variety of comparisons of different upwind TVD schemes, namely CUBISTA [4] , ADBQUICKEST [5], SMART [6], SUPERBEE [7], van Albada [8], and van Leer [9, 10, 11], with the TOPUS scheme for several computational fluid dynamics benchmark cases, such as advection of scalars, gas dynamics and simple flows. However, even with the relative success obtained with the (original) TOPUS scheme [2], there is still the need for improved it to enable fluid flow computations to be used routinely in engineering practice.

In the present article, we aim to extend the TOPUS scheme to simulate fluid flow problems of increasing complexity, and to answer the basic question: what reliable, accurate and easy-to-program upwinding scheme should be employed for fluid dynamics with jump discontinuities? The emphasis of the study is not to trace the evolution of upwind-biased schemes, nor to provide rigorous analysis for them, but rather to present the versatility of the TOPUS scheme for resolving more complicated PDEs than those reported in authors's earlier paper [2]. In addition, the paper intends to supply results of simulations for both (representative) compressible and incompressible flows. These computations aim to acquaint the researcher in computational fluid dynamics with the virtues of the TOPUS scheme.

The remainder of the paper is organized as follows. In the next section, it is presented a summary of the TOPUS scheme and its modification used in this work. In Section 3, computational results and applications are made by means of a series of numerical experiments for a variety of PDEs. Finally, some concluding remarks are drawn in Section 4.

2 Summary of the TOPUS scheme and its modification

In this section, we review the (original) TOPUS scheme [2], an upwinding technique for approximation of cell interface values in reconstruction formulas, and its modification on Cartesian meshes in the framework of the finite difference method.

The normalized generic form (see Leonard [3]) of a high resolution upwind scheme for advection term discretizations is given by

where and are, respectively, the normalized values of the convected variable ϕ at the boundary interface f between two control volumes and at the neighboring upwind node U. According to Leonard [3, 12], a bounded high resolution second and/or third order accurate scheme (in general, non-linear) within the CBC region must pass through points O(0,0), Q(0.5,0.75), P(1,1) and with inclination of 0.75 at Q. Passing through Q will provide second order accuracy and passing through Q with a slope of 0.75 will give third order accuracy.

A modification of the TOPUS scheme is constructed by assuming that thevariable is related to by a fourth degree polynomial function for 0 < < 1, and by the NECBC1 (New CBC) scheme of Jian et al. [13] for < 0 and > 1. In the original TOPUS scheme (see [2]), the NECBC1 scheme is replaced by the FOU (First Order Upwind) scheme. By imposingthe conditions of Leonard presented above, plus the condition that is a continuously differentiable function at P, one obtains

where

is the normalized variable of Leonard and α is a free parameter. The notations ϕ_D, ϕ_U and ϕ_R represent, respectively, the value of the variable ϕ at Downstream, Upstream and Remote-upstream locations, which are selected according to the sign of the advection velocity (upwind/downwind direction) at the interface f. The motivation for choosing this compact stencil for TOPUS scheme is its computational efficiency, reducing the amount of storage and exchange of information, and simplifying the implementation of the boundary conditions.

The corresponding flux limiter for the original TOPUS scheme, which was used for the simulation of conservation laws (see [2]), can be written, in a commonly used notation, as

where r_fis a local smoothness measure satisfying Sweby's monotonicity preservation condition (see, for example, [1]) when it tends to zero, and it is given by

where Δ and ∇ are forward and backward difference operators, respectively. The original TOPUS scheme inside the TVD region of Harten can be found in [2]. The results and applications presented in the next section were generated by using the free parameter α equal to 2, since for this value of α the scheme has guaranteed to be oscillation-free and has provided satisfactory results for standard problems (see [2]).

It should be stated that, in the most simulations presented in this work, divergence has been observed with the use of the original (frozen TOPUS) version of the TOPUS scheme (2) (for compressible flows). And, from the insight gained of the modified van Albada limiter (see [14]), the following improvement to the TOPUS limiter is proposed

where a = ϕ_D - ϕ_U and b = ϕ_U - ϕ_R; ε is a small parameter that prevents indeterminacy in regions of zero gradients and is taken to be of O(δ_x³), δ_x being the mesh spacing. This modification for the TOPUS limiter is similar in many respects to that of [8, 14] and [15], and is being exploited in this work to prevent spurious numerical oscillations at large flow gradients and discontinuities (e.g. shocks), improve numerical (global) convergence, and handle the clipping and squaring effects of smooth extrema.

The reader is referred to Ferreira et al. [2] to see how TOPUS may be incorporated into the discretized form of a number of PDEs. In addition, in this reference, it is provided a discussion concerning the stability of the computations and the choice of the CFL (Courant-Friedrichs-Lewy) parameter. The issue of non-linear stability for the TOPUS scheme is also addressed in [2] by checking the numerical time dependent total variation on progressively refined meshes. We should, however, remark that in practice it seems that for many problems spatial accuracy is more crucial than temporal accuracy; hence in the calculations presented in this paper, for simplicity, the first order Euler method hasbeen used for marching in time.

Computations of the inviscid flows are performed by using the CLAWPACK software package of LeVeque [16], incremented with the TOPUS and other upwinding schemes presented in this work. CLAWPACK is a general purpose and open-source software developed at the University of Washington for resolving hyperbolic PDEs. This software, in the context of finite volume methodology, uses the Godunov's method with a correction term; one can provide additional code for initial and boundary conditions and add new limiters. Incompressible viscous flows are simulated by using the primitive variable (front-tracking/finite difference) Marker-And-Cell (MAC) technique: this is a special case of the projection method of Chorin [17] described by Harlow and Welch [18] (see also McKee et al. [19]). This method, associated with the upwind TOPUS scheme (and other upwinding schemes), has been incorporated into the Freeflow code of Castelo et al. [20] to solve complex incompressible moving free surface flows.

3 Simulation results and applications

In authors's earlier paper [2], representative 1D/2D test cases have been presented so that the methods discussed could be compared. In this section, other 2D validation cases, as well as verification tests, will be presented. The flexibility and robustness of the TOPUS scheme are considered by solving complex flows in two or three space dimension. In both cases, a comparison with well recognizedhigh-resolution schemes is also performed. The objective here is to investigate whether the TOPUS scheme could effectively solve a real engineering problem. Seven test cases have been selected in order to assess the behaviour of the TOPUS scheme, namely three inviscid compressible flows and four viscous incompressible flows with moving free surfaces.

3.1 2D inviscid compressible flows

In this section, the TOPUS scheme is used to compute 2D non-linear system of hyperbolic PDEs of the form

where U is the conservative state vector, and F(U) and G(U) are the convective flux vectors along the x- and y-directions, respectively. The specific flows simulated here are: (i) the circular dam-break, modelled by the inviscid shallow water equations; (ii) the steady transonic flow around the NACA 0012 airfoil, modelled by Euler's equations; and (iii) the compressible Orszag-Tang MHD vortex, modelled by ideal magnetohydrodynamics (MHD) equations.

Test case 1 (Circular dam-break problem). The TOPUS scheme is initially tested for the simulation of the collapse of a circular dam, a free surface shallow flow described in details by Stecca et al. [21] and modelled by Eq. (5) with the conservative state vector and convective flux vectors given by

where h is the water depth; u and v are the x and y velocities, respectively; and g is the acceleration of the gravity. The aim in this test is to demonstrate the ability of the TOPUS scheme of accurately reproducing shock and rarefaction waves. As presented by Stecca and co-authors, the problem consists of the instantaneous breaking of a cylindrical tank initially filled with 2.5 m deep waterat rest. The circular column of water is suddenly released and, then, a shockwave propagates in the radial direction (outwards) while a rarefaction wave moves inwards. The wave generated by the breaking of the tank propagates into still water with an initial depth of 0.5 m. The computed results for the water depth h(x, y, t) (sliced on the x-axis) are compared with the results provided by Stecca et al. [21]. The mesh used was 100 × 100 computational cells and the reference solution was obtained by using the SUPERBEE scheme on a fine mesh of 1000 × 1000 cells and at CFL number 0.9. Figure 1 displays comparisons, using two different values of the CFL number (CFL=0.10 and 0.45), between the results of [21] and those obtained with the TOPUS scheme. From this figure, one can see that the numerical method equipped with the TOPUS captures the essential physical mechanism of the problem and provides the best numerical solution. Moreover, the scheme shown to be less dissipative than other schemes. In addition, one can observe some oscillations (almost imperceptible) appearing in Figure 1(a); this can easily be removed by refining the mesh.

Test case 2 (Transonic flow around the NACA 0012 airfoil). The second problem is that of a steady inviscid compressible flow over a NACA 0012 airfoil at freestream Mach number M_∞ = 0.8 and angle-of-attack α_at = 1.25 deg. This problem is modelled by Eq. (5) with the conservative state and convective flux vectors given by

where E is the total energy, ρ is the density, and p is the pressure; other variables have been defined previously. This test case is computed using a mesh size of 251 points over the airfoil surface, 151 points in radial direction and the farfield boundary is set at 70 chords of radius. The CFL number is set as a constant value of 0.8 and the maximum density residual for accepting convergence is chosen to be 10^-7. The pressure coefficient distributions, C_p, on the upper and lower surfaces of the airfoil obtained with TOPUS and van Albada limiters are plotted in Figure 2. From this figure, it is seen that both TOPUS and van Albada limiters provide similar results, showing that the strength of the shock is in good agreement with the ones given in literature. These data also indicate that TOPUS is slightly less dissipative than the van Albada limiter at the shock. Away from the shock waves, both TOPUS and van Albada schemes produce almost identical results.

Further investigation of these results can be achieved by inspecting the entropy generated by the numerical solutions. Hence, Figure 3 presents the entropy generated at the airfoil surface by the two schemes for the same flight condition.

The clear conclusion from this figure is that the entropy generated by the two schemes is quite comparable. One can see that TOPUS creates slightly more entropy at the airfoil surface than the van Albada limiter. Again, these results emphasize that TOPUS has essentially the same shock capturing characteristics as the widely used van Albada limiter for such inviscid transonic applications.

Finally, drag and lift coefficients (C_d and C_l) are summarized in Table 1; besides the comparison between TOPUS and van Albada schemes, we have included results for the present test case obtained by Zhou et al. [22], Caughey [23], Rizzi [24], Jameson and Martinelli [25], Pulliam and Barton [25], Du et al. [27] and Venkatakrishnan [14]. Such data provide for a more quantitative comparison of the presently proposed scheme. One can see in Table 1 that the present results for lift and drag coefficients are between those provided by the van Albada limiter and those provided by the centered schemes. Again, the current results are very close to those provided by the van Albada limiter, except that we obtain a slightly higher value of lift coefficient, which is probably a consequence of the less dissipative behavior at the shock, as previously discussed, and also a somewhat higher drag coefficient. We believe that the higher drag coefficient is associated with the fact that TOPUS is generating slightly more entropy at the airfoil surface than the van Albada limiter, as indicated in Figure 3. Hence, TOPUS produces more spurious drag than the van Albada limiter, explaining the higher C_d values. However, one should notice that, clearly, such additional spurious drag is quite lower than what is generated by the other schemes compared in Table 1. Furthermore, the current results for both lift and drag coefficients are well within the ranges reported in literature.

Test case 3 (The Orszag-Tang MHD vortex). The interaction of a moving plasma with a magnetic field which produces shocks, vortices and other smooth structures is simulated here. This vortical flow field contains many significant features of MHD turbulence and has been a challenging benchmark test to check the accuracy of upwinding schemes (see [28] and references within). When compared to the Navier-Stokes equations, the MHD equations are more complicated since they support a family of waves that propagates at different speeds in an anisotropic manner. This problem is modelled by Eq. (5) with the conservative state vector and convective flux vectors given by

where u = (u, v, w)^T and B = (B_x, B_z, B_z)^T represent, respectively, the velocity and magnetic fields, and p^* = p + (1/2)|B|² is the total pressure.For the simulation of this complex flow, the computational domain is set as [0.2, π] × [0.2, π], with double-periodic boundary conditions and initial conditions given by

Figure 4 shows this property along the line z = 0.625π obtained with the ADBQUICKEST [5] and TOPUS schemes at CFL of 0.8, and the scheme of Balbas et al. [28] at the more restrictive CFL condition of 0.4. From this figure, one can see that the results with TOPUS are comparable to the published results in above reference and show its ability to capture shocks sharply as well as resolving the central vortex.

In the following, the observed accuracy of the TOPUS scheme on this complex flow is assessed by marching to a fixed time of t = 5 at CFL number of 0.75. Table 2 gives the L₁ errors and the corresponding orders of convergence for the TOPUS, ARORA-ROE [7], ADBQUICKEST, SUPERBEE and MC schemes. One can see that, practically, the same order (≈ 2.56 as the mesh is refined) is observed for all schemes.

3.2 2D and 2.5D viscous incompressible flows

From now on, computational results for 2D and axisymmetric (2.5D) viscous incompressible flows involving free surfaces are presented. The basic equations for the simulation of incompressible fluid flows are the Navier-Stokes and mass conservation equations which describe the conservation of momentum and mass, respectively. In Cartesian or cylindrical coordinates, these equations are given by

where u = u(r, z, t) and v = v(r, z, t) are, respectively, the components in the r- and z-directions of the local velocity vector field of the fluid; p is the ratio of scalar pressure field to constant density. The non-dimensional parameters Re = U₀L₀/ ν and Fr = denote the associated Reynolds and Froude numbers, respectively, in which U₀ is a characteristic velocity scale, L₀ is a characteristic length scale, and ν is the kinematic molecular viscosity coefficient (constant) and g = [g_r, g_z]^T is the gravitational acceleration. The parameter τ in Eqs. (6) and (8) is used to specify the coordinate system, namely: when τ = 0, Cartesian coordinates are considered (r is interpreted as x and z as y);and when τ = 1, cylindrical coordinates are assumed.

Test case 4 (A 2D liquid jet impinging onto a solid smooth surface at high Reynolds number). This problem concerns a 2D smooth fluid jet impinging normally onto a horizontal surface at high Reynolds number. This free surface flow (in a laminar regime) has been selected because there is (see [29]) an approximate analytical solution for the thickness of the fluid layer flowing on the horizontal (rigid) surface. It is difficult to simulate this problem because the free surface boundary conditions must be specified on an arbitrarily moving boundary. The calculations reported below were obtained by using the 2D version of the Freeflow code [20]. This code, equipped with TOPUS, ADBQUICKEST and CUBISTA schemes, ran this problem at a Reynolds number of 2.0 ×10³, which was based on the maximum velocity U₀ = 1.0 m/s and diameter of the inlet 2a = 0.02 m. The distance between the inflow section and the rigid wall (the inflow-to-plate distance) was 0.037 m. The boundary conditions were the usual no-slip at the solid surface and no-shear stress at the free surface. A mesh size of 200 × 50 (δ_x = δ_y = 0.001 m) computational cells was employed. By using this mesh, a comparison was made between the free surface height (the total thickness of the layer h), obtained from the numerical method and from the analytical viscous solution of Watson [29]. This is displayed in Figure 5 using the ADBQUICKEST, CUBISTA and TOPUS schemes. One can see from this figure that the numerical results using the TOPUS scheme are generally in good agreement with the analytical solution, displaying small differences in some regions of the flow. It can also be observed, from this figure, that ADBQUICKEST and CUBISTA gave similar results, with TOPUS providing marginally better ones.

Test case 5 (A 2.5D liquid jet impinging onto a solid smooth surface at moderate Reynolds number: a stationary circular hydraulic jump). When a 2.5D jet of liquid impinges on a flat (horizontal) plate it can, for certain moderate values of the Reynolds number, create a (circular) hydraulic jump. This occurs at a critical radius, where there is a sudden transition from shallow rapidly flowing fluid to deep, much slower flowing fluid. A better understanding of this phenomenon and the instabilities when it is turbulent is of commercial interest, since jet impingement is often used in cooling systems and the flow of the fluid beyond the jump can degrade the efficiency of the system. Probably, the first author to study the influence of fluid viscosity on the jump radius was Watson [29]. The purpose of this test is three-fold. Firstly, we wish to demonstrate that the TOPUS scheme is capable of simulating this complex moving free surface flow, with emphasis on the total thickness of the fluid layer h of the liquid and on the position of the jump. Secondly, we wish to compare the numerical solutions obtained with the approximate analytic results of Watson [29] where appropriate. Finally, we wish to check the effect of the numerical parameters δ_x = δ_y (spatial resolution) and the δ_t (time step) on the numerical results. The study was carried out by varying one parameter while keeping the other constants.

We begin by verifying that the TOPUS scheme provides good estimates for the position of the jump. For this, the scaling relations for the radius of the jump

of Brechet and Néda [30], and r_jump = Q^5/8ν^-3/8g^-1/8 of Bohr et al. [31] were used for comparison. The radius of the inlet a₀ = 0.008 m and the velocity of the fluid at this boundary U₀ = 3.75 ×10^-1 ms^-1 have been used as the scaling parameters. The jet flow rate Q = πU₀a₀² = νRea₀ = 0.75 × 10^-5 m³s^-1, producing a Reynolds number of 250, and a constant inflow-to-plate distance of H = 0.03 m were employed in the simulations. The jump was identified as the location where the derivative possesses its maximum (see [32]). The 2.5D version of the Freeflow code [20] equipped with the TOPUS scheme was run on this problem using three meshes, namely 200 × 126 (δ_x = δ_y = 0.00025 m); 400 × 252 (δ_x = δ_y = 0.000125 m) and 800 × 504 (δ_x = δ_y = 0.000625 m) computational cells (known hereafter as Mesh I, Mesh II and Mesh III, respectively). Table 3 shows the jump radii obtained from the simulation results and the theoretical scaling laws. One can see that the calculated estimates for the jump with the TOPUS scheme on the three meshes, particularly the one on the fine mesh (Mesh III), are in reasonable agreement with the theoretical scaling law of Brechet and Néda [30], a result that has later been experimentally verified by Hansen et al. [33].

A comparison was then performed between the fluid layer h, obtained from the numerical results and the viscous/inviscid solution of Watson [29]; this is displayed in Figure 6. The numerical solutions were calculated by using Meshes I, II and III at a time step of 1.3 ×10^-4s. We restricted the analysis to the region 0.2 < < 0.8 because Watson's analysis is only valid under the restriction r >> a₀ and the presence of the outflow-boundary (see [29]). One can clearly see that the numerical solutions do not show oscillatory behavior, and as the mesh size is decreased, the solution converges, indicating the convergence of the method for this complex non-linear free surface flow. Watson's approximate solution is only valid over a restricted range of r and the results presented in Figure 6 are extremely good over that range. In fact, when the mesh was refined once more (Mesh IV = 1600 × 1008) the numerical solution (not shown) converged to a solution very close to that obtained in Mesh III.

Finally, in order to check the effect of the time step on the numerical solution we compute, on the Mesh I (200 ×126), the fluid layer h using different time steps (from 10^-3s to 10^-6s). In Figure 7, the numerical results of three simulations using the time steps 1.3 ×10^-4s, 6.5 ×10^-5s and 2.7 ×10^-6s are presented. It can be seen that no significant effect was detected in the numerical solutions by reducing or increasing the time step. The independence of results with such time step variations shows that the time step of the order of 10^-4 used for the computation of the circular hydraulic jump is appropriate for obtaining accurate results.

3.3 Applications: Full 3D moving free surface flows

We conclude this paper by demonstrating the applicability of the TOPUS to more realistic engineering problems, namely: collapse of a liquid block and circular hydraulic jump. These flows are of significant industrial and environmental importance but are difficult to simulate because the boundary conditions must be specified on an arbitrarily moving surface. The governing equations for simulating unsteady incompressible Newtonian free surface flow in three space dimensions are the momentum equations and the continuity equation. In index notation they are, respectively, given by

where all quantities have previously been defined. The 3D version of the Freeflow simulation system (see [20]), equipped with the MAC methodology and TOPUS scheme, was used in a similar way to previous sections. Details of the free surface boundary conditions can be found in [19].

Test case 6 (Collapse of a liquid portion of fluid). Results are presented nowfor the collapse of a column of water onto a horizontal wall for 2D and 3D cases. This free surface flow problem was first studied experimentally in detail by Martin and Moyce [34], and more recently by Koshizuka and Oka [35] to investigate the spreading velocity and the falling rate of water columns. By using the TOPUS scheme, we performed a simulation of this unsteady free surface flow. The geometry used is a fluid column (a = 0.05 m wide and 2a = 0.1 m high) in the 2D case and a fluid block (b = 0.05 m length, a = 2b width and c = 2b height) in the 3D case, both in hydrostatic equilibrium and confined between walls. Initially, a wall is instantaneously removed and the fluid is subject to gravity and free to flow out along a rigid horizontal wall. Our transient 3D numerical simulation of this free surface problem is illustrated in Figure 8.

In order to compare with the experimental data given by Martin and Moyce [34] and Koshizuka and Oka [35], the free-slip boundary condition was used to model the flow at the walls. The Reynolds number based on the characteristic length D = 2a and the characteristic velocity U = was chosen to be Re ≈ 99 ×10³ (|g| = 9.81 ms^-2). The meshes used in this problem were: 150 × 75 (δ_x = δ_y = 0.002 m) computational cells in the 2D case; and 150 × 50 × 75 (δ_x = δ_y = δ_z = 0.002 m) computational cells in the 3D case. Figure 9 shows the 2D/3D numerical results and the experimental data for the position of the fluid front Xmax versus time (the 3D numerical results were obtained by a cutting plane at position y = 0.05m). As shown in this figure, both 2D and 3D calculations with TOPUS agree fairly well with the experimental data especially in comparison with the results of Martin and Moyce [34] and Koshizuka and Oka [35]. In order to provide a stiffer test for the performance of TOPUS we compared the calculated surge front position (as a function of non-dimensional time) against other sophisticated techniques, for example, smoothed particle hydrodynamics (SPH), boundary element method (BEM), level set and an approach by Ritter (see Colagrossi and Landrini [36]). Figure 10 displays this comparison and, once more, TOPUS compares very favorably with these more recent results, giving us confidence in the numerical solution.

Test case 7 (Circular hydraulic jump). In a similar manner to the 2.5D simulation case (see Test case 5.), a 3D jet of viscous fluid at high Reynolds number was projected onto a horizontal rigid wall with an appropriate prescribed velocity U₀, so that a hydraulic jump would occur. The Reynolds number, based on the maximum velocity U₀ = U_max = 1.0 m/s and diameter of the inlet D = 0.05 m, was 1.0 ×10³. The mesh used was 120×120×10 computational cells. Figure 11 shows a qualitative comparison between the experimental results of Ellegaard et al. [37] and the results obtained with our numerical method. One can clearly see from this figure that the numerical method captured, at least qualitatively, the essential physical mechanism (e.g. the circular hydraulic jump and surface waves on subcritical region) of this complex free surface flow.

4 Conclusions

In this paper, an alternative practical upwinding scheme (TOPUS) for advection term discretization has been presented. Several numerical experiments have been performed to verify the accuracy and non-oscillatory shock resolution of this approach to more complicated fluid dynamics PDEs than those presentedby Ferreira et al. [2]. Applications of the method to a number of free surface flow problems with increasing complexity have also been presented.

The main conclusions that can be drawn are:

(i) TOPUS scheme is simple to implement in multidimensional problems. An additional advantage of the scheme is that it produces physical solutions for both hyperbolic and parabolic systems;

(ii) the performance of the TOPUS scheme performed well on different numerical tests, providing good comparisons to experiment, especially considering high Reynolds numbers and complex flow physics; and

(iii) the advantages of the improvement in the TOPUS scheme is apparent - convergent computation and wide applicability for both compressible and incompressible flows.

Acknowledgments. The research reported here has been supported by the FAPESP (grants 09/16954-8, 09/15892-9 and 10/16865-2), CNPq (grants 305447/ 2010-6, 306808/2011-0 and 472945/2011-4) and CAPES (grant PECPG1462/08-3).

Received: 22/IX/12.

Accepted: 28/IX/12.

#CAM-703/12.

• [1] N.P. Waterson and H. Deconinck, Design principles for bounded higher-order convection schemes - a unified approach Journal of Computational Physics, 224 (2007), 182-207.
• [2] V.G. Ferreira, R.A.B. Queiroz, G.A.B. Lima, R.G. Cuenca, C.M. Oishi, J.L.F. Azevedo and S. McKee, A bounded upwinding scheme for computing convection-dominated transport problems Computers & Fluids, 57 (2012), 208-224.
• [3] B.P. Leonard, Simple high-accuracy resolution program for convective modeling of discontinuities International Journal for Numerical Methods in Fluids, 8 (1988), 1291-1318.
• [4] M.A. Alves, P.J. Oliveira and F.T. Pinho, A convergent and universally bounded interpolation scheme for the treatment of advection International Journal for Numerical Methods in Fluids, 41 (2003), 47-75.
• [5] V.G. Ferreira, F.A. Kurokawa, R.A.B. Queiroz, M.K. Kaibara, C.M. Oishi, J.A. Cuminato, A. Castelo, M.F. Tomé and S. McKee, Assessment of a high-orderfinite difference upwind scheme for the simulation of convection-diffusion problems International Journal for Numerical Methods in Fluids, 60 (2009) 1-26.
• [6] P.H. Gaskell and A.K.C. Lau, Curvature-compensated convective transport: SMART, a new boundedness-preserving transport algorithm International Journal for Numerical Methods in Fluids, 8 (1988), 617-641.
• [7] M. Arora and P.L. Roe, A well-behaved TVD limiter for high-resolution calculations of unsteady flow Journal of Computational Physics, 132 (1997), 3-11.
• [8] G.D. van Albada, B. van Leer and W.W. Roberts, A comparative study of computational methods in cosmic gas dynamics Astronomy & Astrophysics, 108 (1982), 76-84.
• [9] B. van Leer, Towards the ultimate conservative difference scheme III. Upstream-centered finite-difference schemes for ideal compressible flow Journal of Computational Physics, 135 (1997), 229-248.
• [10] B. van Leer, Towards the ultimate conservative difference scheme V. A second-order sequel to Godunov's method Journal of Computational Physics, 32 (1979), 101-136.
• [11] B. van Leer, Upwind and high-resolution methods for compressible flow: from donor cell to residual-distribution schemes Communications in Computational Physics, 1 (2006), 192-206.
• [12] B.P. Leonard, Universal limiter for transient interpolation modeling of the advective transport equations: The ULTIMATE conservative difference scheme NASA Technical Memorandum 100916, ICOMP-88-11 (1988).
• [13] W. Jian, P. Traoré and R. Hubert, The extended of convective Boundedness criterion, CP1281, ICNAAM, Numerical Analysis and Applied Mathematics, International Conference 2010, Vol. 1, Ed. T.E. Simons, G. Psihoyios and Ch. Tsitouras.
• [14] V. Venkatakrishnan, On the accuracy of limiters and convergence to steady state solutions AIAA paper 93-0880, (1993), Reno, NV.
• [15] M.J. Kermani, A.G. Gerber and J.M. Stockie, Thermodynamically based moisture prediction using Roe's scheme 4th Conference of Iranian Aerospace Society, Amir Kabir University of Technology, Tehran, Iran, 27-29 (2003).
• [16] R.J. LeVeque, Finite Volume Methods for Hyperbolic Problems Cambridge University Press, New York (2004).
• [17] A.J. Chorin, Numerical solution of the Navier-Stokes equations Mathematics of Computation, 22 (1968), 745-762.
• [18] F.H. Harlow and J.E. Welch, Numerical calculation of time-dependent viscous incompressible flow of fluid with free surface Physics of Fluids, 8 (1965), 2182-2189.
• [19] S. McKee, M.F. Tomé, V.G. Ferreira, J.A. Cuminato, A. Castelo, F.S. Sousa and N. Mangiavacchi, The MAC method Computers & Fluids, 37 (2008), 907-930.
• [20] A. Castelo, M.F. Tomé, C.N.L. Cesar, S. McKee and J.A. Cuminato, Freeflow: an integrated simulation system for three-dimensional free surface flows Journal of Computing and Visualization in Science, 2 (2000), 1-12.
• [21] G. Stecca, A. Siviglia and E.F. Toro, Upwind-biased FORCE schemes withapplications to free-surface shallow flows Journal of Computational Physics, 229 (2010), 6362-6380.
• [22] G. Zhou, L. Davidson and E. Olsson, Transonic inviscid/turbulent airfoil flow simulations using a pressure based method with high order schemes Lecture notes in Physics, Springer-Verlag, 453 (1995), 372-377.
• [23] D.A. Caughey, Diagonal implicit multigrid algorithm for Euler equations AIAA Journal, 26 (1988), 841-851.
• [24] A. Rizzi, Spurious entropy and very accurate solutions to Euler equations AIAA Paper, 84 (1984), 1644-1651.
• [25] A. Jameson and L. Martinelli, Mesh refinement and modeling errors in flow simulation AIAA Journal, 36 (1998), 676-686.
• [26] T.H. Pulliam and J.T. Barton, Euler computations of AGARD working 07 airfoil test cases AIAA Paper, 85-0018, (1985).
• [27] X. Du, C. Corre and A. Lerat, A third-order finite-volume residual-based scheme for the 2D Euler equations on unstructured grids Journal of Computational Physics, 230 (2011), 4201-4215.
• [28] J. Balbás, E. Tadmor and C.C. Wu, Non-oscillatory central schemes for one- and two-dimensional MHD equations: I Journal of Computational Physics, 201 (2004), 261-285.
• [29] E.J. Watson, The radial spread of a liquid jet over a horizontal plane Journal of Fluid Mechanics, 20 (1964), 481-499.
• [30] Y. Brechet and Z. Néda, On the circular hydraulic jump American Journal of Physics, 67 (1999), 723-731.
• [31] T. Bohr, P. Dimon and V. Putkaradze, Shallow-water approach to the circular hydraulic jumps Journal of Fluid Mechanics, 254 (1993), 635-648.
• [32] J.W.M. Bush and J.M. Aristoff, The influence of surface tension on the circular hydraulic jump Journal of Fluid Mechanics, 489 (2003), 229-238.
• [33] S.H. Hansen, S. Horlück, D. Zauner, P. Dimon, C. Ellegard and S.C. Creagh, Geometric orbits of surface waves from a circular hydraulic jump Physical Review E, 55 (1997), 7048.
• [34] J.C. Martin and W.J. Moyce, An experimental study of the collapse of liquid columns on a rigid horizontal plate Philosophical Transactions of Royal Society of London, series A. Mathematical, Physical and Engineering Sciences, 244 (1952), 312-324.
• [35] S. Koshizuka and Y. Oka, Moving-particle semi-implicit method for fragmentation of incompressible fluids Nuclear Science and Engineering, 123 (1996), 421-434.
• [36] A. Colagrossi and M. Landrini, Numerical simulation of interfacial flows by smoothed particle hydrodynamics Journal of Computational Physics, 191 (2003), 448-475.
• [37] C. Ellegaard, A.E. Hansen, A. Haaning and T. Bohr, Experimental results on flow separation and transitions in the circular hydraulic jump Physica Scripta, 67 (1996), 105-110.

Publication Dates

• Publication in this collection
  28 Nov 2012
• Date of issue
  2012

History

• Received
  22 Sept 2012
• Accepted
  28 Sept 2012