5. Numerical Examples 107
5.5. Non-Isothermal Laminar Flow: 2D Heated Cavity
The heated cavity is an example for flow in a rectangular enclosure, that is relevant for many industrial applications like cooling of electronic devices or venting in buildings, see e.g. [MH96] for analysis of such a problem. The effects of convection resulting from a heated wall in presence of a gravitational force can be examined.
5.5.1. Features of the Test Case
We consider laminar, non-isothermal flow in a cavity. The left vertical wall of the two-dimensional domain Ω = (0,1)2 is heated, the right one is cooled. The flow is then driven by this difference, no external forces occur. Due to the gravityg≡(0,−1)T, the fluid rises at the warm wall and sinks at the cold one.
Dirichlet boundary conditions for the temperature θleft = 0.5, θright = −0.5 and homo-geneous Dirichlet boundary conditions for the velocity are applied. We use the Prandtl number Pr= 0.71 of air at room temperature and vary the Rayleigh numberRa.
We solve the following dimensionless formulation of the Oberbeck-Boussinesq equations:
∂tu−Pr∆u+ (u· ∇)u+∇p−Pr Ra θey = 0 in (0, T)×Ω,
∇ ·u= 0 in (0, T)×Ω,
∂tθ−∆θ+ (u· ∇)θ= 0 in (0, T)×Ω
with Prandtl and Rayleigh numbers as introduced in Definition 2.1.1. The variables are of the order of the respective characteristic quantities.
The mesh is adapted to resolve the boundary layer and is randomly distorted by 1% (see appendix, Figure B.14). Here, the equidistant grid points are transformed by mappings Tab: [0,1]2 →[0,1]2 of the form
(x, y)T 7→Tab((x, y)T) :=
x− 1
2π(1−a) sin(2πx), y− 1
2π(1−b) sin(2πy) T
with parameters 0< a, b <1 chosen asa≈Nu−1 andb≈Nu−1/3, as suggested in [Löw11].
The so-called Nusselt numberNuis introduced in the section below; for the grid generation, we use reference values from the literature. We point out that in case of Ra = 104, the maximal aspect ratio of the cells is approximately 1 : 2.8. For the highest Rayleigh number considered, the maximal aspect ratio is nearly 1 : 50.8. The aspect ratio of neighboring cells does not differ much. Note that the analytical results from Chapters3and 4hold for isotropic meshes. This numerical test case goes beyond this. We refer the reader to [Ape99], where the approximation properties of finite element spaces on anisotropic meshes are
studied and local error estimates are derived.
Throughout this section, we use (Q2/Q1)∧Q1∧(Q2/Q1) elements.
5.5.2. Description of Benchmark Quantities
In order to validate the numerical results, we evaluate certain benchmarks and compare them with results from the literature.
The heat flux is defined asq:=uθ−α∇θ. In this example, we are interested in the heat transport from the warm left wall to the cold one, so we consider the horizontal heat flux qx(x, y, t) =uxθ−α∂θ∂x. From the total horizontal heat flux
hqxiy(x, t) :=
Z 1 0
qx(x, y, t)dy, (5.6)
the dimensionless Nusselt number Nu is calculated as the ratio between the convective heat transport h and the diffusive heat transport for motionless fluid in an areaA:
Nu(x, t) = h(x, t)L
α withh(x, t) =hqxiy(x, t)
A∆θ . (5.7)
In the dimensionless context, this reduces toNu(x, t) =hqxiy(x, t). It is an easy calculation that Nu(x, t) is independent of x when a steady state is reached. Using the stationary equation for θand the fact that u is solenoidal gives
∂hqxiy(x) due to the Dirichlet boundary conditions ofu andθ.
For two-dimensional flow, the streamfunction ∇Φ := u⊥ = (−uy, ux)T can be defined.
Local extrema indicate centers of a vortex and can be used as benchmarks. We solve for Φ numerically through the associated Poisson problem
∆Φ =∇ ·(∇Φ) =∇ ·u⊥ in (0, T)×Ω, Φ = 0 in (0, T)×∂Ω,
as described in [Löw11]. We compare the absolute values of the global extrema, denoted by Φmax, with values from the literature.
5.5.3. Numerical Experiments
The velocity magnitude and temperature for different Rayleigh numbersRa∈ {104,106,107} are presented in Figure 5.16; in the appendix, Figure B.15 all Ra ∈ {104,105,106,107} are depicted as an overview. No stabilization and N = 642 cells are used. In all of these cases, the flow is laminar and the solution reaches a stationary state that is shown in the pictures.
Figure 5.16.: Heated cavity, temperature (left) and velocity magnitude (right) with stream-lines, without stabilization,Pr= 0.71,Ra∈ {104,106,107},N = 642. For small Rayleigh numbers, a central main vortex of highest velocity occurs. For higher Ra, the boundary layer in the temperature near the heated and cooled walls becomes
thinner and the fluid in the middle of the cavity becomes almost motionless. Boundary layers within the velocity also emerge with growingRaas the flow near the walls becomes faster.
Not only this behavior but also the benchmark quantities in Table5.1are in good agree-ment with the literature [Löw11,LQ91]. Even for Ra = 107, no stabilization for neither velocity nor temperature is needed. Note that for our simulations, we use a grid that is transformed in a way consistent with [Löw11] but distorted randomly. We observe grid convergence of the Nusselt number and of Φmax. Furthermore, Nuavg is approximated better than Nu(0.5), as expected, since Nu(x) for some fixed x ∈ [0,1] depends on how well the temperature gradient near the walls is resolved. However,Nuavg is averaged over x∈[0,1] and is subject to the temperature difference only:
Z 1 0
∂θh
∂x dx=θh(1, y)−θh(0, y) =θleft−θright.
Nonetheless, we observe that Nuavg and Nu(0.5) are very similar. In case of the exact solution, they would coincide as the solution reaches a steady state.
Ra √
N (Nuavg)ref Nuavg Nu(0.5)ref Nu(0.5) (√
RaΦmax)ref √
RaΦmax
104 16 - 2.24478 - 2.23529 - 5.07359
32 2.24481 2.24481 2.24195 2.24436 5.07367 5.07367 64 2.24482 2.24482 2.24410 2.24656 5.07367 5.07367 128 2.24482 2.24482 2.24464 2.24719 5.07367 5.07367
105 16 - 4.52124 - 4.52350 - 9.60613
32 4.52162 4.52162 4.52283 4.52194 9.61490 9.6151 64 4.52163 4.52163 4.52192 4.52095 9.61570 9.61572 128 4.52164 4.52164 4.52170 4.52074 9.61637 9.61638
106 16 - 8.81573 - 8.85490 - 16.8347
32 8.82502 8.82502 8.82484 8.82173 16.79189 16.7934 64 8.82519 8.82519 8.82546 8.82189 16.81011 16.8103 128 8.82520 8.82520 8.82530 8.82174 16.81013 16.8102
107 16 - 15.3718 - 16.0458 - 43.1554
32 16.51578 16.5156 16.4921 16.4710 30.18100 30.1805 64 16.52302 16.5230 16.5211 16.5017 30.16094 30.1612 128 16.52309 16.5231 16.5228 16.5033 30.16377 30.1638 Table 5.1.: Nusselt numbersNuand global extrema of the streamfunction Φmax for
differ-ent Ra and different numbers of cellsN (at stationary point), compared with simulations by [Löw11].
For our purpose to study stabilization, the case of even higher Rayleigh numbers is more interesting. [LQB98] investigates the transition point, when unsteady and even turbulent flow occurs inside a differentially heated cavity. The critical Rayleigh number is found to be Racrit∈[1.81·108,1.83·108], where the flow becomes unsteady.
We consider a flow withRa= 108, which is slightly below this critical value. In Figure5.17, temperature and velocity are displayed without and with grad-div stabilization. N = 642 cells are used within the cavity. Although a stationary solution is expected, the unstabilized results are unstable and show unphysical oscillations and vortices in both the velocity and the temperature; no steady state is reached. Grad-div stabilization withγM = 1 cures this situation.
Figure 5.17.: Heated cavity with Pr = 0.71 and Ra = 108, temperature (left) and ve-locity (right) with streamlines without (top) and with grad-div stabilization (bottom),N = 642.
Table 5.2 compares the resulting benchmark values for different stabilization scenarios with [LQ91]. The presented values are obtained via averaging in time. Again, in all cases, we observe that refinement of the grid yields better agreement with the reference data.
There is a considerable difference in the quality of the benchmarks if no stabilization is
√
N 16 32 64 reference
Nuavg 18.7694 27.3830 30.0965 30.225 Nuavg GD 24.5473 27.9058 30.2237
Nuavg SUt 24.7862 27.9774 30.2238 Nuavg SUut 26.5130 28.2333 30.2239
Nu(0.5) 23.4069 29.2566 29.5245 30.225 Nu(0.5) GD 30.9351 27.1668 30.1839
Nu(0.5) SUt 29.5851 27.5024 30.1838 Nu(0.5) SUut 27.4597 28.7397 30.1825
√
RaΦmax 113.7240 96.2385 55.3487 53.8475
√
RaΦmax GD 108.2694 88.4510 53.8126
√RaΦmax SUt 109.1343 88.3693 53.8135
√
RaΦmax SUut 104.7587 87.1177 53.8136
Table 5.2.: Nusselt numbers Nu and extrema of the streamfunction Φmax for Ra = 108 and different numbers of cells N (averaged over t ∈[0.07; 0.1]) with different stabilizations, compared with benchmarks by [LQ91], that are obtained by using N = 1282 cells. GD denotes (γM, τMu, τLθ) = (1,0,0), SUt (γM, τMu, τLθ) = (1,0,kuhk∞,M−2), SUut (γM, τMu, τLθ) = (1,12h/kuhk∞,M,kuhk−2∞,M).
applied compared to the stabilized case. Adding LPS SU stabilization for temperature and velocity yields slightly better results than grad-div stabilization only. Note that anisotropic grids are used that are adapted to resolve the boundary layers.