how to find two-phase mass flow rate (not mass flux) 1

[CJC0117]

Hi all, I am trying to figure out how to size a PRV for two-phase flow. All the sizing calculations I've come across give methods like HEM or the Omega method in order to find the maximum two-phase mass flux, G. Then the formula A=W/(K*G) (where A is required orifice area, W is the mass flow rate, K is the discharge coefficient, and G is the mass flux) is supposed to be used to find the required orifice area. My problem is that I cannot figure out the standard procedure for calculating the mass flow rate W in the case of two-phase flow. I found this thread,

http://www.eng-tips.com/viewthread.cfm?qid=4187,

but I didn't really understand the explanation given in the first reply. I looked at the calculations that Digitial Solutions Tech. uses in their Mach II program. It seems that they use W=Q/L (where W is the mass flow rate, Q is the heat input, and L is the latent heat of vaporization) to find the mass flow rate. However, I am hesitant to use this formula because I thought it applied only to the all-vapor case. Could someone please clarify this for me, or else pexplain to me the standard method, if any, for finding W? I know that API 520 says there is not "any recognized procedure for certifying the capacity of PRVs in two-phase flow service." How can the orifice areas be calculated then? Thanks.

 

[apetri]

you got some good information in a very recent thread,
there is a procedure for HEM discussed here

http://eng-tips.com/viewthread.cfm?qid=329093

The model for the nozzle applied in the procedure discussed above is

hin+1/2*vin^2 = ho+1/2*vo^2

(in = inlet, o = orifice)
where vo (for a critical flow) is the speed of sound
ho, vo calculated at vena contracta conditions

the vena contracta conditions are determined with a series of constant entropy flash operations,
to solve the constant entropy flash (i.e. given initial t,p and final p (or p) the procedure finds tout (or pout) which gives Sout=Sin).

There are several models for two phase flow,
I utilize Prode Properties which includes

1) HEM Homogeneous Equilibrium (Solution of Mass Flux integral)
2) HNE Homogeneous Non-equilibrium (HEM with Boling Delay and Gas-Liquid Slip Contributes)
3) HNE-DS , Homogeneous Non-equilibrium
4) NHNE Non-homogeneous Non-equilibrium

and every model offers some advantages in specific areas of application.

To my knowledge there are not simple models for two-phase,
even API Omega requires to solve a entropy (or enthalpy) flash operation,
from your posts it seems that you wish to write a code,
in my opinion you need a very good knowledge of thermodynamics to complete this task,
or you may use a commercial product for solving entropy/enthalpy/volume as suggested in previous thread.

 

[Latexman]

Are you working on a specific scenario, like a fire case, or on a general procedure?

Good luck,
Latexman

Need help writing a question or understanding a reply? forum1529

 

[CJC0117]

I would preferably like to formulate a general procedure, but if you know a way for finding W in the fire case, I'd be glad to hear it, as that's often the worst case emergency. The information in the last thread I started indeed was very helpful. I've since been able to complete my excel sheet that uses numerical integration to find the maximum mass flux (I've only tested it with once example though so it may need more work). I've taken a look at the Prode Properties program, namely the ISPF() function. However, if I'm not wrong, it requires the temperature, pressure, and flow rate of the nozzle inlet stream before it can output the area. So then the mass flow rate must be known beforehand. Or can the program calculate that as well? It'll probably come down to me having to use this program or one like it, but I'd just like to know how the program works and how it finds the two-phase mass flow rate W (if it does), or if the mass flow rate for the inlet stream must be known beforehand. If so, how can I calculate it for a two-phase inlet stream?

 

[apetri]

the function ISPF() in Prode Properties given (for a stream) the flow and conditions t,p can calculate the required discharging area (solving the isentropic model with HEM, HNE, NHNE models),
it's a general method good for vapor, liquid and vapor+liquid flows

for a PSV and fire case you need to define the discharging flow which depends from many variables, see API 521 or equivalent standards,
note that properties as latent heat depend from operating conditions and compositions, also wetted surfaces, intensity of radiation etc. are not constant, you need to define some average values or simulate the whole process (not easy since you do not know a priori what is going to happen)
in these cases a good practice is to select conservative values.

 

[Latexman]

The vapor only flow rate combined with the flow area (superficial vapor velocity) is the driving force behind entraining the liquid to give two phase flow. Thinking in terms of fire case, you need the heat input due to fire, the latent heat of vaporization, and the diameter of the tank/vessel. The character of the liquid is also important, like viscosity, foamy, surface tension. Have you looked at the AICHE book "Emergency Relief System Design Using DIERS Technology"? I think Appendix I, especially I-B, would be a big help to you.

Good luck,
Latexman

Need help writing a question or understanding a reply? forum1529

 

