firemodels / fds

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

H_FIXED and H_LOGLAW #999

Closed gforney closed 9 years ago

gforney commented 9 years ago
Application Version: FDS 5.4.3
SVN Revision Number: 5629
Compile Date: 22 feb 2010
Operating System: Windows 7

This is a follow up of issue 987 which was closed believing all was fine.

I ran the validation suite SP2009-test with all the different convection
models availiable. Two things looks strange.

1. H_FIXED for the PT front gives h = 25 for PT ABC-1 and -2 but not for PT
ABC-3 and -4.

2. H_LOGLAW gives a heat transfer coefficient with a magnitude of 0.5 which
seems very low. All the other models (except H_EDDY) gives h with a
magnitude of at least 5.

Joakim

Original issue reported on code.google.com by jsandstrom42 on 2010-03-02 16:14:36


gforney commented 9 years ago
We can re-open closed cases. But we'll  just work this new issue here.

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

gforney commented 9 years ago
Joakim,

I appreciate the want to "give it a try" with new features.  But we are still very
much in the development phase for both of these options (H_EDDY and H_LOGLAW) and so
I would like to stick to simpler problems.

In particular, I am interested in recreating Nusselt number correlations using FDS
the same way we recreated the Moody diagram for the friction factor (attached).  I
very much liked your vertical wall test case.  But even that can be simplified.  So,
to get us started I have created a case called 'vertical_convection.fds' (attached).
This case mimics section 7-11 of J.P. Holman on Free Convection in Enclosed Spaces.

The left and right walls are held at a fixed temperature.  The enclosure spacing and
the temperature difference together with the default fluid (air) determine the
Rayleigh number (Ra).  For the L/d of 10 in this case, the Nusselt correlation is

Nu = .046*Ra^{1/3}

The DEVCs measure the heat flux at both walls.  When they are at steady state, then
we determine h = q/(T_hot-T_cold), and from there Nu = h*d/k, where k is the thermal
conductivity.  This becomes one point on our Nu vs. Ra map.

It may take quite a long run to reach steady state.  I have attached results from a
quick test run to 100 sec.  As can be seen from the plot, the heat flux is not close
to steady state.

So, the first step is to verify that we reach a steady state with this case.  Else
there is a bug.  Next, we can pick a typical Ra for fire cases and test the models
for that case.  The log law model admittedly still requires tuning.  Then we need to
map out the remaining parameter space and study grid dependence.  If you would like
to help with these details, that would be great.  If not, I understand you probably
have other priorities.  But I just want to make sure we stay simple for now.

Thanks for your help.
Randy

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 18:37:32


gforney commented 9 years ago
(No text was entered with this change)

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 18:38:13

gforney commented 9 years ago
There indeed seems to be a bug in H_LOGLAW.  It does not attain a steady state.  But
the default FDS correlations do.  I have attached updated input files and results.

FDS gives h = 3.8 W/m*C and the correlation gives h = 2.8 W/m*C.  So they are in the
ballpark.  Obviously need to fine tune.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 19:15:21


gforney commented 9 years ago
Status: H_EDDY flies off the rails, like there is some amplifying feedback loop.  Log
law model becomes steady but q does not match on the hot wall and cold wall...
clearly a bug somewhere.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 19:53:48

gforney commented 9 years ago
Here's something interesting:  With H_FIXED (hot) = 10 and H_FIXED (cold) = 5, the
system reaches a steady state and both values of q match (so far so good).  But the
global value of h is 3.1 based on the output from the DEVCs.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 20:04:40

gforney commented 9 years ago
I forgot my thermal resistance formula!

q = (T1-T2)/sum(R) = dT/(1/h1 + 1/h2).

When this applied, the global h based on the H_FIXED values comes out correctly.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 20:53:00

gforney commented 9 years ago
Attached is the latest input file.  I played with this case for quite a while today
and came to the conclusion that, as currently implemented, H_EDDY and H_LOGLAW do not
work at all -- H_EDDY is at least an order of magnitude too high and H_LOGLAW does
not correctly treat the cold side wall and never achieves an energy balance.  At one
point I worried that the periodicity was having an effect.  So I just made the
y-direction walls adiabatic.

