NBSIR 77-1394

A Computer Code for the

Simulation of Turbulent

Swirling Flow Through

Closed-Conduit Flow Meters

R. W. Davis and E. F. Moore

Fluid Meters Section

Mechanics Division

Institute for Basic Standards

National Bureau of Standards

Washington, D.C. 20234

December 1977

Final

U S. DEPARTMENT OF COMMERCE

NATIONAL BUREAU OF STANDARDS

NBSIR 77-1394

A COMPUTER CODE FOR THE

SIMULATION OF TURBULENT

SWIRLING FLOW THROUGH

CLOSED-CONDUIT FLOW METERS

R. W. Davis and E. F. Moore

Fluid Meters Section

Mechanics Division

Institute for Basic Standards

National Bureau of Standards

Washington, D.C. 20234

December 1977

Final

U.S. DEPARTMENT OF COMMERCE, Juanita M. Kreps, Secretary

Dr. Sidney Harman, Under Secretary

Jordan J. Baruch, Assistant Secretary for Science and Technology

NATIONAL BUREAU OF STANDARDS, Ernest Ambler, Acting Director

ABSTRACT

A computer code is presented for the simulation of turbulent swirling

flow through closed-conduit flow meters. The code is reasonably fast and

accurate. The response of a flow meter to changes in its operating

environment is easily assessed. Through the use of computer plotting routine

a complete picture of the flowfield through a meter can be obtained. An

example flowfield for a target meter is presented.

1 . INTRODUCTION

A computer code currently in use at NBS has proven highly useful in

determining the characteristics of closed-conduit flow meters operating

under realistic conditions [1]. This code, called TEACH, has been developed

at the Department of Mechanical Engineering, Imperial College of Science

and Technology, London, England. It is in wide use throughout the United

States to model such turbulent flows as those found in combustors and in

nuclear reactors [2, 3]. The TEACH code has generally been found to give

reasonably accurate results within reasonable computing times. It was thus

the logical candidate to pick in testing the feasibility of using turbulence

modeling to study the flow through various types of closed-conduit flow

meters. The results have demonstrated this feasibility, and computer simulation

is now regarded as a vital tool in NBS' effort to improve flow metering

practice

The TEACH code can be used to obtain steady solutions for axisymmetric

flows with or without swirl, a swirl subroutine having been added at NBS. The

ke-model is employed to determine the turbulent shear stresses [4]. The resulting

solution is in terms of velocities, pressure, and turbulence intensity. These

can be computer plotted in such a way as to give a detailed view of the flowfield

through a metering device. Changes in this flowfield due to variations in such

parameters as upstream velocity profile, Reynolds number, or swirl number are

easily assessed. A target meter flowfield is presented as an example.

2. BASIC EQUATIONS

The momentum and continuity equations for turbulent incompressible

stationary flow employed in the ker-model are

p(q * V)q = - VP - Vx[p

eff

(Vxq) ] (1)

V • q = 0 ( 2 )
in which q = (U, V, W) , where U, V, and W are the mean velocity components in

the radial, axial, and tangential directions, respectively, in an axisymmetric

cylindrical reference frame; P is the mean pressure; p is the constant density;

and p cc is the sum of the laminar and turbulent viscosities. The turbulent

etr

viscosity, p , is assumed to be a function of the turbulent kinetic energy,

k, and dissipation rate, e, of k; i.e.,

y

t = C

y

p k

2 /e (3)

where the empirical constant C is 0.09. The values of k and e are determined

y

by the following transport equations written in axisymmetric cylindrical

coordinates (r, z) : 3k . „ 3k, 9 /eff 9k, , 1 9 ,ry eff 9k, ^ „ „ „ P ( v + u

= 77 — Tr) + - — (-— T~) + vG - C PE (4) 9z 9r 9z a, 9z r 9r a, 9r t D k k

