CFD Simulations AC3-02
Induced flow in a T-junction
Application Challenge 3-02 © copyright ERCOFTAC 2004
Overview of CFD Simulations
Simulations of this test case with the code ESTET Version 3.2 (supplied by EDF) have been carried out for a range of the test data. These include a study of the influence of the bulk velocity in the auxiliary pipe as indicated in Table 1, a comparison of two turbulence models and grid sensitivity studies.
The turbulence models considered here are :
- the standard high Reynolds number k-epsilon model by Launder and Spalding (1974)
- a high Reynolds number Reynolds Stress Model (or SMC for « Second Moment Closure ») with isotropic dissipation, turbulent self-transport modeled by the usual gradient diffusion following Daly and Harlow (1970) and the pressure-strain correlation term recommended by Launder (1989), i.e. consisting of the sum of Rotta's return to isotropy and isotropization of production terms.
The grid dependence study has been carried out with meshes of 100 000, 400 000, and 1 500 000 nodes for each turbulence model and each auxiliary mass flow rate considered. Aspect ratios (ratio of the longest to the smallest edge in each cell) have been kept approximately between 1 and 2.
Simulation CFD1 is detailed hereafter. It corresponds to EXP1, with k-epsilon model and a coarse mesh. The 17 remaining simulations differ from CFD1 by the turbulence model, the grid and the auxiliary pipe velocity only : Table 3 CFD B provides a summary.
Remark : with the Reynolds Stress Model all conditions are identical to those used for the k-epsilon model, except for the two differences reported hereafter :
- the inlet conditions for the Reynolds stresses are derived from the k-epsilon inlet conditions assuming isotropic turbulence : the diagonal components are taken equal to 2/3 k whereas the extra-diagonal components are set to zero. The inlet conditions for the other variables remain unchanged.
- the dissipation and Reynolds stresses are cell-centered whereas with the k-epsilon model, k and epsilon are located at the vertices.
Simulation Case CFD1
Solution strategy
The flow investigated being unsteady, the convergence had to be determined from time averages. Hence, once the vortical motion established, time average values of the variables were estimated in the dead leg. The vortex penetration was then determined by visualizing secondary velocities. The calculations were stopped once the value for the mean vortex penetration was stabilized.
Computational Domain
The computational domain represents the geometry described in Figure 1 with the same sharp-edged junction but with different pipe lengths to minimize the size of the numerical model. The inlet of the main pipe has been placed relatively near from the junction to reduce the size of the numerical model. Hence, special attention has been devoted to determination of inlet boundary conditions. The outlet is sufficiently far away to prevent any recirculating structure from crossing it. The height of the auxiliary pipe has been chosen large enough to prevent the corkscrew pattern from reaching its top.
The mesh is single block structured in cartesian orthogonal coordinates. Curved boundaries are represented by special features of the code which enables the treatment of slanted boundary cells (prisms, tetrahedrons...).
Figure 3 : Computational domain. © ERCOFTAC 2004 Boundary Conditions
Default initial conditions implemented in the code have been used for each simulation.
As for the boundary conditions at the inlets (main pipe and auxiliary pipe), Dirichlet conditions are used for all transported variables. In the main pipe, the incoming flow is supposed to be fully developed. Hence, the mean velocity and turbulent quantities were obtained from preliminary periodic pipe computations with the same bulk velocity and the same cross-section. For the auxiliary pipe, due to the very small flow-rates, it is not necessary to describe so accurately the incoming flow. The mean velocity was taken constant and turbulent quantities were estimated from usual experimental laws for the friction velocity. In fully developed pipe flows with zero-roughness, and for Reynolds numbers Re (based on the bulk velocity VA and the hydraulic diameter D) ranging between 5,000 and 30,000, the friction velocity ut can be determined from ut = VA ((0.3164/Re0.25)/8)0.5. The turbulent kinetic energy k can then be estimated from ut as k = ut2/0.3 and the dissipation as ε = ut3/(0.42 D 0.1).
At the outlet (main pipe), zero gradient conditions were applied for all transported variables.
Wall functions were used at solid walls.
The velocity in the main pipe is 9.2m/s, given in Table 1. The auxiliary pipe velocity can be found in Table 3 CFD B. © ERCOFTAC 2004 Application of Physical Models
The standard high Reynolds number k-epsilon model proposed by Launder and Spalding (1974) has been used for CFD1. See Table 3 CFD B for the other simulations. © ERCOFTAC 2004 Numerical Accuracy
The numerical techniques are based on finite difference and finite volume discretizations. The equations are discretized on a 3D semi-staggered grid. Velocity is located at the vertices and pressure is cell-centered. When the k-epsilon model is used, the turbulent variables are located at the vertices (whereas they are cell-centered with the Reynolds stress model).
The algorithm for the solution of transient Navier-Stokes equations relies on a segregated velocity-pressure formulation.The advective terms are treated by a method of characteristics. The trajectory is approximated by a second order Runge Kutta scheme with a third order 3D Hermitian polynomial for interpolation. Diffusion with explicit and implicit source terms for dynamic variables is solved implicitly. For the computation of the velocity components, a third step is required in order to prescribe the mass conservation, leading to a Poisson equation for the pressure increment. In order to avoid non physical oscillations of the pressure field and the associated difficulties in obtaining a converged solution, a variant of the Rhie and Chow interpolation is used as in Méchitoua et al., (1994). The overall precision of the discretization is first order in time and second order in space. © ERCOFTAC 2004 CFD Results
Figure 4 compares the length of vortex penetration obtained by numerical simulation with experimental measurements.
NAME
GNDPs
PDPs (problem definition parameters)
SPs (simulated parameters)
Sensitivity analysis
Inlet Re
Main pipe
Bulk velocity VM (m/s)
Auxiliary pipe
Bulk velocity VA (m/s)
Velocity ratio VA / VM
DOAPs
Grid (number of nodes)
Turbulence model
CFD 1
895 000.
9.2
0.092
1%
vortex penetration
100 000
k-epsilon
CFD 2
895 000.
9.2
0.092
1%
vortex penetration
400 000
k-epsilon
CFD 3
895 000.
9.2
0.092
1%
vortex penetration
1 500 000
k-epsilon
CFD 4
895 000.
9.2
0.046
0.5%
vortex penetration
100 000
k-epsilon
CFD 5
895 000.
9.2
0.046
0.5%
vortex penetration
400 000
k-epsilon
CFD 6
895 000.
9.2
0.046
0.5%
vortex penetration
1 500 000
k-epsilon
CFD 7
895 000.
9.2
0.023
0.25%
vortex penetration
100 000
k-epsilon
CFD 8
895 000.
9.2
0.023
0.25%
vortex penetration
400 000
k-epsilon
CFD 9
895 000.
9.2
0.023
0.25%
vortex penetration
1 500 000
k-epsilon
CFD 10
895 000.
9.2
0.092
1%
vortex penetration
100 000
SMC
CFD 11
895 000.
9.2
0.092
1%
vortex penetration
400 000
SMC
CFD 12
895 000.
9.2
0.092
1%
vortex penetration
1 500 000
SMC
CFD 13
895 000.
9.2
0.046
0.5%
vortex penetration
100 000
SMC
CFD 14
895 000.
9.2
0.046
0.5%
vortex penetration
400 000
SMC
CFD 15
895 000.
9.2
0.046
0.5%
vortex penetration
1 500 000
SMC
CFD 16
895 000.
9.2
0.023
0.25%
vortex penetration
100 000
SMC
CFD 17
895 000.
9.2
0.023
0.25%
vortex penetration
400 000
SMC
CFD 18
895 000.
9.2
0.023
0.25%
vortex penetration
1 500 000
SMC
Table 3 CFD-B Summary description of all test cases © copyright ERCOFTAC 2004
Contributors: Frederic Archambeau - EDF - R&D Division
Site Design and Implementation: Atkins and UniS
Top Next