ES_CriFlow  Compressible Flow Analysis Engine of Steam, Rev. 1 









4.1 Nozzle 














Steam shows different behaviors in compressible flow from those of ideal gas because of its inconsistent characteristics depending on the conditions. In this concern, different approaches have been used in the analyses of compressible flow in power plant engineering, depending on their conditions of steam, i.e. superheated steam, saturated steam and saturated water.
This category of analysis in power plant engineering includes cascade heater drain analysis, safety valve vent stack analysis and steam blowout analysis. There are many reference papers for these subjects as listed in the References at the end of this paper. Ref. No. 1 through 3 are for cascade heater drain analysis; Ref. No. 4 through 6 are for safety valve vent stack analysis; and Ref. No. 7 is for steam blowout analysis.
The reference papers use basically the equations for ideal gas with various experimental coefficients and experimental equations in order to supplement the inadequacy of ideal gas laws and equations for the compressible flow analysis of steam. And they use different experimental equations and coefficients depending on the conditions of steam. In case of Ref. NO 1, for example, a pseudoisentropic exponent is used to apply Fannoline equation to the twophase flow of saturated water. Also they use several fundamental laws additionally to check whether the analysis results are correct or not.
Anyway, the target of the analysis method presented by these papers is to provide results with positive margins, and the system so designed should have conservative margins in any case.
ENGSoft Inc. has developed a core program, ES_CriFlow, for the compressible flow analysis of steam, using steam table, the basic equations of flow dynamics and the speed of current desktop computers. The outputs of ES_CriFlow runs had been compared with the data of Reference papers, and presented in a separate page titled as "ES_StmNzl and ES_StmPipe Program Test Results". Through the comparison, ES_CriFlow has been proved as effective for the compressible flow analysis of steam regardless of the conditions of steam. Since the ES_CriFlow outputs had not been proved by actual experimental test, the term of "effective" has been used. Meanwhile, it was found that the ES_CriFlow run outputs are well in line with the measurements performed in a actual cascade heater drain system and presented in Ref. NO. 3, and therefore ENGSoft Inc. wants to say that the ES_CriFlow run output is correct.
The program "StmNzl" and "StmPipe" have been developed for testing the effectiveness of the method and its algorithm of ES_CriFlow. Using ES_CriFlow, ENGSoft Inc. developed also several engineering programs of ES_CasDrain for cascade heater drain line analysis, ES_SVVent for safety valve vent stack analysis, and ES_BlowOut for steam blowout analysis.
The page presents the basic equation of flow dynamics for compressible flow analysis and the equation of ideal gas flow analysis, and then finally presents the method and algorithm used in ES_CriFlow.
2. Basic Equations of Flow Dynamics (TOC)
Basic equations which applies to both ideal gas and steam are as follows.
 Energy Equation (The Fist Law of Thermodynamics)
 Continuity Equation
 Equation of Motion
The equations presented here are for onedimensional, steady flow.
The general energy equation in infinitesimal form is as below.
q + w = du + d(p * v) + d(V^2 / 2/ g) + d(h)
where, 
q 
: Heat added per mass of flowing fluid 

w 
: Work added per mass of flowing fluid 

u 
: Internal Energy 

p 
: Static Pressure 

v 
: Specific Volume 

V 
: Fluid Velocity 

h 
: Elevation Head 
Since the process of nozzle and short length pipe can be assumed as a adiabatic process with no external work, q and w become zero in the equation above. According to the definition of enthalpy, du + d(pv) = dH, and in case of gas and steam the elevation head(dh) can be neglected.
Therefore, the energy equation to be applied to nozzle and short length pipe is as below.
d(H) + d(V^2 / 2 / g / J) = 0
where, 
H 
: Enthalpy, kcal/kg 

V 
: Fluid velocity, m/sec 

g 
: Gravity acceleration = 9.81 m/sec2 

J 
: Joule constant, 427 kgm/kcal 
The above equation means that the fluid velocity is generated by the consumption of enthalpy and the fluid velocity is another form of energy the enthalpy transformed. This equation is effective as long as the process is adiabatic even though there is friction in the process.
The enthalpy at zero velocity is called as total enthalpy which is the maximum value of enthalpy in the process, and the total enthalpy is calculated as below.
Ht = H + V^2 / 2 / g / J
where, 
Ht 
: Total Enthalpy, kcal/kg 
If the precess is adiabatic, then the total enthalpy value through the process is constant.
The continuity equation below is applied to all fluid flow process.
W / A = V / v
where, 
W 
: Mass flow rate, kg/sec 

v 
: Specific volume, m3/kg 