/lT 9e , TT 9e, 9 , M eff 3e, , 1 3 , ry eff 9e, p(v + ^

H (—^ + 7 37^ +

e e y

t e

2 C

l

e k

G C

2

P k (5)

where the empirical constants are a, = 1.0, o = 1.3, C = 1.0, C = 1.44, and

K. 0 U JL r 1 no. , r _ 0 r ^

/U . 2 , ,3V. 2, , r 9V 3U.2 .9W.2 .9 W. . 2 C„ - 1.92; and G = 2 (—) + (-) + (—) ] + [-r- + -r~] + (—) + rr- (-) ] . I dr r 9z dr dz dz dr r The origin of the values for the empirical constants is discussed in [4].

3. NUMERICAL SCHEME

The TEACH code employs the SIMPLE finite difference scheme. This scheme is

described in [5, 6], while the code itself is discussed in [2, 7]. Only a brief

review of TEACH will thus be given here.

The TEACH code employs a variably-spaced, rectangular grid. The radial cell

dimension is only a function of radial position, and the axial cell dimension

depends only on axial position. Normal velocities are specified along cell

boundaries, and pressure, tangential velocity, turbulent kinetic energy, and

dissipation rate are defined at cell centers. The system of equations

(1) - (5) is finite differenced and solved iteratively by underrelaxation. The

2
nonlinear convective terms in these equations are differenced using a hybrid

scheme comprised of central and conservative upwind differencing [8]. At

each iterative step, the finite difference analogs of equation (1) are solved

implicitly using a line-by-line procedure involving a tridiagonal matrix solver.

The velocity field resulting from this solution of the momentum equations for a

given iteration will probably not satisfy the continuity equation (2). Using

this known (but probably incorrect) velocity field, equation (2) can be transformed

into an equation for the pressure correction required to make the velocity field

solenoidal. This pressure correction is calculated in the same implicit

manner and is used to determine the final velocity field for a given iteration.

These velocities are then used in the solution (by the usual implicit scheme)

of equations (4) and (5) for k and e. Finally, equation (3) for y^_ is updated

and the next iteration begins. The iterations proceed until the sum

of the absolute values of the residual mass and momentum sources is less than

some preset small value. This summation is performed separately for the mass

flow and each component of momentum, with the resulting sums to each be less

than some small percentage of the mass flow and largest component of

momentum, respectively, at the inlet to the computed region.

Boundary conditions must be specified at the inlet to the region of

flow containing the metering device. The required conditions are radial

distributions of U, V, W, k, and e. The geometry of the flow meter itself

is handled by not allowing fluid to cross the appropriate grid lines. The

condition at the exit of the computed region is that axial gradients are zero.

Shear stresses along walls are based on the logarithmic law; i.e., wallfunction method of [4]. Axisymmetry is enforced at the centerline (r = 0) of

the computed section.

3
4. COMPUTER PLOTTING

The plotting program makes extensive use of four CalComp (California

Computer Products) subroutines. Three of these, PLOT, SYMBOL, and NUMBER,

are part of CalComp 's basic software package that adapts program output

for plotting. The fourth subroutine, FLINE, is classified as a "functional"

subroutine, in which basic subroutines are applied to more specified needs.

Both basic and functional software is available from CalComp.

The PLOT subroutine converts the output data from TEACH (U, V, W, P, k, e)

into hardware codes for a pen plotter. SYMBOL enables printing of alphanumeric

data, such as titles for plots. This subroutine also creates the symbols used

to plot the data points. NUMBER is used to provide axis value labels on plots.

The FLINE subroutine has the capability of plotting an entire array of data

points with symbols, and then drawing straight lines or smooth curves between

successive data points. A parabolic curve-fit is generally used at NBS

.

The plotting program generally runs in less than three minutes and uses

about 40K words of storage on the NBS UNIVAC 1108. About 1.5 hours is usually

required to make the graphs from the resulting output. When a complete set

of plots is generated, approximately 9 meters of 30 cm paper is used. Some

examples of these plots for a target meter are presented in the next section.

5. EXAMPLE: TARGET METER

As an example of what can be done with TEACH and the plotting program,

