Following the announcement of the v6 beta release, we
thought it might be helpful expand a little bit on some of the topics in the
release notes. Our aim here is to
explain why certain changes were necessary or why certain choices were made in
a way that might provide more insight than can be gleaned from just reading the
formal documentation (though reading the guides is still highly encouraged!).
The focus in this article will be hydrodynamics and
turbulence. Hydrodynamics is the
foundation of any CFD code, especially a reacting flow code like FDS. Chemical species cannot react until they have
mixed, and (in a diffusion flame) they cannot mix until they have been
transported. The velocity field,
therefore, plays a critical role in fire dynamics. In a model like FDS, other things being
equal, the model for the turbulent viscosity determines the behavior of the
velocity field and hence to a large extent the dynamics of the fire.
In v5, the turbulent viscosity was modeled by constant
coefficient Smagorinsky (1963) (csmag), which has the drawback being overly
dissipative (have you ever seen a 10-foot fire look like a candle flame?). Thus, finer grid resolution was required to
achieve accurate results. But even
highly resolved flames did not look qualitatively correct. The reason: csmag is not convergent. That is, the modeled term does not go away
(as it should) when the grid is refined to the level of a DNS.
The original candidate to replace csmag was the dynamic
Smagorinsky model (dsmag). While this
model added about 20% to the cost of computing the velocity field, in theory
this cost could be recovered by achieving more accurate results (statistically)
with coarser grid resolution. So, dsmag
was implemented and run through a battery of verification and validation tests
with more or less satisfactory results.
On the plus side, we were getting very realistic-looking flames at
moderate grid resolution. But on the
minus side, if the resolution was too coarse (often inevitable in engineering
calculations), the flames were erratic---instead of a 10-foot flame looking
like a candle, dsmag might produce an unstable ball of fire.
The solution to this dilemma came in the form of Deardorff’s
model (1972), with a twist. In the model
Deardorff originally proposed, he solved a transport equation for the subgrid
kinetic energy (sgs ke). This strategy
is expensive, but it allows for the inclusion of complex subgrid physics, like
unresolved buoyancy production. To avoid
this cost, we employ a simple algebraic model for the sgs ke based on the scale
similarity model of Bardina (1986). The
result is a model that is cheap and performs reasonably well at both coarse and
fine resolution. As an added benefit,
the framework is now in place to utilize a transported sgs ke, a topic that
might be worth exploring as FDS gets pushed to model large-scale outdoor
flows in both wildfire and wind engineering applications.
While the new turbulence model is certainly significant, the
change that really defines v6 is the new scalar transport scheme. The ripple effects of this modification were
not fully appreciated until a couple of years after the initial
implementation. What follows is the “saga
of removing spurious temperature wiggles.”
Take FDS 5 and simulate a simple fire plume with open
boundaries and 20 C ambient air temperature.
Now look closely at the temperature field near the corner of the
burner. You will see unphysical
excursions of the local temperature to well below ambient, near 0 C. Clearly this is undesirable. What is not as obvious is that there are also
unphysical excursions well above ambient, within the fire. Of course, 1020 C versus 1000 C is not as
dramatic a problem.
So why does this happen?
The root cause is central differencing in the solution of the density
equation. Purely centered schemes are
notorious for generating dispersion errors---wiggles. To combat this problem, v5 uses a boundedness
correction which tries to prevent the scalar fields (density and species
concentration) from going above or below defined limits. With mass fractions the limits are easy to
set (0 and 1). With density there is
only one hard limit (0).
In v6, we have decided to try another approach: total
variation diminishing (TVD) scalar transport.
This is a fancy way of saying “we use just the right amount of upwinding.” Pure upwind schemes are too dissipative
(translation: inaccurate). But TVD
schemes are specially designed to track scalar discontinuities with minimal
dispersion error and minimal dissipation.
This effectively solved the temperature wiggles problem.
Several TVD schemes are implemented in v6. The default for LES is Superbee (Roe, 1986),
so chosen because this scheme does the best job preserving the scalar variance
in highly turbulent flows with coarse grid resolution. The default scheme for DNS is Charm (Zhou,
1996) because the gradient steepening used in Superbee forces a stair step
pattern at high resolution, while Charm is convergent. A few other schemes (including Godunov and
central differencing) are included for completeness; more details can be found
in the Tech Guide.
These modifications were completed by the summer of 2009
(yes, you are reading that right). It
seemed at this point that we were very close to a v6 release. Then… the wheels came off the bus. In what we thought would be a routine pass
through the verification suite, we hit a snag: our energy budget cases were
off, by about 10%. In addition, several
of the validation cases now showed low upper layer temperatures---yep, by about
10%.
Hmmm. A lot of head
scratching ensued---two years, in fact, until Jason Floyd finally developed a
test case that pinpointed the problem.
Imagine a T-mixer in which equal mass fluxes of a hot gas with a low cp
(heat capcity) and a cold gas with a high cp mix and exit together. The critical aspect of this case is that the
flow paths are one cell thick and we artificially set the diffusivity and conductivity
to zero. Thus, any mixing is purely
numerical. What we found was that the
outlet stream temperature was mass weighted instead of enthalpy weighted---there
was no accounting for the variation in cp.
The pieces of the puzzle finally started to come together because
another relatively recent development (actually in v5) was the inclusion of
variable temperature specific heats.
The aftermath of this discovery was a torturous dissection
of the energy equation (see appendix of the Tech Guide) illustrating that basic
assumptions of the low-Mach number approximation can hide the effects of the
numerical mixing from TVD transport schemes. The bright side of the story is that through
this analysis we were able to derive corrections for numerical mixing that
force FDS to satisfy the discrete, conservative form of the sensible enthalpy
equation, ultimately solving the temperature issues.
The benefit of the delayed release is that many other parts
of the code were improved in the mean time.
So stay tuned for our next installment.