V 
: Fluid velocity, m/sec 

A 
: Flow path area, m2 
The equation of motion is one of the energy equation which represents the energy balance among pressure head, velocity head, elevation head and friction head. The momentum equation says the same thing with the equation of motion, but in force terms instead of head terms in the equation of motion.
If we neglect the elevation head for gas and steam flow, the momentum equation in force terms is described as below in infinitesimal forms.
A * dP = T * 3.14 * D * dL + W * dV / g
where, 
A 
: Flow path area 

dP 
: Pressure difference 

T 
: Wall Shear Stress 

dL 
: Duct length 

D 
: Hydraulic diameter of duct 

g 
: Gravity acceleration 

W 
: Mass flow rate 

dV 
: Velocity difference 
The first term in the equation represents the force for flow by pressure difference, the second represents the resistance against flow by friction, and the third represents the momentum change by velocity due to specific volume change.
The equation means that the force by pressure difference is used to overcome the friction and to make momentum change.
In the equation, the friction term written in wall shear stress can be replaced with the term in the friction factor commercially used, i.e. T = f * V^2 / v / 8. Under the assumption that the duct flow path area is constant along the duct length, when dividing the equation by area(A or 3.14 * D^2 / 4) and multiplying by specific volume(v), we can write the equation as below.
v * dP = f * dL / D * V^2 / 2 / g + W / A * v * dV / g
(?) 
A 
: Flow path area, m2 

dP 
: Pressure difference, kg/m2 

f 
: Friction Factor 

dL 
: Duct length, m 

D 
: Hydraulic diameter of duct, m 

g 
: Gravity acceleration = 9.81 m/sec2 

dV 
: Velocity difference, m/sec 

W 
: Mass flow rate, kg/sec 

v 
: Specific volume, m3/kg 
In the compressible flow analysis, critical pressure should be calculated in order to find out the possibility of choked flow, and in order to calculate the critical pressure the calculation of sonic velocity is required.
For compressible fluid, the sonic velocity is calculated as below.
Vc = (dP / dRo * g)^0.5 @ isentropic infinitesimal pressure change
where, 
Vc 
: Sonic Velocity, m/sec 

dP 
: Infinitesimal pressure difference, kg/m2 

dRo 
: Density change by dP in isentropic process, kg/m3 

g 
: Gravity acceleration, 9.81 m/sec2 
Critical pressure exists in compressible fluid flow, where the fluid velocity equals the sonic velocity of the fluid.
Before the downstream pressure reaches the critical pressure, the compressible mass flow rate increases as much as the downstream pressure is decreased. However, once the downstream pressure reaches or is less than the critical pressure, the compressible mass flow rate does not increase even though the downstream pressure is decreased further.
This specific phenomenon is called as "choked flow", and the energy difference between the choked exit and ambient conditions is dissipated by shock wave and/or turbulence.
For the case of compressible nozzle flow, the pressure, mass flow rate and velocity have the relationship described below.
Type of Flow 
Pressure 
Mass Flow Rate 
Velocity 
Subcritical Flow 
P2 = P3 > Pc 
F < Fc 
V2 < Vc 
Critical Flow 
P2 = P3 = Pc 
F = Fc 
V2 = Vc 
Choked Flow 
P2 = Pc > P3 
F = Fc 
V2 = Vc 
where, 
Pc 
: Critical Pressure 

Fc 
: Critical Flow Rate 

Vc 
: Critical Velocity = Sonic Velocity 
If a convergentdivergent de Laval nozzle is used, then a supersonic flow can be made. However, such a supersonic flow is normally not used in power plant engineering.
The first thing in analysis of compressible flow of nozzle and piping is to calculate the critical pressure(Pc) and then to check whether the discharge pressure(P3) is higher than the critical pressure. If P3 >= Pc, the nozzle or piping exit pressure(P2) is selected as P3, and if P3 < Pc then P2 is selected as Pc.
In case of orifice, actually there is no critical pressure. The mass flow rate increases as much as the discharge pressure is decreased till zero absolute pressure. This is because, as the fluid velocity increases, the friction generated at the edge of the orifice is absorbed into the flowing fluid itself and increase the fluid temperature, the mass flow rate per unit area decreases and more mass flow can pass the orifice. However, the increase rate of mass flow rate is very slow when the discharge pressure is less than approximately 50% of the inlet pressure which is the critical pressure ratio of the nozzle of general gas.
4. Compressible Flow Analysis of Ideal Gas (TOC)
4.1 Nozzle
4.1.1 Mass Flow Rate per Unit Area
The mass flow rate per unit area of nozzle expanding isentropically (P * v^k = Const.), is calculated by the following equation, which is derived from the energy equation, continuity equation and BoyleCharles equation (P * v = R * T).
W / A = ((2 * g * k) / (k  1))^(0.5) * (P1 / v1)^(0.5) * (r^(2/k)  r^((k+1)/k))^(0.5)
where, 
r 
: Pressure ratio = P2 / P1 

