Elsevier

Journal of Computational Physics

Volume 284, 1 March 2015, Pages 230-260
Journal of Computational Physics

A new velocity–vorticity formulation for direct numerical simulation of 3D transitional and turbulent flows

https://doi.org/10.1016/j.jcp.2014.12.030Get rights and content

Abstract

Accuracy of velocity–vorticity (V,ω)-formulations over other formulations in solving Navier–Stokes equation has been established in recent times. However, the issue of non-satisfaction of solenoidality conditions on vorticity is not addressed in the literature which can possibly lead to non-physical solution. In this respect, here, we have developed and reported conservative rotational form of the (V,ω)-formulation which preserves the solenoidality condition on vorticity in a much simpler way compared to other formulations. Superiority of rotational form over the conventional Laplacian form of (V,ω)-formulation is also shown [by comparing the results for flows inside cubical lid driven cavity (LDC)]. For solving the 3D Navier–Stokes equation using a staggered grid, we use optimized compact schemes for (a) interpolation and (b) evaluation of first and second derivatives. As illustrations, we have solved problems of (i) flow inside a 3D lid driven cavity (LDC), whose solutions are compared with experimental results reported by Koseff and Street [19] and (ii) 3D transitional flow of an equilibrium zero pressure gradient (ZPG) boundary layer over a flat plate as demonstration of the effectiveness of the rotational form of (V,ω)-formulation.

Introduction

Direct numerical simulations (DNS) and large eddy simulations (LES) of incompressible flows are generally attempted using either primitive or derived variable formulations. Fractional step methods with three time level Adams–Bashforth scheme for time integration are generally used in primitive variable formulation [1] for time-accurate solution of Navier–Stokes equation (NSE). It has been shown in [2] that use of three time level Adams–Bashforth algorithm introduces spurious numerical mode, which removes high frequency components of the initial solution, invalidating any claim of performing DNS by this approach. Another difficulty with this approach is the specification of wall boundary conditions at the intermediate stage. This is usually achieved via extrapolation from values obtained in previous time steps. Perot [3] has shown that while this fix may provide numerical stability, but it gives rise to large time splitting error. To minimize this error, one is again forced to take small time steps. An approach that avoids this problem by developing correct conditions is given in [4], but is very memory intensive [5]. In [6], fourth-order accurate explicit and compact finite-difference methods for incompressible NSE in (p,V)-formulation is presented where four-stage, fourth-order Runge–Kutta (RK4) time integration scheme is used. In this formulation, the solenoidality constraint on velocity is enforced for each Runge–Kutta stage by solving respective Poisson equations for pressure correction and is computationally very intensive.

Thus, it is natural to seek formulations which eliminate pressure altogether from the governing equations. In 2D flows, we note the advantage of stream function-vorticity (ψ,ω)-formulation over the primitive variable formulation, for reasons of lesser number of unknowns as compared to three variables in (p,V)-formulation; exact satisfaction of mass conservation everywhere in the flow field and obtaining vorticity directly from the governing equation and not by numerical differentiation of computed quantities. However, this formulation's extension to 3D flows by using vector potential-vorticity formulation is not straightforward due to complex and coupled boundary conditions. In view of this, (V,ω)-formulation is more amenable, as we show the equivalent governing equation for Dv=V in (V,ω)-formulation to be entirely different, which allows for a simple satisfaction of Dv=0 condition, as compared to (p,V)-formulation.

(V,ω)-formulations have been used by several researchers to solve external [7] and internal flow problems [8]. Different references report different forms of (V,ω)-formulation. Mostly, two forms of this formulation have been used: (i) a non-conservative form [9], [10], [11] and (i) the Laplacian form [8], [12]. The difference between these two formulations lie in the way nonlinear convective and vortex-stretching terms in the vorticity transport equations are handled. In non-conservative form, the vorticity transport equation (VTE) is given asωt+(V)ω=(ω)V+1Re2ω

Note that the first term on the right hand side of Eq. (1) is the vortex stretching term, which is absent for 2D flows. In contrast, the Laplacian form of VTE is given asωt+×(ω×V)=1Re2ω

