firemodels / fds

Fire Dynamics Simulator
https://pages.nist.gov/fds-smv/
Other
663 stars 622 forks source link

Final energy in closed adiabatic box not same for N2 injection mass which varies by rate only #1023

Closed gforney closed 9 years ago

gforney commented 9 years ago
Based on discussion thread 
http://groups.google.com/group/fds-smv/browse_thread/thread/10c8146fa45a3c28#

Using FDS 5.4.3 SVN 5210, Win7 32bit

I created two simple models of a closed room with adiabatic walls.  An
inlet on the floor injects 2000degC N2 into the room.  One file
injects over 5 seconds (fast inject). The other over 20 seconds (slow
inject). Both inject the same total mass of N2 (verified using mass
dump).

Radiation is False, Stratification is False, Radiative_Fraction is 0,
Background_species=Nitrogen.

I integrated the CONV_LOSS in the hrr file and both integrate to the
same total energy added. RAD was 0 and COND was negligible.

However, integration of the energy in the room using
&DEVC ID='intEnthalpy', XB=0,30,0,30,0,40, QUANTITY='SPECIFIC
ENTHALPY', STATISTICS='MASS INTEGRAL' /
shows a 11% difference in energy at 30 seconds (TEND). Pressure is
also about 11% different.

I would have expected that the final state of the box would be
independent of the path (injection rate) since no energy is lost and
the same total mass and energy is added. I would appreciate any
comments on this regarding the physics I am not seeing or the code
behavior.

I have upload to the files section the zip file "fds_test.zip" which
has the input files and some png files showing mass, conv_loss, and
integrated specific energy plots. The input files are less than 20
lines and run in 1 minute. Thanks. 

You can also find the files here
http://drop.io/8bp12sv
or see attached.

Original issue reported on code.google.com by grgrgrknght on 2010-03-31 14:49:35


gforney commented 9 years ago
I'll take a look at it.

Original issue reported on code.google.com by mcgratta on 2010-03-31 15:06:17

gforney commented 9 years ago
The box is not completely adiabatic. The injection vent will create or drain energy

from the system, depending on its temperature. In your case, the vent is a different

size. This may not explain everything, but have you checked how much the heat losses

from the VENT contribute?

Original issue reported on code.google.com by mcgratta on 2010-03-31 16:58:04

gforney commented 9 years ago
I was aware of that and that is why I set the vent surface temp to 20C right after
the leak ends. If you don't do that, conv_loss shows a clear slow rise/fall due to
conv. after the leak ends. Even if vent size allows some extra conv energy during the
leak, it would be captured in the conv_loss integration that shows both models energy
input as equal.

Original issue reported on code.google.com by grgrgrknght on 2010-03-31 17:02:55

gforney commented 9 years ago
I created Fillbox_r0a.fds which is the slow inject with a large vent similar to fast
inject. I decreased the MASS_FLUX to account for larger vent. Total mass injected is
still the same between r0a and r1.

Convection energy increased as did conduction as expected, but more that I though it
would. See attached. However, integrated specific enthalpy DECREASED for r0a? See plot.

Original issue reported on code.google.com by grgrgrknght on 2010-03-31 17:32:46


gforney commented 9 years ago
I wonder if this is an issue of how we compute the fixed T BC.  We compute CP_TERM
using TMP_F and RHOWAL (the average of RHO_W and RHO_G) and the quantities
U_W,RHO_W,and TMP_W are not computed simultaneously even though each is dependent
upon the value of the other.  Since we are no longer using constant cp, do we still
have the correct implementation of this boundary condition (i.e. given the user
specified T and U at the inlet, do we input the correct enthalpy?)

Original issue reported on code.google.com by drjfloyd on 2010-03-31 18:11:44

gforney commented 9 years ago
But the two cases have nearly the same CONV_HRR from the _hrr.csv file. That is just

the flow of energy in or out of the domain. I'm not sure we can rule out the effect

of the vent itself not being adiabatic. 

Original issue reported on code.google.com by mcgratta on 2010-03-31 18:21:50

gforney commented 9 years ago
If I set H_FIXED = 0 on the vent I see similar behavior (with a slightly lower
enthalpy change).

In dump we compute H_G using the average of wall and gas temperatures.  In wall,
compute using only the wall temperature.  In wall we update the wall BC's in serial
and not in parallel.  Pressure is computed based on USUM (volume flow) and DSUM.  If
there is an inconsistency between dump and wall in how the boundary is being
evaluated that could be reflected in the USUM.