P1 
: Nozzle inlet pressure, kg/m2 abs. 

v1 
: Nozzle inlet specific volume, m3/kg 

W 
: Mass flow rate, kg/sec 

A 
: Nozzle throat area, m2 

g 
: Gravity acceleration = 9.81 m/sec2 
4.1.2 Critical Pressure and Nozzle Throat Pressure
The critical pressure of ideal gas expanding isentropically in nozzle, is calculated by the following equation, which is derived from the equation above where the mass flow rate per unit area has the maximum value.
Pc = (2 / (k + 1))^(k / (k  1)) * P1
where, 
Pc 
: Critical pressure, kg/m2 abs. 

k 
: Specific heat ratio of ideal gas = cp/cv 

P1 
: Nozzle inlet total pressure at V1 = 0, kg/m2 abs. 
The nozzle exit pressure(P2) is selected considering the discharge pressure(P3) according to the method described in the Clause 3 above.
As described above, there is no choked flow in orifice flow, and the mass flow rate is calculated by the following equation, which is applicable for turbulent flow.
Y = 1  0.41 * (P1  P2) / P1 / k
W / A = 0.598 * Y * (2 * g * (P1  P2) / v1)^(0.5)
4.3 Adiabatic Pipe with Friction (TOC)
The difference of pipe analysis from that of nozzle is the involvement of friction and heat loss. The existence of friction means that the process is not isentropic and the existence of heat loss means that the process is not adiabatic.
In pipe flow analysis, friction is always there. However, the heat loss is not.. When analyzing long distance lines such as cross country gas pipe lines, the heat loss consideration is must. However, when analyzing relatively short pipe lines, the heat loss can be neglected. Almost all pipe lines in power plant engineering is short enough to be analyzed without the consideration of heat loss.
The equation for ideal gas flowing through constant crosssectional area duct with friction and heat loss is called as Rayleigh Line equation, while the equation for ideal gas flowing through constant area duct with friction but in adiabatic process, is called as Fanno Line equation. Herein, the Fanno Line equation is introduced for its applicability to power plant engineering.
The assumption given to Fanno Line equation is as below.
Assumption of Fanno Line Equation 
Ideal gas (constant specific heat) 
Steady, onedimensional flow 
Constant friction factor over the length of duct 
Adiabatic flow (no heat transfer through wall) 
Effective conduit diameter D is four times hydraulic radius (crosssectional area divided by wetted perimeter) 
Elevation changes are unimportant compared with friction effects 
No work added to or extracted from the flow 
The Fanno Line equation is for the flow which has sonic velocity(Mach no. = 1) at the duct(or pipe) exit plane, i.e. choked flow.
4.3.1 Choked Flow
As depicted in the diagram below for choked flow, the pipe exit pressure(P2) is same with the critical pressure(Pc) and the fluid velocity at pipe exit is same with sonic velocity, i.e. Mach no. = 1. And, the discharge pressure is equal or less than the critical pressure.
The Fanno Line equation defines the relationship among pipe inlet Mach no.(M1), flow resistance coefficient(K) and specific heat ratio of the ideal gas(k). The Fanno Line equation written below is derived from energy equation, continuity equation, equation of motion, BoyleCharles equation (P * v = R * T) and sonic velocity equation (Vc = (k * g * R * T)^(0.5)).
K = 1 / k * (1 / M1 ^ 2  1) + (k + 1) / 2 / k * ln(M1 ^ 2 * (k + 1) / ((k  1) * M1 ^ 2 + 2))
where, 
K 
: Flow resistance coefficient, K = f * L / D 

k 
: Specific heat ratio = cp/cv 

M1 
: Pipe inlet Mach no. 
The equations for critical ratios of pressure, temperature and velocity are as blow, which are shown in term of pipe inlet Mach no.(M1).
Pc / P1 = M1 * (((k  1) * M1 ^ 2 + 2) / (k + 1))^(0.5)
Vc / V1 = 1 / M1 * (((k  1) * M1 ^ 2 + 2) / (k + 1))^(0.5)
Tc / T1 = ((k  1) * M1 ^ 2 + 2) / (k + 1)
where, 
Pc 
: Critical pressure 

