Gas Pipe (Fanno)

The gas pipe element of type Fanno is a pipe element with constant cross section (Figure [*]), for which the Fanno formulae are applied [85].

Figure: Geometry of the Gas Pipe element
\begin{figure}\epsfig{file=Gaspipe.eps,width=11cm}\end{figure}

The friction parameter is determined as

f=\frac{64}{\text{Re}} (62)

for laminar flow ( Re<2000) and

\frac{1}{\sqrt{f}} = -2.03 \log \left( \frac{2.51}{\text{Re}\sqrt{f}} +\frac{k_s}{3.7 D} \right). (63)

for turbulent flow. Here, k_s is the diameter of the material grains at the surface of the pipe and Re is the Reynolds number defined by

Re \text{Re} = \frac{U D}{\nu}, (64)

where U is the fluid velocity and \nu is the kinematic viscosity.

It is described by the following parameters (to be specified in that order on the line beneath the *FLUID SECTION,TYPE=GAS PIPE FANNO ADIABATIC or *FLUID SECTION,TYPE=GAS PIPE FANNO ISOTHERMAL card):

The default gas pipe is adiabatic, i.e. there is no heat exchange with the pipe. Alternatively, the user may specify that the pipe element is isothermal. This means that the static temperature does not change within the pipe. In that case the energy equation in one of the two end nodes of the element is replaced by an isothermal condition.

The form factor \varphi is only used to modify the friction expression for non-circular cross sections in the laminar regime as follows:

f=\varphi \frac{64}{\text{Re}}. (65)

Values for \varphi for several cross sections can be found in [13]. For a square cross section its value is 0.88, for a rectangle with a height to width ratio of 2 its value is 0.97.

Instead of specifying a fixed diameter and length, these measures may also be calculated from the actual position of given nodes. The version in which the radius is calculated from the actual distance between two nodes a and b is described by the following parameters (to be specified in that order on the line beneath the *FLUID SECTION,TYPE=GAS PIPE FANNO ADIABATIC FLEXIBLE RADIUS or *FLUID SECTION,TYPE=GAS PIPE FANNO ISOTHERMAL FLEXIBLE RADIUS card):

The version in which the radius is calculated from the actual distance between two nodes a and b and the length from the actual distance between nodes a and c is described by the following parameters (to be specified in that order on the line beneath the *FLUID SECTION,TYPE=GAS PIPE FANNO ADIABATIC FLEXIBLE RADIUS AND LENGTH or *FLUID SECTION,TYPE=GAS PIPE FANNO ISOTHERMAL FLEXIBLE RADIUS AND LENGTH card):

Although a gas pipe looks simple, the equations for compressible flow in a pipe are quite complicated. Here, the derivation is first presented for the adiabatic case. Starting from the energy equation ([*]) and using the relation dh = c_p dT for an ideal gas one arrives at:

c_p dT + v dv = 0. (66)

By means of the definition of the Mach number ([*]) one gets

\frac{dT}{T} = - (\kappa -1) M^2 \frac{dv}{v}. (67)

Because of the ideal gas equation p=\rho r T this can also be written as:

\frac{dp}{p} = -[1+(\kappa -1) M^2] \frac{dv}{v}. (68)

Figure: Differential pipe element
\begin{figure}\epsfig{file=momentumpipe.eps,width=11cm}\end{figure}

Looking at Figure ([*]) the momentum equation can be derived by applying Newton's second law, which states that the sum of the forces is the change of momentum (D is the diameter of the pipe, A its cross section):

A p - A (p + dp) - \tau \pi D dx = \frac {\rho A v dt (v + dv)- \rho A v dt (v)}{dt}, (69)

or, using Darcy's law (\lambda is the Darcy friction factor)

\tau = \frac {\lambda}{4} \frac{\rho}{2} v^2, (70)

\rho \frac{dv^2}{2} + dp + \frac{\lambda}{D}\frac{\rho}{2} v^2 dx =0. (71)

One can remove the density by means of the gas equation arriving at:

v^2 \frac{dv}{v} + r T \frac{dp}{p} + \frac{\lambda}{2} v^2 \frac{dx}{D}=0. (72)

Combining this with what was obtained through the energy equation ([*]) leads to (removing p in this process):

\frac{dv}{v}-\frac{\kappa \lambda}{2} \left( \frac{M^2}{1-M^2} \right )\frac{dx}{D}=0. (73)

This equation contains both M and v. We would like to get an equation with only one of these parameters. To this end the equation defining the Mach number ([*]) is differentiated and the energy equation in the form ([*]) is used to remove T, leading to:

\frac{dM}{M} = \frac{dv}{v} \left( 1 + \frac{1}{2} ( \kappa -1) M^2 \right ). (74)

In that way, the previous equation can be modified its final form:

\frac{dM}{M} = \frac{\kappa M^2}{2(1-M^2)} \left ( 1 + \frac{\kappa -1}{2} M^2\right ) \lambda \frac{dx}{D}, (75)

expressing the Mach number as a function of the distance along the pipe. This equation tells us that for subcritical flow (M < 1) the Mach number increases along the pipe whereas for supercritical flow (M > 1) the Mach number decreases. Consequently, the flows tends to M=1 along the pipe. Notice that by assigning the sign of v to \lambda the above equation also applies to negative velocities. Substituting Z=M^2 and integrating both sides yields:

\int_{Z_1}^{Z_2} \frac{1-Z}{\kappa Z^2} \frac{1}{(1+\frac{\kappa-1}{2}Z)} dZ =\int_0^L \lambda \frac{dx}{D}. (76)

Since (partial fractions)

\frac{1-Z}{Z^2(1+\frac{\kappa-1}{2}Z)}=-\frac{\kappa+1}{2} \frac{1}{Z} +\frac{1}{Z^2} + \left( \frac{\kappa-1}{2} \right )\left( \frac{\kappa+1}{2}\right ) \frac{1}{(1+\frac{\kappa-1}{2} Z)}, (77)

one obtains finally

\frac{1}{\kappa} \left( \frac{1}{M_1^2} - \frac{1}{M_2^2} \right ) +\frac{\kappa +1}{2 \kappa} \ln \left[ \left( \frac{1+\frac{\kappa-1}{2} M_2^2}{1+\frac{\kappa-1}{2} M_1^2} \right) \left( \frac{M_1^2}{M_2^2} \right ) \right] = \lambda \frac{L}{D}, (78)

linking the Mach number M_1 at the start of the pipe with the Mach number M_2 at the end of the pipe, the pipe length L and the Darcy friction coefficient \lambda.

Notice that Equation ([*]) can be used to calculate an estimate of the mass flow due to a given pressure gradient by assuming incompressibility. For an incompressibile medium in a pipe with constant cross section the velocity is constant too (mass conservation) and Equation ([*]) reduces to:

\frac{dp}{\rho} + \frac{\lambda}{2} v^2 \frac{dx}{D}=0. (79)

Integrating yields:

v=\sqrt{\frac{2 D}{\lambda L} \frac{(p_1-p_2)}{\rho}}, (80)

or

\dot{m}=A \sqrt{\frac{2 D}{\lambda L} \rho (p_1-p_2)}, (81)

which can finally also be written as:

\dot{m}=A \sqrt{\frac{2 D}{\lambda L} \frac{{p}_1}{r {T}_1} (p_1-p_2)}. (82)

For an estimate of the mass flow in the gas pipe the above static variables p and T are replaced by the total variables p_t and T_t, respectively. Equation ([*]) is the governing equation for an adiabatic gas pipe. In order to apply the Newton-Raphson procedure (linearization of the equation) with respect to the variables {T_t}_1, {p_t}_1, \dot m, {T_t}_2 and {p_t}_2, this equation is first derived w.r.t M_1 and M_2 (denoting the equation by f; the derivation is slightly tedious but straightforward):

\frac {\partial f}{\partial M_1} = \frac{2}{\kappa M_1} \left [ \frac{M_1^2-1}{M_1^2 ( 1 + b M_1^2)} \right], (83)