[CJC0117]

apetri, I revisited API 521, namely section 5.15.3 for how to find flow rates in the fire case. For the all liquid case a hydraulic expansion equation is used and for the all vapor case the heat input is divided by the latent heat of vaporization (W=Q/L). However, as far as I can tell it is not explicitly stated whether or not the all vapor flow rate can be used for two-phase flow.

Latexman, I have read that section of the DIERS project manual and I have been using that method to predict if two-phase flow is likely. But I don't know how or why the superficial vapor velocity would be used for the two-phase case, if the superficial vapor velocity is just the velocity of the vapor as if it were all vapor flow. Won't the actual velocity differ if liquid is entrained with the vapor? And if the liquid is displacing volume that would otherwise be occupied by vapor, won't the mass flow rate be much higher than if it were all vapor flow? I've been looking through the standards (API 520, 521, 2000) to see if it is explicitly mentioned or heavily implied that the all vapor flow mass flow rate can be used for two-phase sizing, but I've had no luck so far. And it just doesn't make sense to me intuitively that that would be accurate.

I've also gone ahead and attached the excel sheet I have so far for sizing, if anyone's interested in taking a look. It's only good for one component and pseudo-one-component mixtures. Blue boxes are input and red boxes are outputs. It uses the numerical integration HEM method explained in API 520 for finding the maximum mass flux. Equations used include the Antoine equation, Peng-Robinson EOS, and Clapeyron equation. The inputs I currently have in there are from an example I saw in an article entitled "Protection of Storage Tanks from Two-Phase Flow Due to Fire Exposure" by Harold Fisher and Harry Forrest. They assumed the mass flux was equal to the all saturated liquid mass flux, 210.52 lb/s*ft[sup]2[/sup] (which is close to the mass flux calculated from my excel sheet, 220.93 lb/s*ft[sup]2[/sup]. The resulting vent area assuming a perfect nozzle came out to be 4.845 ft[sup]2[/sup] (which means the mass flow rate is 1020 lb/s). However, when I use the all vapor mass flow rate at relieving conditions in my excel sheet, I calculate the area to be 0.09212 ft[sup]2[/sup], which is obviously much less. It's even less than the area I calculate for all vapor venting (0.5037 ft[sup]2[/sup]). Meanwhile, I tried using the ISPF() function in Prode Properties to find the vent area using the all vapor mass flow rate at relieving conditions, and I'm getting 0.43143 ft[sup]2[/sup] (though I won't claim I know how to properly use that program).

 

[apetri]

>I revisited API 521, for the all vapor case the heat input is divided by the latent heat of vaporization (W=Q/L).

with API procedure you can calculate the fluid flow (vapor state) to discharge dividing heat input / latent heat, see also the comments about how to calculate latent heat in recent threads, for that case there is only vapor (no liquid phase)

>I tried using the ISPF() function in Prode Properties to find the vent area using the all vapor mass flow rate at relieving conditions, and I'm getting 0.43143 ft2 (though I won't claim I know how to properly use that program).

I use regularly the Excel page nozzle.xls distributed with Prode Properties which allows to size and rate a PSV,
there are only a few steps required to size (or rate) a PSV,
as a test case we can calculate the area required for discharging a flow of 15000 Kg/h , two-phase mixture with composition C1 0.7 C6 0.3 from 20 Bar.a 300 K to atmospheric pressure.

1) as first step you can define the units of measurement (see image Prodeunits.jpg) for this case I require Bar , K, Kg/h , cm2 , you can select from many US and SI units

2) as second step you should define the composition and the thermodynamic models to apply, (see image ProdeStream.jpg) for stream 1 we define as composition C1 0.7 and C6 0.3 , Soave Redlich Kwong as model

3) finally in Excel page nozzle.xls we set HEM as model, 20 Bar.a and 300 K as Pin, Tin 1.0 Bar.a as tout, 15000 Kg/h as flow then we click on Calculate button to find the solution

the program returns almost immediately the calculated values of 6.82 cm2 for calculated area and 272 K as estimated outlet temperature

note that this is a two phase (gas+liquid) test case but the procedure works for both vapor and vapor+liquid phases

 

[CJC0117]

Thanks for the reply apetri. I hadn't even realized there were all these excel sheets like nozzle.xls saved in the Prode installation folder until I read this post. It makes it much easier to use the program. Just to be sure, the Pout entered is the backpressure, right? I'm not sure because in your example you said you were discharging to atmosphere but then you used Pout=1 bar abs.

>with API procedure you can calculate the fluid flow (vapor state) to discharge dividing heat input / latent heat, see also the comments about how to calculate latent heat in recent threads, for that case there is only vapor (no liquid phase)