Vc 
: Critical velocity = Sonic velocity 

Tc 
: Critical temperature 
When the flow resistance coefficient(K) and specific heat ratio(k) are known, the pipe inlet Mach no.(M1) can be calculated by the Fanno Line equation above. With M1 calculated, the critical ratios can be calculated. However, in order to get the conditions of pipe inlet and outlet from the critical ratios, either of pipe inlet or outlet condition should be known.
Furthermore, since the Fanno Line equation is based on choked flow, the critical pressure at the pipe exit should be known in order to select the pipe exit pressure. Therefore, the Fanno Line equation with the critical ratio equations is not enough to analyze the compressible pipe flow of ideal gas.
In case of compressible nozzle flow, since the flow conditions along the nozzle are determined by isentropic process, knowing of the nozzle inlet condition with discharge pressure is good enough for analysis. However, in case of the adiabatic pipe flow with friction, the flow conditions is indeterminate and vary depending on the mass flow rate per unit area(W/A) and pipe resistance coefficient(K).
Therefore, in order to analyze the adiabatic pipe flow with friction the mass flow rate per unit area and pipe resistance coefficient should be known in addition to pipe inlet condition and pipe discharge pressure.
Pipe inlet condition here means the total condition at zero velocity. Choked flow of adiabatic pipe with friction is depicted below.
Critical pressure of adiabatic pipe flow with friction for ideal gas can be derived as below.
The critical temperature ratio equation above may be rewritten as below by substituting T1 with T0 where the velocity equals zero and M0 = 0.
Tc / T0 = 2 / (k + 1)
Applying the BoyleCharles equation (P * v = R * T), the equation above can be written,
(Pc * vc) / (P0 * v0) = 2 / (k + 1)
Pc = 2 / (k + 1) * (P0 * v0) / vc (Eq. 4.3.2  1)
P0 and v0 are the known, while vc can be calculated from the continuity equation (W / A = Vc / vc)and sonic velocity equation (Vc = (k * g * P * v)^(0.5)), wherein the mass flow rate(W) and pipe area(A) are constant through the pipe length.
Two equations can be reduced for vc as below by eliminating Vc.
vc = k * g * Pc / (W / A)^2 (Eq. 4.3.2  2)
Substituting vc with (Eq. 4.3.22), (Eq. 4.3.21) can be written,
Pc = (W / A) * (2 * P0 * v0 / k / (k + 1) / g)^(0.5)
The critical pressure calculation equation above is applicable to all adiabatic pipes, whatever processes the flow go through in the pipe.
4.3.3 Subcritical Flow
Subcritical pipe flow has the pipe exit pressure(P2) same with discharge pressure(P3) and subsonic velocity at the pipe exit. The subcritical pipe flow can be analyzed using a imaginary pipe between P2 and Pc as depicted below.
The pressure ratio of the imaginary pipe which is in choked flow, is as below by Fanno Line equation.
Pc / P2 = M2 * (((k  1) * M2 ^ 2 + 2) / (k + 1))^(0.5)
Since the Pc and P2 are known, the above equation can be rearranged for M2 as below.
M2 = (((1 + (k  1) * (k + 1) * (Pc/P2)^(2))^(0.5)  1) / (k  1))^(0.5)
With M2 known, the resistance coefficient of the imaginary pipe(Ki) can be calculated by Fanno Line equation as below.
Ki = 1 / k * (1 / M2 ^ 2  1) + (k + 1) / 2 / k * ln(M2 ^ 2 * (k + 1) / ((k  1) * M2 ^ 2 + 2))
Since total pipe including actual pipe(K) and imaginary pipe(Ki) is in choked flow again, the pipe inlet condition can be calculated according to the method described in Clause.4.3.1.
The general sonic velocity equation described in Clause 2.4 may be rewritten in isothermal process as below.
Vc = (dP / dRo * g)^0.5 @ isentropic infinitesimal pressure change
= (k * dP / dRo * g)^(0.5) @ isothermal infinitesimal pressure change
Integrating the isothermal equation by using BoyleCharles equation (P * v = R * T), we can get the sonic velocity equation of ideal gas as below.
Vc = (k * g * R * T)^(0.5) = (k * g * P * v)^(0.5)
where, 
Vc 
: Sonic velocity, m/sec 


k 
: Specific heat ratio, cp/cv 


R 
: Gas constant, R = R0 / M



g 
: Gravity acceleration, 9.81 m/sec2 


T 
: Gas absolute temperature, degree Kelvin 


P 
: Gas absolute pressure, kg/m2 abs. 