a target meter flowfield is presented in Fig. 1. The target meter is merely

a disc of zero thickness (T = 0) suspended in the center of a circular pipe at a

distance of four pipe radii from the inlet. The beta ratio (ratio of disc radius

to pipe radius) is 0.5, and the Reynolds number based on pipe radius and mean

inlet velocity is 2 r 10“*. The inlet conditions are based on fully-developed pipe

4
flow at this Reynolds number, with the axial velocity profile conforming to a

standard l/9-power law. A solid-body rotation is superimposed on these

inlet conditions in order to assess the effects of swirl on the meter. The

inlet swirl number, which represents the ratio of axial flux of angular

momentum to axial flux of axial momentum, is 0.10. All quantities are

nondimensional , with the characteristic length being pipe radius and the

characteristic velocity being mean inlet velocity. There are 100 variablyspaced grid points axially and 20 variably-spaced ones radially. It required

about 400 iterations to reach the final solution, which took about 40 minutes

and 55K words of storage on the NBS UNIVAC 1108.

Figure (la) shows pressure and axial velocity variations through the meter.

The pressure difference across the disc can be integrated to give the force

on the disc, which is related to flow rate. Figure (lb) gives wall shear

stress (x) and swirl number (S) anywhere in the meter. Figure (lc) is a

streakline (segmented streamline) plot, with each line segment being in the

local flow direction at a given point. The recirculation zone along the pipe

centerline behind the disc shows clearly, as does the concentration of grid

points near the disc. Figure (Id) presents turbulence intensity based on local

mean velocity, I, along the pipe wall and pipe centerline. The flattened peak

along the centerline results from cutting off a spike which occurs at the end

of the recirculation zone where the local mean velocity goes to zero (I -* °°) .

Figure (le) shows axial velocity profiles at several critical points in the meter.

Note the power-law inlet profile. The outlet profile has not yet returned to

the fully-developed case. Figure (If) presents tangential velocity profiles at

the same locations as for the axial ones. The inlet profile is linear,

denoting solid-body rotation.

5
6.

CONCLUSION

The feasibility of performing computer simulations of turbulent swirling

flow through closed-conduit flow meters has been demonstrated. A turbulence

code, such as the TEACH code presented here, can, along with proper plotting

routines, define entire flowfields quickly and cheaply. It is thus expected

that this sort of simulation, in conjunction with well-chosen validation

experiments, will become more common in the determination of the performance

of flow meters and associated pipeline elements, such as upstream flow conditioners.

7 . ACKNOWLEDGEMENT

The authors are grateful to Dr. G. E. Mattingly for his assistance and

encouragement in the task of setting up these simulations at NBS

.

8. REFERENCES

1. Davis, R. W. and Mattingly, G. E., Numerical Modeling of Turbulent

Flow through Thin Orifice Plates, Flow Measurement in Open Channels

and Closed Conduits , ed. by Lafayette K. Irwin, National Bureau of

Standards Special Publication 484, p. 491, October 1977.

2. Lilley, David G., Primitive Pressure-Velocity Code for the Computation

of Strongly Swirling Flows, AIAA Journal , Vol. 14, p. 749, June 1976.

3. Chen, Yi Bin and Golay, Michael W. , Coolant Mixing in the LMFBR Outlet

Plenum , COO-2245-44TR, Dept, of Nuclear Energy, MIT, June 1977.

4. Launder, B. E. and Spalding, D. B., The Numerical Computation of

Turbulent Flows, Comp. Meth. in Appl. Mech. and Egg. , Vol. 3, p. 269, 1974.

5. Caretto, L. S., Gosman, A. D. , Patankar, S. V., and Spalding, D. B.

,

Two Calculation Procedures for Steady, Three-Dimensional Flows with

Recirculation, Proc. Third Int . Conf. on Nuro. Meth. in Fluid Mech. ,

ed. by Henri Cabannes and Roger Temam, Vol. II, Spr inger-Verlag , New York,