Some of the references on (V,ω)-formulations mentioned above, also discuss about the solenoidality conditions on velocity and vorticity. For example, Laplacian form of (V,ω)-formulation is used in [13], where the solenoidality condition on vorticity is enforced based on evolution equation for divergence of vorticity (Dω). In this reference, the authors also claim that their formulation implicitly satisfies solenoidality condition on velocity due to the way the governing equations are derived. Laplacian form has also been used in [12], [14]. In [14] solenoidality of vorticity is satisfied by integrating the wall-normal component of vorticity (ωy) from free-stream to the wall for flows over a rotating circular disk for which no evolution equation is solved. However, such procedure is only applicable for the external flows for which the component of vorticities at the free-stream are zero. In [12], no such attempt is made, although the computations are claimed to exactly satisfy divergence-free conditions for both velocity and vorticity. Non-conservative forms of the (V,ω)-formulation are used in [10], [11], [12]. Evolution equation for Dω is reported in [8], where schemes are constructed to satisfy the solenoidality condition on ω. In [15], the computed vorticity field ω is made divergence-free by adding to it the gradient of a scalar function ϕ˜ asωtr=ω+ϕ˜ where ωtr is the unknown divergence-free vorticity vector. Taking divergence of Eq. (3), one obtains the governing equation for the scalar ϕ˜ as2ϕ˜=ω Hence in [15], four Poisson equations are solved at each time-step, which makes the procedure computationally more expensive and time consuming. While adopting such procedure of obtaining divergence-free vorticity field, one encounters the problem of specifying the boundary condition for ϕ˜ which is addressed in [8]. However, all the 3D LDC results reported in [8], [12], [15], [16] are only for the subcritical Reynolds number for which divergence error in vorticity does not affect the final flow structures. This aspect is also shown in the present work via results obtained for the cubical LDC.

In the present work, we report a new rotational form of (V,ω)-formulation in conservative form. We have also derived the evolution equations for divergence of vorticity for the rotational, Laplacian and non-conservative forms of this formulation and shown that enforcing solenoidality condition on vorticity is much simpler for rotational form compared to the other two forms. In the present work, we have used optimized compact schemes for spatial discretization, which have been shown in [17] to have superior spectral accuracy compared to traditional schemes and specially central difference schemes. Also, an optimized Runge–Kutta scheme [18] is used for time integration which provides high accuracy solution to the problems solved here. This formulation is used to solve flow inside 3D rectangular LDC comparing the obtained results with experimental results reported in [19]. Results obtained for cubical LDC are also used to compare the Laplacian and rotational forms. The method used in [14] to solve the disturbance evolution for flows over a rotating disc is applicable only for external flows, as Eqs. (12) to (14) in this reference can be applied to external flows only. On the other hand, our method reported in the present work is general and can be applied to both external and internal flows.

DNS of 3D NSE for the receptivity of flows past a flat plate using Laplacian form of (V,ω)-formulation are also reported in [7], to simulate the evolution of disturbances in 3D flow field. However, this references does not comment about the effects or satisfaction of solenoidality conditions on vorticity. The enforcement of wall-boundary condition on vorticity in [7] is not obtained in the conventional way of specifying Ω=×V as given by Eqs. (12a)–(12c) in this reference.

The other test example considered in the present work involves simulation of transition of an equilibrium boundary layer excited by a Gaussian bump and spanwise modulated wall-exciter in the 3D framework by DNS using (V,ω)-formulation in rotational form. Corresponding DNS of 2D Navier–Stokes equation have been reported in [20], [21]. A numerical investigation of velocity and vorticity field in boundary layer transition process are reported in [22]. In [23], the boundary layer was shown to support growing spatio-temporal wave fronts, which have been shown to grow and cause transition to turbulence in [20]. The deterministic disturbances were shown from the receptivity stage to the fully developed turbulent stage by solving the 2D Navier–Stokes equation without any limiting assumptions. An attempt is made to see that if this is also the case for 3D transitional flow. Solution of such flow in 3D framework has been reported in [24].

The governing equations for (V,ω)-formulation and the evolution equations for Dv and Dω have been developed in the next section. Computational details for: (i) flow inside 3D rectangular LDC and (ii) 3D transition of an equilibrium boundary layer over a flat plate caused by a Gaussian bump and spanwise modulated exciter have been provided in this section. In Section 3, results and discussion for these two test cases have been provided. The paper closes with summary and conclusions in Section 4.

Section snippets

Governing equations and numerical formulation

The governing equations for rotational form of VTE are described next. Additionally, the associated velocity Poisson equations are also derived. In the following, numerical methods adopted to solve the problem are discussed in details.

Computational details, boundary conditions, results and discussion

To show the efficacy of the proposed rotational formulation, using the high accuracy compact schemes of interpolation and derivative evaluation, here we have studied test problems, which simulate flows from onset to transitional and turbulent states exhibiting typical power spectrum.

Summary and conclusions