v 
: Gas specific volume, m3/kg 
4.5 Compressible Flow through Increaser and Reducer (TOC)
When fluid flows through increaser or reducer, there is flow resistance by vortex in addition to that by friction due to wall shear stress.
When fluid flows through a increaser having the crosssectional area changed suddenly, the pressure just behind of main stream exiting from the smaller crosssectional area is lower than that of main stream by vortex, which disturbs the free expansion of main stream and acts as flow resistance.
When analyzing the compressible flow through increaser and reducer, the flow resistance by vortex should be considered, even though the crosssectional area changes suddenly and so there is no friction by wall shear stress.
Meanwhile, even if there is no friction by wall shear stress and there is only flow resistance by vortex, the flow process is not insentropic process, i.e. not reversible adiabatic process, but a kind of polytropic process without external work having the polytropic exponent normally higher than isentropic exponent(k). This is true even though the flow process is adiabatic process. Reversible adiabatic process has no increase in entropy. But adiabatic process has the increase of entropy, even though there is no heat transfer through boundary. Of course, the adiabatic process includes the reversible adiabatic process.
In order to understand the flow mechanism and K value calculation of increaser and reducer, the theory of incompressible flow is presented first.
4.5.1 Analysis in Incompressible Flow (TOC)
A. Increaser
In a increaser having crosssectional area suddenly enlarged as depicted above, Bernoulli equation between the location 1 and 2 can be written as below neglecting the elevation head for gas. Bernoulli equation is the equation of motion for incompressible flow. The flow resistance by vortex is called as hf in the equation below.
P1 * v + Vel1^2 / 2 / g = P2 * v + Vel2^2 / 2 / g + hf
Rearranging the equation for hf,
hf = (Vel1^2  Vel2^2) / 2 / g  (P2  P1) * v (Eq. 4.5.1  1)
Equation of momentum between location 1 and 2 is as below.
(P1  P2) * A2 = W / g * (Vel2  Vel1)
Using continuity equation (W = A2 * V2 / v), the equation above may be arranged as below.
(P2  P1) * v = Vel2 * (Vel1  Vel2) / g (Eq. 4.5.1  2)
Eliminating pressure terms from (Eq. 4.5.1  1] and (Eq. 4.5.2  2) and rearranging for hf, we may have,
hf = (Vel1  Vel2)^2 / 2 / g = (1  A1 / A2)^2 * Vel1^2 / 2 / g (Eq. 4.5.1  3)
The equation above for hf is theoretical equation. However, it was revealed by experiments that the equation represents the actual value very well for the increaser having the crosssectional area enlarged suddenly.
Comparing (Eq. 4.5.1  3) with Darcy equation, the value corresponding to the resistance coefficient(K) is (1  A1 / A2)^2 which is based on the velocity of location 1. Conversion to the velocity of location 2 can be done by dividing the value by (A1 / A2)^2.
In case of the increaser having the crosssectional area enlarged gradually, there is the friction by wall shear stress in addition to the flow resistance due to vortex, because the fluid flows in contact with the increaser wall. While there are many books providing with the experimental equations for the resistance coefficient of increaser, herein the experimental equations provided in Ref. No. 8 is introduced.
If ? <= 45 o (based on Aup velocity)
K = 2.6 * Sin(?/2) * (1  Aup / Adown)^2
if 45 o < ? <= 180 o (based on Aup velocity)
K = (1  Aup / Adown)^2
We can see in the equations above that only the flow resistance by vortex exists if the increaser angle is larger than 45o and the value is same with the theoretical value. If the increaser angle is smaller than 45o, the flow resistance is increased due to additional friction by wall shear stress.
B. Reducer (TOC)
The flow pattern of reducer having the crosssectional area contracted suddenly is as depicted below.
Due to fluid inertia the crosssectional area of flow path is contracted smaller than the physical area of the smaller end of the reducer passing over the location, and then enlarged to the physical area of downstream pipe. The flow resistance by the contraction itself is so small to be negligible, and the flow resistance of reducer mainly consists of the flow resistance by enlargement as in increaser, plus the friction by wall shear stress in case of the graduallycontracted reducer.
The experimental equations for the resistance coefficient of reducer presented in Ref. No. 8 is as below.
If ? <= 45 o (based on Adown velocity)
K = 0.8 * Sin(?/2) * (1  Adown / Aup)
if 45 o < ? <= 180 o (based on Adown velocity)
K = 0.5 * Sqr(Sin(?/2)) * (1  Adown / Aup)
Since reducer has always both flow resistance by enlargement and the friction by wall shear stress except the reducing angle of 180o, the equations above include both also. Meanwhile, the flow resistance only by enlargement can be calculated as below, which is the case of 180o reducing angle.
K = 0.5 * (1  Adown / Aup)
4.5.2 Analysis of Ideal Gas
In order to apply the resistance coefficient for incompressible flow to compressible flow, the specific volume change due to pressure change should be considered.
If the flow process is isentropic process, then the analysis can be done straight using several equations. However, since the process of increaser and reducer is not isentropic, the analysis should be done by tryanderror method using the basic equations.
The ideal gas equations which can be used for irreversible adiabatic process are as below.
1) BoyleCharles equation
P * v = R * T
2) Continuity equation
W / A = V / v
3) Energy equation
H0 = H + V^2 / 2 / g
4) Equation of Motion
v * dP + (V1^2  V2^2) / 2 / g  hf = 0
5) Enthalpy equation
dH = cp * dT
6) Constant pressure specific heat equation
cp = k * R / J / (k  1)
6) Entropy equation
dS = cp * Ln(T2 / T1)  R / J * Ln(P2 / P1)
For enthalpy and entropy calculation of ideal gas, the normal ambient condition(0 oC and 10332 kgf/m2 abs.) shall be used as reference condition.
In the analysis, it is assumed that the downstream condition of increaser or reducer is known by previous analysis of downstream pipe, and the upstream condition of increaser or reducer shall be calculated using tryanderror method.
It is noted that the upstreamend flow velocity of increaser is higher than that of downstreamend, and so the upstreamend pressure of increaser is lower than that of downstreamend. Meantime, the lower limit of the upstreamend pressure of increaser is the critical pressure at the upstreamend.
Viceversa, the upstreamend pressure of reducer is higher than that of downstreamend, and the higher limit is the inlet pressure of total pipe. The upstreamend pressure of reducer can not be lowered below the critical pressure in any case.
Hereinafter, the upstream end of increaser and reducer will be called as Loc 1 and the downstream end as Loc 2.
A. Calculation Method of Increaser
Since the lower limit of Loc 1 pressure is the critical pressure, the pressure of Loc 1 is the value between Loc 2 pressure and the critical pressure. Loc 2 pressure is the higher limit and the critical pressure is the lower limit.
If Loc 2 pressure is equal or lower than the critical pressure, Loc 1 pressure is the critical pressure and the other condition of Loc 1 is the critical condition.
The Loc 1 enthalpy is the value between the total enthalpy of the fluid(the higher limit) and the insentropic expansion enthalpy at Loc 1 pressure(the lower limit).
The condition of Loc 1 is selected by tryanderror method searching the Loc 1 pressure and enthalpy ranges described above, which meets the equation of motion.
B. Calculation Method of Reducer
Calculation method of reducer is same with that of increaser, except that the pressure range to search is different.
5. ES_CriFlow  Compressible Flow Analysis of Steam (TOC)
Herein the tryand error method and algorithm of ES_CriFlow for the analysis of compressible flow of steam is presented.
The key difference of steam analysis from that of ideal gas is the use of steam table instead of BoyleCharles equation, enthalpy and entropy equation, and the use of tryanderror method to get steam conditions from steam table for each iteration.
5.1.1 Input Data
The input data required for the analysis of nozzle are as below.
 