So, I like this case as a place to start. The FDS default models seem to work well
with either FDS5 or FDS6.  Until we see something do a better job for this case, we
will stick to these correlations.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-02 23:27:21


gforney commented 9 years ago
I tried to use the simple convection model that we created to force the error of
H_FIXED not working properly but could not make it work. All the surfaces with
IOR=negative got a variable H_FIXED close to that of the default model.

Added a devc for h just to see. For H_EDDY, the heat flux on the cold wall is
identical over the height but for the hot wall it differs over height. h over the
cold wall differs over hight which is strange since the heat flux don't. See the
attached excel sheet. (data for h is hidden in vertical_convection)

I do have other things to do as well but this is important for my studies as well so
I can just as well spend some time understanding the models and help out the best I
can.

The loglaw model seems to be a surface model developed for weather forecast. Are
there any previous studies on this model for the fire application? Where should I
start if I want to learn more about heat transfer models for LES (or other models
applicable for FDS)?

Do you have an mpi version of FDS, I can only compile the serial version and it takes
forever...

Joakim

Original issue reported on code.google.com by jsandstrom42 on 2010-03-03 15:00:19


gforney commented 9 years ago
Since for H_EDDY, the relation q=h*deltaT doesn't add up. Are they calculated
separately in this method or is it a numerical error of some other kind? Shouldn't
q
be calculated based on h or the other way around and shouldn't that mean that they
would always add upp in this comparison of q=h*deltaT even if h is calculated wrongly?

Joakim

Original issue reported on code.google.com by jsandstrom42 on 2010-03-03 15:26:57

gforney commented 9 years ago
Joakim,

Thanks for looking at this.

I am unclear on some of your comments.

I think the H_FIXED works fine.  Let's back up. Your new case with all the new DEVC
is more complicated than we need at the moment.  Before we get into H_EDDY, convince
yourself that H_FIXED works.  As I said above, you have to think of this in terms of
thermal resistance.  At steady state, the heat flux is

q = (T_hotwall - T_coldwall)/(1/h_hot + 1/h_cold)

So, from the original set of DEVCs let the case get to SS and then check that this
relationship holds.  I find that it does.  If you find otherwise, we need to
reconcile this simplest of cases before we go on.

R

Original issue reported on code.google.com by randy.mcdermott on 2010-03-03 18:01:51

gforney commented 9 years ago
Bad input file, forget my comments... Discovered some errors today. It seems as
H_FIXED works fine, the relationship holds.

One question on the default model. The forced convection is calculated with L=1. It
gives h gets smaller than for L<1. As you might have noticed, I look into the SP2009
case in the validation suite where L on the Plate Thermometers never gets larger than
L=0.1. Is it so that the forced convection coefficient for these plates should be 10
times that calculated as forced convection in FDS where L=1? As of now, the natural
convection is almost always dominant (ending at h=10 in version 5.4.3 and the forced
convection giving h=5).

Assuming that the above is true, should it be that h for the PT's should be 50
instead of 10?

Just curious

Joakim

Original issue reported on code.google.com by jsandstrom42 on 2010-03-04 09:57:58

gforney commented 9 years ago
Joakim,

For the case you described (L=.1), h will increase by roughly 1.6 instead of 10. 
This is because the L dependence on h ~ L^{-1/5}.  Here's why.  The Nusselt number
correlation is (I think from Incropera and DeWitt)

Nu = h*L/k = 0.037*Re^{4/5}*Pr^{1/3}

Now set k = Cp*mu/Pr and Re = rho*v*L/mu.  When you work this out, you will see that
h ~ L^{-1/5}.  This weak dependence on L is what allows us to get away with
neglecting it.  But the justification for these assumptions must rest on a successful
validation case.  I think we need to proceed as we are, testing FDS's ability to
reproduce relevant Nu correlations for a range of Re and Ra.

R

Original issue reported on code.google.com by randy.mcdermott on 2010-03-04 14:42:49

gforney commented 9 years ago
Forgot about the Re, even though I should have remembered.

Thanks for all the help

Joakim

Original issue reported on code.google.com by jsandstrom42 on 2010-03-04 16:03:26

gforney commented 9 years ago
I have a follow up on question 1. I have tried to assign H_FIXED to the plate
thermometers in the SP2009 validation suite calculateion. This works well for half
of
the PT's but for the other half, h is the default value (I assume).