Original issue reported on code.google.com by drjfloyd on 2010-03-31 18:49:36

gforney commented 9 years ago
Good point -- I'm doing an experiment now using just lagrangian particles that do 
nothing but blow N2. In my first trial, 

p2/p1 = (rho2/rho1)^gamma

to about 0.1% error, assuming gamma=1.4 (which we don't use directly in FDS, but it

should nonetheless be close). I believe that this scenario is isentropic.

Original issue reported on code.google.com by mcgratta on 2010-03-31 18:56:35

gforney commented 9 years ago
Attached are three plots I just made where particles within a 1 m cube inject N2 at

a fast and slow rate. The horizontal lines indicate the ideal result, like what I 
wrote in the last comment. This is telling me that the problem is most likely in the

way the gas is injected or heat losses to these non-adiabatic surfaces. One other 
thing to check is to reduce the time step to make sure that this is not just a 
problem with integration accuracy.

Original issue reported on code.google.com by mcgratta on 2010-03-31 20:53:03


gforney commented 9 years ago
The other two files.

Original issue reported on code.google.com by mcgratta on 2010-03-31 20:55:14


gforney commented 9 years ago
I modified the cases to get rid of the linear ramps (ramped up and down over 0.0001
s
rather than 1 s) and set H_FIXED to zero.

The input specifies the injection of 640 kg of N2 at a temperature of 2000 C.  The
total enthalpy addition should be 640 kg * 2274 K * 1.173 kJ/kg-K (cp_bar at 2000 C)
or 1.71E6 kJ.

For 5 s of injection time the HRR file shows a change of 1.68E6 and the DEVC file
when I do H+dP V gives 1.66E6 kJ.   

For 20 s of injection time the HRR file shows a change of 1.65E6 and the DEVC 1.2E6.

For 1 s the values are HRR 1.8E6 and DEVC 2.5E6

Original issue reported on code.google.com by drjfloyd on 2010-04-01 13:57:57

gforney commented 9 years ago
fixed the time step at a small value and the results were basically the same.

Original issue reported on code.google.com by drjfloyd on 2010-04-01 15:51:51

gforney commented 9 years ago
I re-ran my test problem with 5.5 SVN 6004. All results have H_FIXED=0 on SURF line.

1. CONV_LOSS is much cleaner and there is no difference in final integrated value
between slow inject over 64m2 (r0), slow inject over 256m2 (r0a), or fast inject
(r1). The absolute integrated value was close to the 5.4 value. See attached image
hrr.png.
2. Integrated enthalpy in the domain is much closer. About 7-8% difference between
r0, r0a, and r1 with r1 still higher. Note the integrated enthalpy value is more than
double the 5.4 result, and hence pressure is higher also. Is this due to a change in
DEVC logic for the enthalpy integration or something else? 
3. The delta domain energy is significantly higher than the final CONV_LOSS
integrated value. See reference line on devc.png plot.

If there is anything I can do to help in terms of running model configurations let
me
know. Also, is there a detailed FDS change log? I saw the summary release notes.

Thanks for your help guys.

Original issue reported on code.google.com by grgrgrknght on 2010-04-08 15:40:24


gforney commented 9 years ago
On 2) some of this I think may be error in time.  If the time steps are locked to
smaller values the error gets better.

On 3) Two things are happening.  a: Enthalpy is increasing due to the enthalpy
associated with the added mass.  b: Enthalpy is increasing as a result of pressure
work (it takes energy to raise the pressure of a sealed box).  If you subtract (delta
P x Volume = Joules of pressure work) from the enthalpy in the DEVC file, you should
get a value close to the gas enthalpy from the HRR file

Original issue reported on code.google.com by drjfloyd on 2010-04-08 15:46:47

gforney commented 9 years ago
Take a look at the new FDS Verification Guide. There is a case in there that I added

called isentropic2, under Energy Balance. It describes a similar case. 

Original issue reported on code.google.com by mcgratta on 2010-04-08 16:47:06

gforney commented 9 years ago
I'm going to close this case. If there are further developments, just continue to 
add comments and I will re-open it. Thanks.

Original issue reported on code.google.com by mcgratta on 2010-04-20 13:27:26