Stagnated condition of nozzle inlet where steam velocity equals zero. Two properties are required for defining the steam condition, e.g. pressure(P1) and enthalpy(H1) 
 
Nozzle discharge pressure (P3) 
5.1.2 Critical Pressure, Nozzle Throat Pressure and Mass Flow Rate per Unit Area
Nozzle analysis of steam is very simple because the process can be analyzed as isentropic process.
First, get from steam table the nozzle inlet entropy(S1).
Select a nozzle throat pressure(P2) lower than the nozzle inlet pressure(P1), and get steam enthalpy(H2) at nozzle throat from steam table using P2 and S1. The condition is the isentropic expansion condition.
Then calculate the sonic velocity at the nozzle throat condition selected, and then the velocity energy by the sonic velocity.
If the sonic velocity energy is different from the enthalpy difference of H1  H2, then select another nozzle throat pressure(P2) and try again.
If the sonic velocity energy converges on the enthalpy difference of H1  H2, then calculate the mass flow rate per unit area at nozzle throat using the sonic velocity and specific volume selected and finish the calculation.
5.2 Adiabatic Pipe with Friction (TOC)
Nozzle analysis is rather simpler than pipe analysis because the process can be analyzed using isentropic process, as described above. However, the adiabatic pipe with friction is different because the process is polytripic of which process is dependent on the given conditions of pipe. Therefore, the mass flow rate per unit area and pipe resistance coefficient shoul be known for analyzing the adiabatic pipe with friction.
There are two kinds of compressible steam pipe analyses in power plant engineering.
The First is to calculate pipe inlet and outlet conditions from the stagnated pipe inlet condition, pipe crosssectional area, mass flow rate and pipe discharge pressure known. These kinds of analyses include the safety valve vent stack analysis and cascade heater drain pipe analysis.
The second is to calculate the maximum mass flow rate in addition to the pipe inlet and outlet conditions from the stagnated pipe inlet condition, pipe crosssectional area and pipe discharge pressure given. This kind of analysis includes the steam blowout pipe analysis.
Actually the second analysis is the iteration of the first analysis by tryanderror method, in which the mass flow rate is searched for the pipe inlet pressure converging on the stagnated pipe inlet pressure given.
In this Clause, the method of the first analysis is described.
5.2.1 Input Data
The input data required for the analysis of adiabatic pipe with friction are as below.
 