I have tried to see if it could be that the devc are placed on the "wrong side" of
the PT but this does not seem to be the case.

The first of the attached files produced my error. I tried to simplify the case but
it didn't produce the error.

Original issue reported on code.google.com by jsandstrom42 on 2010-03-08 09:17:19

gforney commented 9 years ago
Joakim,

Please attach your input file.  I need to see exactly how you are trying to output
h.

Thanks.
R

Original issue reported on code.google.com by randy.mcdermott on 2010-03-08 13:21:58

gforney commented 9 years ago
I use the line

&DEVC XYZ = ..., IOR = ..., QUANTITY = 'HEAT TRANSFER COEFFICIENT', ID = 'H_whatever'
/

Original issue reported on code.google.com by jsandstrom42 on 2010-03-08 13:28:12


gforney commented 9 years ago
Joakim,

You need to set H_FIXED=5. on the 'PLATE BACK' as well.  Note that in the attached
input file I have removed the DUMP line and changed T_END so you can run a quick test
to confirm.

R

Original issue reported on code.google.com by randy.mcdermott on 2010-03-08 13:53:54


gforney commented 9 years ago
Thanks, is there an answer to why? so I don't do the same mistake again..

Original issue reported on code.google.com by jsandstrom42 on 2010-03-08 14:01:53

gforney commented 9 years ago
Remember that H_FIXED is tied to the particular SURF line.  Notice that the OBST line
used to create the plate thermometers use both 'PLATE BACK' and 'PLATE FRONT'.  They
have zero thickness and so the IOR on the PT DEVC tells the device which side to choose.

Does this answer your question?

Original issue reported on code.google.com by randy.mcdermott on 2010-03-08 14:50:28

gforney commented 9 years ago
Yes it does. I thought that correct IOR would be sufficient. It works now anyhow,
thanks again.

Original issue reported on code.google.com by jsandstrom42 on 2010-03-08 15:56:29

gforney commented 9 years ago
Using the most recently attached input file worked fine but changing H_FIXED from 5
to 50 produced numerical instability. Also setting H_FIXED=5 to the surrounding walls
produces numerical instability but changing H_FIXED on the walls from 5 to 50
increases the lenght of the calculation (still running).

I'll give it a try and see if I can simplify the error. Seems though as some of the
errrors requires large calculations to be produced...

The calculations where performed with svn 5629.

Original issue reported on code.google.com by erikasandstrom@hotmail.com on 2010-03-15 09:59:42

gforney commented 9 years ago
If you are not using H_EDDY, then probably no need to use VAN_DRIEST=T.

Or, you might also want to try CHECK_VN=T on MISC.  Often, these issues end up being
Von Neumann instabilities, which FDS5 does not check for in LES mode by default.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-15 12:39:08

gforney commented 9 years ago
Does CHECK_VN=T solve the problem or tell us that this is the problem?

I'll insert CHECK_VN=T on MISC and let you know.

Original issue reported on code.google.com by erikasandstrom@hotmail.com on 2010-03-15 12:57:04

gforney commented 9 years ago
Well, I'd say both.  If the problem is due to VN instability, this should take care
of it.  But to be safe, you might also set VN_MIN=0.4 and VN_MAX=0.5.

Again, this probably only needed if you keep VAN_DRIEST=T.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-15 12:59:45

gforney commented 9 years ago
If I remove VAN_DRIEST, there should be no instability?

I removed VAN_DRIEST and added CHECK_VN.
I do another one setting VN_MIN and MAX too.

Original issue reported on code.google.com by erikasandstrom@hotmail.com on 2010-03-15 13:05:48

gforney commented 9 years ago
This is my guess.  VAN_DRIEST damps the eddy viscosity at the wall.  This can be
destabilizing.  Regardless, CHECK_VN is a safe thing to do.  It is cheap and catches
many problems when the viscosity is small.  It is default for FDS6 options because
the dynamic model has a range of viscosity from molecular to several times that of
the constant coefficient model.

Original issue reported on code.google.com by randy.mcdermott on 2010-03-15 13:10:13

gforney commented 9 years ago
(No text was entered with this change)

Original issue reported on code.google.com by randy.mcdermott on 2011-08-26 19:56:40