So is it acceptable to use flow = heat input / latent heat even for a two-phase inlet stream? Is this how I would calculate the flow (in your example you used 15000 kg/h) even for a two-phase inlet stream?

 

[CJC0117]

Nevermind the question about the atmospheric pressure being 1 bar. I was thinking in psi.

 

[apetri]

the formula flow = heat input / latent heat gives the amount of vapor generated at certain conditions, usually we assume that PSV discharges vapor phase (being on top of the vessel), since liquids to vaporize require a certain amount of heat in this way you can limit the maximum temperature of the liquid inside vessel (and vessel).
A different case could be a PSV on a line after a control valve, depending from process conditions there you could have a vapor+liquid flow but there are many other cases which you may consider.

 

[CJC0117]

Thanks apetri.

I spoke to an expert in the area of two-phase vent sizing and he suggested that when calculating the two-phase mass flow rate, I use the following equations, which can be derived from general principles:

W[sub]v[/sub] = Q/L

dV[sub]v[/sub]/dt = W[sub]v[/sub](1/ρ[sub]v[/sub] – 1/ρ[sub]l[/sub])

W[sub]tp[/sub] = (dV[sub]v[/sub]/dt)ρ[sub]tp[/sub]

Where:
W[sub]v[/sub] = rate of vapor mass generation due to heat input [lb/s]
Q = rate of heat input [BTU/s]
L = latent heat of vaporization [BTU/lb]
dV[sub]v[/sub]/dt = net volumetric vapor generation rate [ft[sup]3[/sup]/s]
ρ[sub]v[/sub] = density of saturated vapor [lb/ft[sup]3[/sup]]
ρ[sub]l[/sub] = density of saturated liquid [lb/ft[sup]3[/sup]]
W[sub]tp[/sub] = two-phase mass flow rate [lb/s]
ρ[sub]tp[/sub] = two-phase density [lb/ft[sup]3[/sup]]

When I implement this method into my excel sheet (attached), calculating W[sub]tp[/sub] at the vessel conditions where the maximum mass flux G[sub]max[/sub] occurs, and then use the equation A = W[sub]tp[/sub]/(K[sub]D[/sub]G[sub]max[/sub]), I get a required area (for a perfect nozzle, K[sub]D[/sub]=1) of 4.186 ft[sup]2[/sup]. This is close to the area calculated in the article I got the example inputs from: 4.845 ft[sup]2[/sup]. Meanwhile, if I calculate the vent area using Prode Properties and nozzle.xls with a stream consisting of all n-heptane, Pin = 0.8597 psig, Tin = 212.7 F, Pout = 0 psig, and Flow = 924.882 lb/s, the required area comes out to be 19.62 ft[sup]2[/sup]. I wonder why there's such a large discrepancy.

 

[apetri]

I am sorry but I fear you are doing some confusion with basic thermodynamics
With Pin 0.8597 psig, Tin = 212.7 F, Pout = 0 psig n-heptane is all vapor (not two-phase), there is not a critical flow,
in this case a formulation for critical flow IS NOT applicable,
in such cases Prode calculates the area solving a constant energy operation but we do not specify the internal area of pipe (the procedure assumes that value) and there could be some errors,
for some devices (tank protection) you may contact the manufacturer and require the flow/pressure curves for the valve.

 

[CJC0117]

