Time-rate of change of P and T in a pressurized rigid volume with a hole in it

I am trying to model the time-rate of change of pressure and temperature of a gas being held within a pressurized rigid volume as it escapes: an unsteady, open system with no heat transfer. Most examples of an open system provided in-class and online typically presume steady-state conditions inside the control volume of interest. I have yet to find any example that would adequately model the change of state for a gas leaving a rigid container (say, for example, a cold gas thruster where the propellant is entirely in a gaseous state, or an air compressor tank being depressurized).

Take this situation for example: an insulated pressurized rigid volume is connected to a nozzle downstream and is allowed to evacuate. I would like to know what the time-rate of change of pressure and temperature of the gas both within the volume AND seen by the nozzle inlet (assuming they're different).

Intuition tells me that the temperature of the gas inside the volume will decrease since it is doing work to push out. I know the nozzle process is isentropic (that part I have done -- it's everything upstream of it that I wish to know more about), but I don't know what kind of process the evacuation is. I've thought that it might be Joule-Thompson expansion, but I think that would only apply if I take the system as a whole (including the nozzle). What I can't figure out is how to get the state of the gas remaining in the volume.

My first attempt at modeling this went like so:

I began with the general energy balance equation for a control volume:

$$\frac{dE_{cv}}{dt}=\dot{Q_{cv}}-\dot{W_{cv}}+\dot{m}e_{in}-\dot{m}e_{out}$$

There is no inward mass or heat flow (but it's not adiabatic since mass is leaving the system), so this reduces to the following:

$$\frac{dE_{cv}}{dt}=-\dot{W_{cv}}-\dot{m}e_{out}$$

The energy associated with the outward flowing gas ($$e_{out}$$) can be rewritten in terms of enthalpy ($$h$$) and kinetic energy ($$\frac{c_{out}^2}{2}$$):

$$\frac{dE_{cv}}{dt}=-\dot{W_{cv}}-\dot{m}(h_{out} + \frac{c_{out}^2}{2})$$

Then by dividing through by the mass flow rate, we can put everything in terms of specific energy:

$$\frac{de_{cv}}{dt}=- \dot{w}_{cv} - h_{out} - \frac{c_{out}^2}{2}$$

So what does this tell me? Exactly what my intuition said: that the total (specific) energy of the gas within the control volume ($$cv$$) will change over time as a function of the work done by the gas in the cv, the enthalpy of the gas which has left the cv, and the velocity of said gas.

So here's where I get stuck:

1. Work is defined as changes in volume at a constant pressure on a closed system, but this is practically the exact opposite! It's a change in pressure at constant volume for an open system. Formally speaking, that would mean that $$\delta{W}=PdV=0$$, but I think that only applies to closed systems. Since this is an open system, does that mean I need to take into account the full work derivative? Does $$dW=PdV+VdP$$? And since the specific volume in the cv is changing, would that imply that $$dW=Pd \nu + \nu dP$$? If that's the case, then I would jump to the following relationship:

$$\dot{w}_{cv} = P_{cv} \frac{d\nu_{cv}}{dt} + \nu_{cv} \frac{dP_{cv}}{dt}$$

Is this a correct assessment for this situation?

1. What is the rate of specific energy change inside the cv ($$\frac{de_{cv}}{dt}$$)? Is that just simply the enthalpy of the gas inside of it (since enthalpy includes internal energy and $$PV$$ energy, this tells me that the energy stored in the compressed gas gets taken into account in the enthalpy term).

Here's what I do know:

The nozzle downstream is controlling the mass flow rate. It's going to be an implicit function since mass flow rate through a nozzle is partially a function of the gas temperature -- but a change in mass flow rate will illicit a change in the gas temperature as it is leaving the volume. I can solve that iteratively though, I think (I already have an extensive program to model the nozzle flow, so I should be able to tack this onto it).

My initial assumption for the gas leaving the volume was that it was isentropic, but that is definitely wrong seeing as it predicts a nearly -80$$^\circ$$C drop in temperature (and this supposition is supported by one very quick blurb in this video). To put it succinctly: I don't know what kind of process models the gas leaving the volume. It's not isentropic, and it doesn't seem to be isenthalpic, it's certainly not isothermal. It IS isochoric but I don't know what equations are associated with that.

I would love to get some insight, especially if any of my initial assumptions are invalid. I'm very rusty at this. I took thermodynamics over a decade ago :| Thank you!

After brushing off the old textbook, I found where my knowledge gap was. Every example given online, in books, etc. always seems to assume steady flow, even in open systems such as this. That's easy enough because when you do the mass-energy balance on the system, the $$\frac{dE}{dt}$$ (which really should be $$\frac{dU}{dt}$$ because it's a change in internal energy) side of the equation will go to zero. But what if it isn't steady-state?

$$\frac{dU}{dt}$$ merely represents the rate of change of internal energy in the c.v., but if you want to know the pressure and temperature at any given time, you have to know what the absolute value of that energy is. Therefore, I think the solution is to track the total internal energy ($$U_{cv}$$ -- that is, $$um_{cv}$$).

From there, at each time step of your simulation, you would use the specific volume $$v$$ and specific internal energy $$u$$ to determine the pressure and temperature inside the c.v. Since the mass flow rate is controlled by the nozzle downstream, you therefore know $$\dot{m}$$ and can use that knowledge to determine the change in mass inside the control volume.

So for my formulation, I did the following:

$$\delta{U} = \delta{Q} + \delta{W} + \delta{M}$$

For an open, unsteady system, our rate of change of total energy ($$E$$) is the sum of all the heating ($$Q$$), work ($$W$$), and mass flows ($$M$$).

$$\dot{U}_{cv} = \dot{Q}_{cv} + \dot{W}_{net} + \dot{M}_{cv}$$

So now let's break it down a little....

1. $$\dot{M}_{cv} = \dot{m}(u + \frac{c^2}{2})_{in} + \dot{m}(u + \frac{c^2}{2})_{out}$$

This is the energy of the fluid that is moving in/out of the system. It has some internal energy ($$u$$) and some velocity ($$c$$). I am neglecting gravitational potential energy ($$gz$$) and chemical potential energy (??). The internal energy ($$u$$) is the first half of the enthalpy term which I'll introduce in a sec. There is no flow in, so only the out term remains.

1. $$\dot{W}_{net} = \dot{W}_{shaft} + \dot{W}_{flow}$$

This is the work done by or on the system, including the work done on the stuff that is inside the c.v. (Note that you'll find slightly different vocabulary to describe all the different kinds of work that can be done on a system. For this, I'm lumping any work done by pistons, propellers, etc. into "shaft", and the work done by moving a fluid that has some momentum associated with it as "flow"). Now let's break out the work terms a little:

2a. $$\dot{W}_{shaft} = 0$$

No work is being done on the fluid in the form of any external device used to stir or otherwise change the energy of the fluid once it is inside the control volume.

2b. $$\dot{W}_{flow} = \dot{m}(Pv)_{in} + \dot{m}(Pv)_{out}$$

This is where "flow work" is taken into account, but i find the name to be a little misleading because "flow" (in my mind, at least) has more to do with velocity than anything else. At any rate, this is all the energy associated with compressing the fluid and moving it into/out of the system. It is the second part of the enthalpy term which I'll introduce in a sec. There's no flow in, so the only term is the outward term.

1. $$\dot{Q}_{cv} = \dot{Q}_{in} + \dot{Q}_{out} = 0$$

This is all the heat that is being added/rejected from the system. It's adiabatic, so no heat is added or rejected and thus it equals 0.

Now put it all together!

$$\dot{U}_{cv} = \dot{m}(Pv)_{out} + \dot{m}(u + \frac{c^2}{2})_{out}$$

Recall that enthalpy ($$h$$) is a term used to combine the internal energy ($$u$$) and flow work ($$Pv$$) strictly as a matter of convenience. It has relations to measurable properties, but is otherwise difficult to intuit all on its own. However, since $$Pv$$ and $$u$$ appear so often together, scientists and engineers have started tabulating enthalpy in addition to internal energy and specific volume. So what we end up with is:

$$\dot{U}_{cv} = \dot{m}(u + Pv + \frac{c^2}{2})_{out} = \dot{m}(h + \frac{c^2}{2})_{out}$$

Now that I have my time-rate of change of energy in the system ($$\frac{dU}{dt} = \dot{m}(h + \frac{c^2}{2})_{out}$$), I can use my downstream isentropic flow properties to determine $$\dot{m}$$. Furthermore, assuming I know something about the geometry of the flow channel, I can use the same isentropic relations to determine the Mach number -- and therefore velocity -- by assuming that the stagnation pressure ($$P_0$$) is constant (this is where I become uncertain about my assumptions...I feel like I'm mixing metaphors, in the sense that I'm mixing isen- assumptions).

Once I have my $$\frac{dU}{dt}$$, I loop it through my code by doing the following:

1. Calculate $$\dot{m}$$ using my nozzle code (not covered here) based on initial stagnation conditions $$P_t$$ and $$T_t$$

2. Use $$\dot{m} \Delta{t}$$ to calculate the remaining gas mass in the c.v. and the specific volume $$v$$

3. Calculate the change in energy $$\frac{dU}{dt}$$ using the aforementioned mass-energy balance equation and properties of CO2 to find $$h$$, and isentropic Mach-Area & Mach-Density relations to find $$c$$ (turns out it's pretty negligible though).

4. Calculate the total internal energy of the remaining gas in the c.v. ($$u m_{cv}$$). For this, I needed specific fluid properties which I acquired from NIST and wrote some functions to perform linear interpolation between points (i.e. input $$P$$, $$T$$, and get $$u$$ returned).

5. Subtract $$\frac{dU}{dt} \Delta{t}$$ from that total internal energy and recalculate specific internal energy at the new mass.

6. Use the same NIST data to interpolate and find a new $$P$$, $$T$$ based on $$u$$, $$v$$.

7. Repeat

What I got as a result was that my small, 30 $$cm^3$$ pressure vessel that was pressurized to 100 psig drops by 60 K down to 145 K in temperature. What is interesting is that the isentropically-calculated rate of pressure drop was very similar to this rate, but the isentropically-calculated temperature was WAY LOWER than this but this is EVEN LOWER than the isentropic temperature. What's more, this seems much more realistic based on my experimental results (wherein I did not observe any phase changes or extreme temperature swings from a single expulsion of gas) This cannot be realistic, can it?

• Your approach looks correct. I ran a quick calculation venting a 30cc chamber full of air through a 1mm hole and got a similar temperature. I would also do a time-step sensitivity study to make sure $\Delta t$ is small enough. Also, with such a large temperature gradient there is bound to be some noticeable heat transfer into the vessel, so the adiabatic assumption might start to break down if the discharge time is long enough. – Carlton Aug 26 '20 at 20:26
• Thank you for double-checking, but unfortunately I realized that my original simulation had a few small errors. When I fixed them, I found that the temperature dropped even greater than the isentropic temperature :( (I've already done a $\Delta t$ study to make sure the results converge at least). Using isentropic assumption, $T$ drops to ~127 K. Using the conservation of energy, $T$ drops to ~143 K. What values did you end up seeing? And what was your approach? – J. Day Aug 27 '20 at 6:58
• I used a 30 cc vessel venting through a 1 mm diameter hole, using a choked flow assumption. You have a nozzle, so the mass flow rate will be a bit different than mine. Then, I solved the ODEs for $\dot m$ and $\dot U$ using the forward Euler method. Knowing $m$ and $U$, I then used the ideal gas equation of state and the equation $U=m \cdot c_v \cdot T$ to calculate the pressure. I solved this system of equations at discrete time steps using $\Delta t = 10^{-2}s$. I assumed the fluid was air, btw. – Carlton Aug 31 '20 at 15:31
• After $1 s$, the pressure in the tank is almost zero, and the temperature has dropped to $47 K$. Although one second is a short period of time, $47 K$ is a very cold temperature, and there would probably be some noticeable heat transfer from the surroundings to the tank (and the air in it). The next thing I would add to this model is a simple heat transfer calculation to see how significant it is. – Carlton Aug 31 '20 at 15:37

Thanks for a nicely presented first question.

First you need a model for the equation of state for the tank. Given the initial conditions, develop formulas that yield static pressure, temperature, work, etc as a function of the mass remaining in the tank. (You need at least one fully defined state point for the tank, I'm assuming it's the initial condition).

Differentiate these wrt time to get rates of change as a function of mass flow rate.

Use your nozzle equation to relate $$P_{in} (static) = P_{tank}$$ to mass flow rate.

You probably need an equation to relate $$\partial T/\partial m$$ to $$\partial P/\partial m$$ as well. That should be all of the pieces.

I don't know why you think an 80 C drop is unreasonable in general. This is how oxygen and CO2 were first made into liquids.

mass $$m = f(P,T)\quad$$ initial mass $$m_1= f(p_1,T_1)$$ for some gas at constant volume.

You might try the Beattie - Bridgeman real gas model.

\begin{align*} p&={\frac {RT}{v^{2}}}\left(1-{\frac {c}{vT^{3}}}\right)(v+B)-{\frac {A}{v^{2}}}\\\text{where}\quad\quad\\A&=A_{0}\left(1-{\frac {a}{v}}\right)&B&=B_{0}\left(1-{\frac {b}{v}}\right)\end{align*}

The five constants for common gasses are widely available. $$v$$ is often in the odd units of $$\frac{m^3}{k\,mol}$$, check this when looking up the constants.

• For this application, I'm using CO2 at 20C and 100 psig as the initial conditions, so I've just been assuming an ideal gas. At that pressure differential, 80C seems excessively high (especially since I've done this experimentally and found nowhere close to that big of a drop). I re-derived the energy balance and found that $h_{cv} = h_{out} + \frac{v^2_{out}}{2}$, so that describes the process (right?), and the ideal gas model describes the state. I'm just getting lost on how to connect everything together in one cohesive model. – J. Day Aug 23 '20 at 20:08

The $$\dot W$$ term in the energy balance usually represents forms of work other than flow work. Since you have a rigid tank ($$dV=0$$), there is no other work being done on/by the system, so you can eliminate $$\dot W$$ from the energy balance. It then simplifies to $$\dot E = \dot m (h + \frac{v^2}{2})$$. I use $$v$$ for velocity btw. Since the internal energy is changing with respect to both mass and temperature, I'd write $$\dot E = \dot m c_v T + mc_v\dot T$$. You can plug that into the energy balance to get $$\dot m c_v T + mc_v\dot T = \dot m (h + \frac{v^2}{2})$$. Note that the rate of change in temperature depends not just on the mass flow rate, but also the instantaneous mass and temperature in the system.

It sounds like you already have an equation for $$\dot m$$ from your nozzle calculations, so it should be straightforward to solve for $$T(t)$$, $$m(t)$$, and $$v(t)$$. You just have to integrate two ODEs ($$\dot T$$ and $$\dot m$$) which you might have to do numerically depending on your nozzle equation. You'll also need initial conditions for $$m$$ and $$T$$.

Edit: On second thought, it might be easier to keep $$\dot E$$ in the energy balance and integrate it to find $$E(t)$$. Then you can figure out $$T(t)$$ directly using $$E=mc_vT$$ (ignoring gravitational and kinetic energy of the fluid).

• One of the things that's throwing me off with this is the work term. Indeed since it's a rigid tank then $dV=0$, but the specific volume IS changing. So even though $\Delta{W}=P\Delta{V}=0$, it can also be said that $dw= P d \nu \neq 0$. Does this not accurately reflect the work being done by the gas inside the control volume? Or am I mixing up something here? – J. Day Aug 23 '20 at 18:31
• To do work requires that something be moved by a force. That can be pressure acting over a change in volume, or force acting over a change in distance, etc. Although the specific volume is changing in your tank, none of its boundaries are moving, which means there's no $PdV$ or $Fdx$ work. There is work being done where the gas escapes into the nozzle however. Picture a small control volume of gas that expands against the atmosphere as it escapes. It's expanding against a pressure, and it's volume will increase as it equilibrates with the atmosphere, so it meets all the criteria to do work. – Carlton Aug 24 '20 at 1:54
• So the $\Delta W = P \Delta V$ term specifically only applies to the control volume? Okay. I think the concept I didn't understand is that work can come in multiple forms, even for a fixed volume. The key concept I was missing was the idea of "flow work" which, for an open system, is expressed as $\dot{W}_{flow} = \dot{m} e_{flow}$ where $e_{flow} = u + \frac{c^2}{2} + gz$ (u is internal specific energy, c is velocity, and gz is gravitational). – J. Day Aug 25 '20 at 1:59