Pressure drop in tunnel fire simulations

[rmcdermo]

See discussion forum thread here

This issue is also related to #6738

Thermal plume cases in tunnels give significantly different results when using variable and constant cp.

[rmcdermo]

For now I want to focus on a thermal plume, without HVAC, and sort out the reason why the variable CP case is not converging. We also need to first be assured that the constant CP case is converging. So far my evidence for this is spotty. I need to put together a verification series.

[ingraban98]

In case it helps.. these are the two cases I tried: V038.fds uses cp(T) V058.fds uses constant cp But some of the other parameters were varied as well.

V038.txt V058.txt

[tkorhon1]

Well, I have (probably) found a thing that is related to this Issue. I have problems with the pressure solver in a really simple one mesh geometry. I have pressure oscillations with a single "tube" that has other end as OPEN. And things get worse, if I put HVAC (to the AMB or to a small volume that has OPEN). The calculation is no good anymore at all. Things were better some years ago. I had a model that modelled just a small piece of the ISO horizontal furnace. I had there on the left and right side HVAC fans that "stirred" the air inside the tube-like furnace model (full sample width, but just 20cm in z dircetion and some 60cm in the third direction => more or less like a tube). This did not work anymore with fds6.7.0. It did work spring 2014.

The readme file (also in the attached zip-file): Problems: FDS 6.7.0 and its pressure solver: fluctuations "Analytical solution" is guite simple for the flow. The pressure should be close to the ambient and the temperature of the air in the tube very close to the tube surface temperature. So, the air temperature is increasing monotonically => density decreases monotonically => calculate the mass of the air in the tube => mass loss => m_dot at the opening of the tube.

Next are just with OPEN boundary on the -x side, single mesh.

v1B: 1 K/s, dt=0.01s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s v1B2: 1 K/s, dt=0.001s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s v1B3: 10 K/s, dt=0.01s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise v1B4: 10 K/s, dt=0.001s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise v1B3G: Gravity -x, 10 K/s, dt=0.01s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise v1B4G: Gravity -x, 10 K/s, dt=0.001s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise

Next are using HVAC duct to connect the heated pipe to the outside. Outside is a small volume within the mesh + OPEN boundary.

v1B: 1 K/s, dt=0.01s, time_ave=F devc, dt_devc=0.1, t_end=60s, t_begin=-1s, rad.=T, h_c=500 v1B3: 10 K/s, dt=0.01s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise v1B4: 10 K/s, dt=0.001s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise v1B3G: Gravity -x, 10 K/s, dt=0.01s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise v1B4G: Gravity -x, 10 K/s, dt=0.001s, time_averaged=.false. devc, dt_devc=0.1, t_end=60s, t_begin=-1s, no noise

TimoK Hot_Tube_Issue.zip

[rmcdermo]

Try with CONSTANT_SPECIFIC_HEAT_RATIO=F and see if that helps.

@tkorhon1 As a developer, I expect more that just "it worked better several years ago". Do a bisect and see if you can pin down where things started to "go bad". These tunnel or tube cases are problematic specifically because we have not buttoned them up tight over the years.

[tkorhon1]

Well,

&MISC TMPA= 20, BNDF_DEFAULT=.FALSE., GVEC= 0.00,0.0,0.0, SIMULATION_MODE='LES' , HVAC_PRES_RELAX=0.1, HVAC_LOCAL_PRESSURE=.TRUE., NOISE=.FALSE., GVEC=-9.81,0.0,0.0, / CONSTANT_SPECIFIC_HEAT_RATIO=.TRUE., / HVAC_PRES_RELAX=0.1, HVAC_LOCAL_PRESSURE=.TRUE., /

So, the constant spec.heat ratio was the default (there are some "comment" lines after the "real" MISC namelist above). I did test the constant_specific_heat_ratio T or F, but this did not help. I also played around with the HVAC pres relax and local pressure, these did not solve the problem either. When I did these tests, my input file was a little bit more complicated (and a little bit larger geometry, i.e., more grid cells). But there the constant spec.heat.ratio did not have any effect (by looking by eye the different graphs in Excel). HVAC_PRES_RELAX has some effect. But it does not help much. Setting HVAC_PRES_RELAX to whatever, I did not get as good ansver than using just OPEN alone without any hvacs. But the plain OPEN case is not nice either.

The constant spec.heat.ratio is not important in my case, probably. The temperatures are moderate (no combustion).