and

\frac {\partial f}{\partial M_2} = \frac{2}{\kappa M_2} \left [ \frac{1-M_2^2}{M_2^2 ( 1 + b M_2^2)} \right], (84)

where b=(\kappa-1)/2. Now, M at position 1 and 2 is linked to T_t, p_t and \dot m at the same location through the general gas equation:

\dot m = \frac{a p_t}{\sqrt{T_t }} M (1+b M^2)^c, (85)

where a=A \sqrt{\kappa}/\sqrt{r} and c=-(\kappa+1)/(2(\kappa-1)). Careful differentiation of this equation leads to the surprisingly simple expression:

dM = e \frac{d \dot m}{\dot m} - e \frac{ d p_t}{p_t} + e \frac{d T_t}{2T_t}, (86)

where

e= \left[ \frac{M(1+bM^2)}{1+bM^2(1+2c)} \right]. (87)

Finally, the chain rule leads to the expressions looked for:

\frac{\partial f}{\partial {T_t}_1} = \frac{\partial f}{\partial M_1} \cdot\frac{\partial M_1}{\partial {T_t}_1}, (88)

etc.

For the isothermal pipe the ideal gas equation leads to:

\frac{d p}{p} = \frac{d \rho}{\rho} (89)

and from the mass conservation, Equation ([*]) one gets:

\frac{d p}{p} = - \frac{d v}{v}. (90)

Furthermore, the definition of the Mach number yields:

\frac{d v}{v} = \frac{d M}{M}, (91)

finally leading to:

\frac{d v}{v} = \frac{d M}{M} = -\frac{d p}{p} = -\frac{d \rho}{\rho}. (92)

By substituting these relationships and the definition of the Mach number one can reduce all variables in Equation ([*]) by the Mach number, leading to:

(1-\kappa M^2) \frac{d M}{M^3} = \frac{\kappa \lambda}{2} \frac{dx}{D}. (93)

This equation shows that for an isothermal gas pipe the flow tends to M=1/\sqrt{\kappa} and not to 1 as for the adiabatic pipe. Substituting Z=M^2 and integrating finally yields:

\frac{1}{\kappa} \left( \frac{1}{M_1^2} - \frac{1}{M_2^2} \right ) + \ln \left( \frac{M_1^2}{M_2^2} \right ) = \lambda \frac{L}{D}. (94)

The above equation constitutes the element equation of the isothermal gaspipe. Applying the Newton-Raphson procedure requires the knowledge of the derivatives w.r.t. the basis variables. The procedure is similar as for the adiabatic gas pipe. The derivative of the element equaton w.r.t. M_1 and M_2 is easily obtained (denoting the left side of the above equation by f):

\frac {\partial f}{\partial M_1} = \frac{2}{\kappa M_1^3} (\kappa M_1^2-1) (95)

and

\frac {\partial f}{\partial M_2} = \frac{2}{\kappa M_2^3} (1-\kappa M_2^2). (96)

The use of an isothermal gas pipe element, however, also requires the change of the energy equation. Indeed, in order for the gas pipe to be isothermal heat has to be added or subtracted in one of the end nodes of the element. The calculation of this heat contribution is avoided by replacing the energy equation in the topologically downstream node (or, if this node has a temperature boundary condition, the topologically upstream node) by the requirement that the static temperature in both end nodes of the element has to be the same. This is again a nonlinear equation in the basic unknowns (total temperature and total pressure in the end nodes, mass flow in the middle node) and has to be linearized. In order to find the derivatives one starts from the relationship between static and total temperature:

T_t = T(1+bM^2), (97)

where b=(\kappa-1)/2. Differentiation yields:

dT_t = dT(1+bM^2)+2bTMdM. (98)

Replacing dM by Equation ([*]) finally yields an expression in which dT is expressed as a function of dT_t, dp_t and d \dot m.


Example files: gaspipe10, gaspipe8-cfd-massflow, gaspipe8-oil.