Stagnated condition of pipe inlet where steam velocity equals zero. Two properties are required for defining the steam condition, e.g. pressure(P0) and enthalpy(H0) 
 
Pipe discharge pressure(P3) 
 
Pipe crosssectional area(A),(to be constant) 
 
Mass flow rate(W) 
 
Resistance coefficient of pipe K = f * L / D 
5.2.2 Critical Pressure(Pc)
The critical pressure is calculated using the facts that the velocity at critical condition is sonic velocity and that the sum of the static enthalpy and the velocity energy at critical condition is same with the stagnated pipe inlet enthalpy.
In the calculation, it should be noted that the sonic velocity is calculated by using insentropic infinitesimal pressure change, even if the pipe process is not isentropic.
The critical pressure must exist in the pressure range below the stagnated pipe inlet pressure. If the critical pressure does not exists below the stagnated pipe inlet pressure(P0), this means that the mass flow rate of input data can not flow through the pipe given even at choked flow condition and the pipe condition given by the input data does not exist.
For each critical pressure selected for tryanderror, the maximum enthalpy is the stagnated pipe inlet enthalpy(H0) and the minimum enthalpy is the enthalpy expanded through isentropic process to the critical pressure selected.
Consequently, the critical pressure is selected by tryanderror method by searching the pressure range below the stagnated pipe inlet pressure and searching the maximum and minimum enthalpies at each pressure described above, in which the mass flow rate calculated by the sonic velocity and the specific volume is converging on the mass flow rate of input data.
5.2.3 Pipe Exit Condition(Location 2)
The pipe exit pressure(P2) is selected by comparing the critical pressure(P2) and the pipe discharge pressure(P3). Selection method is described in Clause 3 above.
The maximum value of the pipe exit enthalpy(H2) is the stagnated pipe inlet enthalpy(H0) and the minimum value is the isentropically expanded enthalpy to pressure P2 selected. The pipe exit enthalpy(H2) is selected by tryanderror method searching between the maximum and minimum values, in which the mass flow rate calculated by the velocity generated by the enthalpy difference converges on the mass flow rate of input data. Other properties of steam at pipe exit can be gotten from steam table using P2 and H2.
5.2.4 Pipe Inlet Condition(Location 1)
The pipe inlet condition is calculated using the momentum equation. The momentum equation for a pipe with friction means that the force by pressure difference between the inlet and outlet equals the flow resistance fore by pipe wall friction plus the momentum increase by velocity increase.
The momentum equation in differential form is as below.
v * dP = f * dL / D * V^2 / 2 / g + W / A * v * dV / g
where, 
A 
: Pipe crosssectional area, m2 

dP 
: Pressure difference, kg/m2 

f 
: Friction factor 

dL 
: Pipe length, m 

D 
: Pipe diameter , m 

g 
: Gravity acceleration = 9.81 m/sec2 

dV 
: Velocity difference, m/sec 

W 
: Mass flow rate, kg/sec 

v 
: Specific volume, m3/kg 
Rewriting the friction term of the above equation using the square of the continuity equation (V^2 = (W / A)^2 * v^2), we have,
v * dP = f * dL / D * (W / A)^2 * v^2 / 2 / g + W / A * v * dV / g
Dividing the equation by v^2, then we have
(1 / v) * dP = f * dL / D * (W / A)^2 / 2 / g + W / A * (1 / v) * dV / g (Eq. 5.2.4  1)
Meanwhile, the differential form of continuity equation is as below.
dV = (W / A) * dv (W / A = constant) (Eq. 5.2.4  2)
Substituting dV of (Eq. 5.2.4  1) by dV of (Eq. 5.2.4  2), integrating the equation (Eq. 5.2.4  1) from pipe inlet to pipe outlet, and then substituting f by K using K = f * L / D, we have,
K = { Integral(dP / v)(from P2 to P1) } / (W/A)^(2) * 2 * g  2 * ln(v2/v1) (Eq. 5.2.4  3)
where, 
K 
: Pipe resistance coefficient, K = f * L / D 