I try to pin down the version, where things start to go not so nice. I was a little bit busy to do these tests before Christmas, and I just returned to work today. But, I'll will pin down when the things started to go bad.

Randy, where could I find the old Windows executables for different releases? It would be easier to start with these old release versions than using Github/Google Code svn to get the old source codes. I can run these examples in my laptop, they are quite fast to calculate.

https://github.com/firemodels/fds/releases: Next is the last one there: FDS-SMV 6.3.0, @rmcdermo rmcdermo released this on Oct 1, 2015

So, if things are not going nice with fds 6.3.0, where can I find older releases? And where can I find older source code? Is it in GitHub, i,e., I can use GitHub to get old source? Or are the older sources only in the Google Code (svn)?

Wbr, Timo

[sbenkorichi]

Hi Timo, Hope this help. https://bintray.com/nist-fire-research/releases/FDS-SMV#files

Best Salah

[ingraban98]

Hi Randy I don't believe that cp(T) vs. cp=const. will be the key to the issue. I ran another case (V070). It is exactly the same as V058.txt (as above). The only parameter I changed was CONSTANT_SPECIFIC_HEAT_RATIO=F

Apparently, there is a major difference if the LES scheme is changed from LES to VLES. The parameter cp(T)/const. has some influence, but it it not the key to a reduction of these oscillations. I hope, your work is not affected by the government shut-down and you'll have some time to spend on this issue. I continued to work on simulations of heat volume sources, but our project is going to run into trouble as the validity of all tunnel simulations appears questionable. What do you think?

[tkorhon1]

Hi Randy,

I found a bug (or feature) that made my old model not to work anymore. Well, the bug/feature was already in the old source, but the thing did not result numerical instability (with some velocity tolerance/pressure iterations parameters in the old days).

The bug/feature is in the wall.f90 and for surfaces, where CASE (SPECIFIED_TEMPERATURE) METHOD_OF_HEAT_TRANSFER

I did following fast fix and tested it (not committed, of course):

  IF (ONE_D%UW<=0._EB) THEN
     IF (SF%TMP_FRONT>0._EB) THEN
        ONE_D%TMP_F = TMP_0(KKG) + EVALUATE_RAMP(TSI,SF%TAU(TIME_TEMP),SF%RAMP_INDEX(TIME_TEMP))*(SF%TMP_FRONT-TMP_0(KKG))
     ELSE
        ONE_D%TMP_F = TMP_0(KKG)
     ENDIF
  ELSE
     IF (SF%TMP_FRONT>0._EB) THEN
        ONE_D%TMP_F = TMP_0(KKG) + EVALUATE_RAMP(TSI,SF%TAU(TIME_TEMP),SF%RAMP_INDEX(TIME_TEMP))*(SF%TMP_FRONT-TMP_0(KKG))
     ELSE
        ONE_D%TMP_F = TMP_0(KKG)
     ENDIF

!!!Timo: ONE_D%TMP_F = ONE_D%TMP_G ! If gas is being drawn from the domain, set the boundary temperature to the gas temperature ENDIF

So, I had TMP_FRONT (and ramp for that also) and these surfaces had also leak_path defined (zone leakage model used). The IF statement with ONE_D%UW<=0._EB is not very robust. These kind of things might be hidden somewhere else and for some other surface type. Below/attached are two figures. One with the original source ("not_corrected") and one with the above correction. And the fds file is also attached (as a txt-ending).

This did not solve my original problem (the pressure fluctuation in a simple tube). There I did not have any leakage, so this was not the problem. Probably.

BUT: Should ONE_D%UW<=0._EB work generally? Are the boundary condition such, that there can be a tiny normal component of velocity on solid boundaries?

Timo Tinysize_insu_v4C.txt

[mcgratta]

Here's an idea that I can try:

WGT = MAX(UW,0)*DT/DX
TMP_F = WGT*TMP_G + (1-WGT)*TMP_FRONT

In this way, cases where there is just a small amount of outflow, like from a leak, will not be treated as they are now. I'll try it.

[mcgratta]

      WGT = MAX(ONE_D%UW,0._EB)*DT*ONE_D%RDN
      IF (SF%TMP_FRONT>0._EB) THEN
         TMP_FRONT = TMP_0(KKG) + EVALUATE_RAMP(TSI,SF%TAU(TIME_TEMP),SF%RAMP_INDEX(TIME_TEMP))*(SF%TMP_FRONT-TMP_0(KKG))
      ELSE
         TMP_FRONT = TMP_0(KKG)
      ENDIF
      ONE_D%TMP_F = WGT*ONE_D%TMP_G + (1._EB-WGT)*TMP_FRONT