p. 60, 1973.

6
6. Patankar, S. V., Numerical Prediction of Three-Dimensional Flows,

Studies in Convection , Vol. 1, ed . by B. E. Launder, Academic Press,

New York, p. 1, 1975.

7. Salimonu, Dauda 0., Mathematical and Programming Details of the TEACH

T and L Recirculating Flow Code , ANL-CT-76-22 , Argonne National

Laboratory, Oct. 1975.

8. Spalding, D. B., A Novel Finite Difference Formulation for Differential

Expressions Involving Both First and Second Derivatives, Int . Journal

for Num. Meth. in Eng. , Vol. 4, p. 551, 1972.

7
2 TOUGCT MCTC* 0ET«=O.S. RE=2«10«»5. T =0

TRftOCT K€TC* SETA=0.5. *C=2»10»»S. TsO

(a)

(b)

TflBOCT HE Tft 8ETO = 0S. H€=2»lO«»5. T=0

TAHOET Hi I EH »tT« = 0.fc. *C*2»lO»«b. TmO

(c)

(d)

TARGET PIETER BETArO.S. REi2»10»«S.

V I M INLET . X UPSTREAM OP PLATE. DOWNSTREAM OF PLATE)

(e)

TARGET PIETER BETArO -5 . RE=2»tO«»S.

TrO

PLATE. DOWNSTREAM OF PLATE)

TARGET PIETER BETArO.S. RE = 2"tO«»5.

O REATTACMMENT . + OUTLET)

(f)

Fig. 1. Target meter flowfield.

8
as-l 14A (REV 7-73)

U S. DEPT. OF COMM.

BIBLIOGRAPHIC DATA

SHEET

1. PUBLICATION OR REPORT NO.

NBSIR 77-1394

2. Gov’t Accession

No.

3. Recipient’s Accession No.

TITLE AND SUBTITLE

A COMPUTER CODE FOR THE SIMULATION OF TURBULENT

SWIRLING FLOW THROUGH CLOSED-CONDUIT FLOW METERS

5. Publication Date December 1977

6. Performing Organization Code

AUTHOR(S)

R. W. Davis and E. F. Moore

8. Performing Organ. Report No.

. PERFORMING ORGANIZATION NAME AND ADDRESS

NATIONAL BUREAU OF STANDARDS

DEPARTMENT OF COMMERCE

WASHINGTON, D.C. 20234

10. Project/Task/Work Unit No.

11. Contract/Grant No.

2. Sponsoring Organization Name and Complete Address (Street , City, State, ZIP)

National Bureau of Standards

Department of Commerce

Washington, D.C. 20234

13. Type of Report & Period

Covered

Final

14. Sponsoring Agency Code

5. SUPPLEMENTARY NOTES

&. ABSTRACT (A 200-word or less (actual summary of most significant information. If document includes a significant

bibliography or literature survey, mention it here.)

A computer code is presented for the simulation of turbulent swirling flow

through closed-conduit flow meters. The code is reasonably fast and

accurate. The response of a flow meter to changes in its operating environment

is easily assessed. Through the use of computer plotting routines, a complete

picture of the flowfield through a meter can be obtained. An example

flowfield for a target meter is presented.

7. KEY WORDS (six to twelve entries; alphabetical order; capitalize only the first letter of the first key word unless a proper name; separated by semicolons)

Flow meter; meter performance; numerical modeling; turbulent flow.

19. SECURITY CLASS 21. NO. OF PAGES

(THIS REPORT)

UNCL ASSIFIED 10

20. SECURITY CLASS 22. Price

(THIS PAGE)

UNCLASSIFIED $4.00

8. AVAILABILITY

I.

1 X Unlimited

For Official Distribution. Do Not Release to NTIS

c J Order From Sup. of Doc., U.S. Government Printing Office Washington, D.C. 20402, SD Cat. No. Cl 3

! X Order From National Technical Information Service (NTIS)

Springfield, Virginia 22151 USCOMM-DC 29042- P 74

I