dP 
: Incremental pressure difference, kg/m2 

P1 
: Pipe inlet pressure, kg/m2 abs. 

P2 
: Pipe exit pressure, kg/m2 abs. 

W 
: Mass flow rate, kg/sec 

A 
: Pipe crosssectional area, m2 

g 
: Gravity acceleration = 9.81 m/sec2 

v1 
: Pipe inlet specific volume, m3/kg 

v2 
: Pipe exit specific volume, m3/kg 
Fanno Line equation for ideal gas is the equation that the pressure integral term of the above equation has been solved using the BoyleCharles equation. However, steam is not simple because the relationship of pressure and specific volume can not be expressed by a equation.
For steam, the integral can be solved by summing up the reciprocal of specific volume for incremental pressure changes with reasonable accuracy. For specific volume, the algebraic average value of the incremental pressures is used.
The pipe inlet pressure(P1) is selected as the pressure at which the resistance coefficient calculated by (Eq. 5.2.4  3) equals the resistance coefficient of input data. The summingup starts from the pipe exit pressure(P2) to the stagnated pipe inlet pressure(P0). Other steam properties at pipe inlet are calculated as described for the pipe exit condition in Clause 5.2.3 above..
Meanwhile, if the resistance coefficient which meets the value of input data can not be found till summing up to P0, that means the pipe condition given by the input data does not exist.
5.2.6 Sensitivities of Program Variables
As described in Clause 1 Introduction, the methods described in Clause 5 have been programmed and then the sensitivities of major variables have been investigated and summarized as below.
1) 
The sonic velocities calculated by either of 100 kg/m2 or 10 kg/m2 infinitesimal pressure changes show no distinct difference, and 100 kg/m2 pressure change has been selected for use in the program. 
2) 
It was found in the integration work of the momentum equation(Eq. 5.2.4  3) that the incremental pressure change by ratio basis is appropriate rather than by algebraic addition. In the program runs it was found also that the pressure increases by either of 1% or 0.1% show no distinct difference, and 1% has been selected for use in the program. 
3) 
The sonic velocity calculation was tried for subcooled water, but the result was not effective because the specific volume change for infinitesimal pressure change got from steam table is too small to make the result meaningful. However, the sonic velocity of flashing saturated water was found effective. 
5.3 Compressible Flow through Increaser and Reducer (TOC)
The analysis method of steam flow through increaser and reducer is same with that of ideal gas described in Clause 4.5 above, except that the steam table is used instead of BoyleCharles equation, enthalpy and entropy equation.
In Ref. No. 7, it was described that the downstream pipe size should be less than the upstream pipe size without any explanation. The truth is that the method given in the paper is not appropriate for the increaser analysis. It does not mean that compressible steam does not flow through the increaser. Using the method described here, the compressible steam flow analysis through increaser can be solved without any problem.
References :
1. Analytical Approach for Determination of Steam/Water Flow Capability in Power Plant Drain Systems by G.S. Liao and J.K. Larson, ASME Publication 76WA/Pwr4
2. Heater Drain Systems by A.L. Cahn, Bechtel Power Corp., presented for Feedwater Heater Workshop of EPRI held in July, 1979.
3. Flow of a Flashing Mixture of Water and Steam through Pipes by M.W. Benjamin and J.G. Miller, Transactions of the ASME, Oct., 1942
4. ASME B31.11992, Appendix II Nonmandatory Rules for The Design of Safety Valve Installations
5. Analysis of Power Plant Safety and Relief Valve Vent Stacks by G.S. Liao, Bechtel Power Corp., Transactions of the ASME, 1974
6. Crosby Pressure Relief Valves Engineering Handbook, Crosby Gage & Valve Company, March 1986
7. Cleaning of Main Steam Piping and Provisions for Hydrostatic Testing of Reheaters (GEK  27065D)", General Electric Co.
8. Crane Technical Paper No. 410, Flow of Fluids, Crane Co., 1977
9. Principles and Practice of Flow Meter Engineering by L. K. Spink, Foxboro
Copyright (c) 2000  2001 ENGSoft Inc., Seoul, Korea, All right reserved.