First off, when I say G[sub]max[/sub], I mean the mass flux at choked conditions, or, if the flow is not choked (and it's not in this case), I mean the mass flux at the backpressure, 0 psig in this case.

I started my procedure by first assuming homogenous flow such that for a liquid full vessel, the flowing void fraction at peak relieving pressure (P[sub]in[/sub]=0.8597 psig in the example) is equal to the vessel void fraction at peak relieving pressure, α=0. So that means the inlet stream consist of all saturated liquid at P[sub]in[/sub]=P[sup]sat[/sup]=0.8597 psig.

I calculated T[sub]in[/sub]=T[sup]sat[/sup]=212.7 F using the Antoine equation (I recognize that Prode Properties may use different Antoine coefficients that result in a Tsat lower than mine - then the inlet stream as I specified it in Prode would be all vapor, as you say).

I also used the Peng-Robinson EOS to calculate the saturated liquid entropy (relative to some reference state) for the inlet stream. This entropy I held constant throughout all subsequent flash calculations (let's call it S[sub]o[/sub]), as per API's method of modeling flow through the relief valve as isentropic flow through a nozzle.

So then I considered different P[sub]out[/sub]'s in fixed increments between the peak relieving pressure and the backpressure. I calculated the entropy of an all saturated liquid (S[sub]l[/sub]) and an all saturated vapor (S[sub]v[/sub]) at each P[sub]out[/sub]. If S[sub]o[/sub] was less than S[sub]l[/sub] at a specific P[sub]out[/sub], then the fluid would be all subcooled liquid at that specific P[sub]out[/sub]. If S[sub]o[/sub] was greater than S[sub]v[/sub] at a specific P[sub]out[/sub], then the fluid would be all superheated vapor at that specific P[sub]out[/sub]. If S[sub]o[/sub] was in between S[sub]l[/sub] and S[sub]v[/sub] then the fluid would be a two-phase mixture at that specific P[sub]out[/sub].

As it turned out for this example, the constant entropy S[sub]o[/sub] never got higher than S[sub]v[/sub] at any P[sub]out[/sub]'s, so aside from the initial inlet stream consisting of all saturated liquid, the flow through the relief valve remained two-phase all throughout. A maximum mass flux was not reached in between P[sub]in[/sub] and the backpressure, 0 psig, so I used the mass flux calculated for P

 

[CJC0117]

...so I used the mass flux calculated for P[sub]out[/sub]=0 psig as my G[sub]max[/sub].

Pages 120-121 of the DIERS project manual show an example of this procedure.

Assuming I did everything right, then Prode must be using some different method or it's because of the difference in ANotine coefficients I mentioned before. Or not. I'm not sure. I'm just speculating.

 

[apetri]

yes, there could be a little difference in calculating saturation temperature,
I see that for n-Heptane at 0 Psig (101327 Pa.a) Prode predicts a saturation temperature of 371.577 K or 209.168 F which means that at 212.7 F for Prode it's all vapor,

however the main problem in your test case is the very low dp,
Pin 0.8597 Psig
Pout 0 Psig

if you transform to absolute values
Pin 15.556 Psia
Pout 14.696 Psia

the ratio Pin/Pout in close to 1.0 which means that you are far from the correct range of application for the isentropic nozzle model in Prode Properties,
as said in my previous post you should not use that model with these conditions,
if you increase pin to 2 or 3 psig (pout being 0 psig) you get more reasonable areas even if you are still very far from critical conditions,
a better approach would be to use the Constant Energy Flash operation where you can specify the internal area of pipe.

However for these very low pressure applications (see API 2000) I believe that the best option is to ask the manufacturer to provide certified flow/pressure curves for a specific model.

 

[CJC0117]

>I see that for n-Heptane at 0 Psig (101327 Pa.a) Prode predicts a saturation temperature of 371.577 K or 209.168 F

Is there a Prode function that can be used in excel that outputs the saturation temperature for a given pressure?

>the ratio Pin/Pout in close to 1.0 which means that you are far from the correct range of application for the isentropic nozzle model in Prode Properties, as said in my previous post you should not use that model with these conditions

Do you mean the isentropic nozzle model itself is invalid at these conditions or that Prode just yields invalid results? I'm not sure I've seen mention of the validity of the isentropic nozzle model relying on pressure drops and ratios. Could you elaborate on this or point me to a reference?

>a better approach would be to use the Constant Energy Flash operation where you can specify the internal area of pipe.

Do you have a reference where I could see an example of how this is done and what the algorithm behind it is?

>However for these very low pressure applications (see API 2000) I believe that the best option is to ask the manufacturer to provide certified flow/pressure curves for a specific model.

This may be a foolish question, but wouldn't such a curve depend on the specific relief fluid, and the specific emergency scenario?

 

[apetri]

pardon me but the title on the Excel page nozzle.xls (which you are testing) is
"Sizing of relief valves (critical and two-phase flow) with isentropic nozzle model"
which I believe should clearly explain the scope and the field of application.
Your test case with vapor, Pin 15.556 Psia Pout 14.696 Psia does require a different procedure,
if it is a low pressure protection device you can contact the manufacturer and ask for the certified flow/dp curves.

 

[PaoloPemi]

::Is there a Prode function that can be used in excel that outputs the saturation temperature for a given pressure?

Properties includes different functions which returns temperature or pressure given a phase fraction (and pressure or temperature) you can specify a vapor, liquid or solid phase fraction,
in the operating manual there are examples provided about how to use these functions from Microsoft Excel, see page 22

:o you have a reference where I could see an example of how this is done and what the algorithm behind it is?

Prode includes functions for solving the constant energy flash operation, see the method EPF() page 59 of operating manual,
these methods permit to model adiabatic, irrversible expansions when the contribute of kinetic energy cannot be neglected,
I use these methods to model the restriction orifice when depressuring a vessel but there are many other applications where they can be useful.

Finally, I agree with apetri about the range of application of the isentropic nozzle model,
for your test case it would be better to use the "performance curves" usually provided by manufacturer

 

[CJC0117]

Thanks PaoloPemi. I didn't know what was meant by "constant energy flash operation" but I get it now after reading the description for the EPF() function.