This code works well in your case, but I need to look at it more closely before committing because this BC has always been troublesome. There are just so many different situations where it matters. However, this seems like a much better approach to what we had before because it is continuous over the range of UW. It's a simple upwind BC.

[tkorhon1]

Hi Kevin!

Yes, your correction seems to work fine enough for the case, where I have leak on the fixed temperature (or actually, tmp ramp) surface.

For the case without a leak (open tube and same tube and hvac duct to amb/open volume) things has not changed. So, I see some (fast) oscillations in pressure and these result oscillations in the flow in/out of the heated tube. I will check this stuff later, right now I'm a little bit busy. I should check if this is related TMP_FRONT boundary or not. So, I should, for example, set some volumetric heat source to heat the gas in the tube.

Timo

[mcgratta]

OK, I've got plenty of time to look at this.

[mcgratta]

Timo, how important is the heat transfer coefficient of 500? This value is usually around 10. This is part of the problem because the temperature in the gas rapidly changes due to the high HTC. A small time step would probably help. Even if I create a tunnel in a single mesh with perfect velocity boundary conditions, the high HTC causes spurious pressure spikes.

[rmcdermo]

CHECK_HT is a flag you can set to automatically obey heat transfer stability constraints

On Mon, Jan 14, 2019 at 12:57 PM Kevin McGrattan notifications@github.com wrote:

Timo, how important is the heat transfer coefficient of 500? This value is usually around 10. This is part of the problem because the temperature in the gas rapidly changes due to the high HTC. A small time step would probably help. Even if I create a tunnel in a single mesh with perfect velocity boundary conditions, the high HTC causes spurious pressure spikes.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/firemodels/fds/issues/7040#issuecomment-454099941, or mute the thread https://github.com/notifications/unsubscribe-auth/AENrwRNtSofiOtYoaVzlpm1mxEU0uBGHks5vDMUKgaJpZM4ZPadc .

-- Sent from my iPhone

[tkorhon1]

Well, high h_c is not important. I wil check the CHECK_HT input.

I used high t_c in my very simple ISO furnace model to heat the gas rapidly. But I could also blow hot gas (and have OPEN holes for the air to get out) and have not so high h_c for the walls. But I need larger than "normal" h_c for the walls to get the gas temperature to be close to the ISO temp-time curve. I do not want to model the ISO furnace gas burners in a simple model.

And I do not want to use EXTERNAL_FLUX on the SURF line either. Sometimes the samples are not planar. For example, they might have stiffeners etc. So, the sample can see itself, so the view factor of the hot furnace walls is less than 1 (or is there also four*pi, but you get my point).

I try CHECK_HT and see, if it helps. I will see, what is the resulting time step and see, if this is the limitting factor for the time step.

Timo

[tkorhon1]

Well, I did some more tests using the simple tube case, same as above. And I set CHECK_HT=T. But this did not help. I was already using DT=0.xxxx so small, that check_ht did not make the time step smaller in "my best, i.e., smallest dt, cases.

Anyway, below are figures for the flows, densities and pressures fo the cases: Tube, one side opens directly to a small volume, where is an OPEN boundary

• DT= 0.01s, 0.002 s, and 0.0001 s Tube, one side has HVAC duct to a small volume, where is an OPEN boundary
• DT= 0.02 s, 0.001 s, and 0.0001 s

The plain open tube (without any HVAC) is not too bad. But there is pressure (and flow) oscillations. But these are more or less same with different DT.

The havc tube case is bad. The DT behaviour is odd. The DT=0.001 s seems to work reasonably, but the smaller DT=0.0001s is really bad.

I attach the input files (I used the same files, just changing DT accordingly)

hot_tube_duct_v1B3NEW3.txt hot_tube_open_v1B3NEW3.txt

And here are the results:

Fire Dynamics Simulator

Current Date : January 15, 2019 15:16:02 Revision : FDS6.7.1-85-g889eccd69 Revision Date : Sat Jan 12 15:17:48 2019 -0500 Compiler : Intel ifort 18.0.1 Compilation Date : Jan 15, 2019 12:56:24

MPI Enabled; Number of MPI Processes: 1 OpenMP Enabled; Number of OpenMP Threads: 1

MPI version: 3.1 MPI library version: Open MPI v2.1.3, package: Open MPI tstopi@espvm5m030.ad.vtt.fi Distribution, ident: 2.1.3, repo rev: v2.1.2-129-gcfd8f3f, Mar 13, 2018