Here, we have derived a new (V,ω)-formulation of Navier–Stokes equation in rotational form of the VTE. This formulation allows simpler satisfaction of the solenoidality condition on vorticity, which makes it superior to the traditionally used primitive variable formulations or the Laplacian and the non-conservative form of the (V,ω)-formulation. In all the computations reported, solenoidality of the velocity is also preserved by exactly satisfying the divergence-free condition on velocity.

References (48)

  • C. Davies et al.

    A novel velocity–vorticity formulation of the Navier–Stokes equations with applications to boundary layer disturbance evolution

    J. Comput. Phys.

    (2001)
  • M.K. Rajpoot et al.

    Solution of linearized rotating shallow water equations by compact schemes with different grid-staggering strategies

    J. Comput. Phys.

    (2012)
  • T.K. Sengupta et al.

    Space–time discretizing optimal DRP schemes for flow and wave propagation problems

    Comput. Fluids

    (2011)
  • L.-P. Wang et al.

    A spurious evolution of turbulence originated from round-off error in pseudo-spectral simulation

    Comput. Fluids

    (2009)
  • S. Nagarajan et al.

    A robust high-order compact method for large eddy simulation

    J. Comput. Phys.

    (2003)
  • P.C. Chu et al.

    A three-point combined compact difference scheme

    J. Comput. Phys.

    (1998)
  • T.K. Sengupta et al.

    Further improvement and analysis of CCD scheme: dissipation discretization and de-aliasing properties

    J. Comput. Phys.

    (2009)
  • S.K. Lele

    Compact finite difference schemes with spectral-like resolution

    J. Comput. Phys.

    (1992)
  • T.K. Sengupta et al.

    A new compact scheme for parallel computing using domain decomposition

    J. Comput. Phys.

    (2007)
  • T.K. Sengupta et al.

    A new compact difference scheme for second derivative in non-uniform grid expressed in self-adjoint form

    J. Comput. Phys.

    (2011)
  • T.B. Gatski

    Review of incompressible fluid flow computations using the vorticity–velocity formulation

    Appl. Numer. Math.

    (1991)
  • T.K. Sengupta et al.

    A comparative study of time-advancement methods for solving Navier–Stokes equation

    J. Sci. Comput.

    (2004)
  • L. Kleiser et al.
  • D.C. Lo et al.

    Numerical solution of three-dimensional velocity–vorticity Navier–Stokes equations by finite difference method

    Int. J. Numer. Methods Fluids

    (2005)
  • Cited by (28)

    • Combined effects of amplitude, frequency and bandwidth on wavepackets in laminar turbulent transition

      2020, Computers and Fluids
      Citation Excerpt :

      A three-wave system known as the Craik triad [15] was discovered to play a major role [14,45] in the nonlinear evolution of wavepackets and these triads are closely associated with N-type transition, which is also known as H-type transition [21,22]. More recently, the mechanism of genesis and evolution of the spatio-temporal wavepacket in a linear framework was explored in [42,43], and the influence of frequency in determining either K-type or H-type transition in wavepackets was described in [3–5,44], who used the synonymous term “spatio-temporal wave front” (STWF) for a wavepacket. This spatio-temporal wavepacket is a modulated and band-limited combination of spatio-temporal eigenmodes of the Orr-Sommerfeld equation as described by Eq. (1) in [3].

    • A single domain velocity–vorticity fast multipole boundary domain element method for two dimensional incompressible fluid flow problems

      2019, Engineering Analysis with Boundary Elements
      Citation Excerpt :

      As representative works someone could mention the [35], that uses a second-order vorticity kinematic boundary condition for the no-slip wall boundary condition, that satisfies continuity, or the work of [3], where the velocity field is corrected with the solution of an additional equation in order to enforce the continuity. In the recent work of 3D velocity–vorticity formulation implemented in [2], optimized compact schemes have been used for the interpolation and the evaluation of the first and second derivatives. In several BEM velocity–vorticity formulations however like in [39], the kinematics equation that satisfies continuity, provides directly the boundary conditions for the vorticity equation, rendering the method into a stable and robust CFD tool.

    • Global spectral analysis of multi-level time integration schemes: Numerical properties for error analysis

      2017, Applied Mathematics and Computation
      Citation Excerpt :

      Taking first order, forward Euler time stepping method turns out to be numerically unstable and for this reason use of four-stage two-time level Runge–Kutta (RK4) schemes [31] are used for higher accuracy. For example, optimal three- and four-stage, time integration methods along with compact schemes for spatial discretization have been proposed in OCRK3 and ORK4 methods, which ensures dispersion relation preservation (DRP) properties [2,36]. The motivation for using this governing equation in an infinite domain is that one does not require any boundary condition and the initial condition is propagated to the right.

    View all citing articles on Scopus
    View full text