Job TITLE : ISO furnace a tiny and fast model Job ID string : hot_tube_open_v1B3NEW3

   Time Step  609900   January 15, 2019  17:53:17
   Step Size:    0.100E-03 s, Total Time:      59.99 s
   Pressure Iterations:      1
   Maximum Velocity Error:  0.73E-06 on Mesh   1 at (   4   6   4)
   Maximum Pressure Error:  0.66E+02 on Mesh   1 at (  64   4   4)
   ---------------------------------------------------------------
   Max CFL number:  0.40E-01 at (   4,   5,   4)
   Max divergence:  0.30E-01 at (   2,   6,   4)
   Min divergence: -0.77E-01 at (   2,   4,   4)
   Max VN number:   0.61E-04 at (   3,   5,   5)

   Time Step  610000   January 15, 2019  17:53:19
   Step Size:    0.100E-03 s, Total Time:      60.00 s
   Pressure Iterations:      1
   Maximum Velocity Error:  0.44E-06 on Mesh   1 at (   4   6   4)
   Maximum Pressure Error:  0.13E+03 on Mesh   1 at (  64   4   5)
   ---------------------------------------------------------------
   Max CFL number:  0.38E-01 at (   4,   5,   4)
   Max divergence:  0.29E-01 at (   2,   6,   4)
   Min divergence: -0.75E-01 at (   2,   4,   4)
   Max VN number:   0.56E-04 at (   2,   5,   5)

Time Stepping Wall Clock Time (s): 9436.691 Total Elapsed Wall Clock Time (s): 9441.909

STOP: FDS completed successfully (CHID: hot_tube_open_v1B3NEW3)

Well, I'm not sure about the revision string, I updated my source on Tuesday, so it had Kevin's fix for the tmp_front boundary in wall.f90. I compile on our Linux cluster, but my Git is on my laptop. So, I change the revision string manually in the makefile that is in Linux.

Timo

[rmcdermo]

Timo,

Given that this is my assigned issue, I'd like to not add heat transfer or HVAC here. It is getting too complicated. This particular issue is focussed on variable Cp differences in a tunnel with a thermal plume.

I suggest starting a new issue and assigning to Jason to have him set up a verification case for HVAC with DT convergence study.

Thanks, Randy

[ingraban98]

Randy The original question was about pressure fluctuations in tunnel simulations with plume. I am not convinced that the variation of cp is the root of the pressure fluctuations. (Although I am not a developer, but just a user in trouble.) Thanks for your help, Ingo

[rmcdermo]

@ingraban98 Fluctuations don't bother me. The flow is quasi-transient, as are all LES calculations. What bothers me is (1) that the solution does not converge and (2) that CP(T) gives drastically different results from CP even though CP(T) is nominally constant. This is a problem, which, if we do not understand, means we cannot hope to have reliable tunnel simulations.

[ingraban98]

Just a shy question: Do you make some progress on this issue? I'd be able to test new versions of the code, although I am not the software developer.

[rmcdermo]

Sorry, I've had no time to devote to this.

[ingraban98]

I am involved in a research project investigating the pressure drop of tunnel fires in longitudinal ventilation. The project involves empirical models and FDS for numerical simulation. Because of this issue, the project will be paused until July. I sincerely hope that the issue can be resolved. Or would you recommend, we should stop using FDS altogether and look for a different model?

[rmcdermo]

I have laid out the path to addressing this issue above. I will work on it as soon as I can. But I have several other issues pressing through June. I can't promise the issue will be resolved by July. Whether you look for another model is up to you.

[ingraban98]

2018_Capabilities_and_Limitations_of_the_FireDynamics.pdf

[msteck90]

Dear developers We are monitoring this issue quite for a while, as we are using FDS for various short- and longterm studies. According to our national tunnel ventilation guideline, we mainly execute our simulations and investigations with fires around 30MW. The solution of this fundamental FDS issue is of outmost importance to us to get reliable results. From the FDS release notes, one can see, that version 6.7.2 is currently under construction. Can you already say, whether the problem will be solved by this version and give a rough statement about the release date? Thanks a lot for your effort.

[rmcdermo]

I am not back in the US until June. I doubt this issue will be resolved by the next maintenance release. I expect another release in the summer sometime. There is a PhD student I met at the Juelich summer school who is going to try to help resolve this. Stay tuned.

[rmcdermo]

This issue was address in this letter to the editor in Fire Technology.