OJFD  Vol.9 No.3 , September 2019
Numerical Analysis of the Transient Process of Flow and Heat Transfer at High Rayleigh Number in a Partially Open Habitat Heated from Above
ABSTRACT
In this work, we studied the thermoconvective instabilities in a pentagonal cavity containing a Newtonian fluid. The cavity provided with a side opening is uniformly heated from above by a constant heat flux. The natural ventilation phenomenon in the classic habitat of the hot climate is thus numerically analyzed with unsteady natural convection equations formulated with vorticity and stream-function variables. The finite volume predictions of two-dimensional laminar natural convection at high Rayleigh number are presented. Results show that the incoming fresh air and the hot air discharge begin with the late start of the convection. The phenomenon intensifies with time and the birth of instabilities improves the homogenisation of temperatures which imply the elimination of very cold and very hot areas. However, the competition between the incoming fresh air and the hot air expansion leads to a perpetual displacement of the thermal front. The cross-sections at the opening of the incoming fresh and outgoing hot air are time-varying and the penetration depth of the fresh air is highlighted by the large convective cells originated from the aperture. The non monotonic variation of the Nusselt number reflects not only the multicell nature of the flow but also expresses the heat lost by the active walls due to the fresh air.

1. Introduction

In the absence of wind, the temperature gradient between the walls of an open cavity and the external environment causes a natural ventilation phenomenon which is the air movement through the openings. In doing so, by the phenomenon of convection, the apertures provide the thermal comfort and/or chemical decontamination of the enclosures. Interests of researchers are mainly focused on the natural ventilation control parameters in partially or fully open cavities.

The square or rectangular shape with a side opening is widely used to model the type of flow that occurs in household refrigerators and ovens under open-door conditions. Thus, works of Mokhtarzadeh-Dehghan et al. [1] in a square cavity opened at the top, showed that the characteristic number of Stanton is significantly influenced by the positions of the active wall and of the opening. Works such as [2] [3] conducted numerically and experimentally in a rectangular cavity with a side opening revealed that convection heat losses are lower than those of a vertical plate with height that is equal to the total length of the cavity. All three walls of the cavity were at a uniform temperature exceeding the ambient. The study [2] [4] [5] of the aperture position in angle ranges from 0˚ - 120˚, relative to the horizontal, revealed that the Nusselt number is a non-linear function of the inclination angle. For electronic cooling, Bilgen and Oztop [5] found that the large aperture at low position and inclined at 60˚ - 90˚ would maximize heat transfer. Heating from above is then found to be the worst case, because heat is mostly transferred by conduction. Several works such as [6] [7] [8] [9] on the natural convection in cavities essentially heated from above, showed that the laminar convection develops for Rayleigh number Ra > 106. Experimental works performed by Wu and Ching [10] in rectangular cavities heated from the top showed that the flow remains laminar for a global Grash of number of approximately 1.8 × 108. Gan [7] and Dubovsky et al. [8] have shown that the ventilation can be performed in a rectangular cavity heated from the horizontal top wall. Thus, Gan [7] has numerically shown that the effective depth for thermal comfort may not coincide with that of the distribution of the fresh air.

The recently numerical study of Kpode et al. [11] predicted the heat transfer rates and flow patterns in pentagonal cavity heated from above. This represented the thermal and dynamic fields of the habitat with a gable roof during a sunny day. The convection developed at high Rayleigh number and the flow became unstable at high time. However, the results of this study cannot be readily adapted to a real habitat of the hot climate that cannot be totally closed. Also the effect of the opening on the fluid dynamics in a habitat with a gable roof is rare. The first attempts are the analytical and experimental works conducted respectively by Hirunlabh et al. [12] and Al-Helal et al. [13]. They explored the two ventilation openings effects and found that the internal temperature that is still greater than the ambient temperature of the outside is correlated with the mass flow.

On the basis of the literature review, it appears that no wok was reported on natural convection in ventilated pentagonal enclosure heated mainly by the top which is nevertheless the most encountered configuration in very hot climate. Indeed, it is the habitat model (residential, commercial building, warehouse, etc.) heated by the roof during a long period of sunshine. Transfers dominated over a long period of the day by conduction, cause a strong thermal stratification, and therefore the stagnant air becomes more and hotter. The birth of the thermoconvective instabilities would have positive effects by improving mixtures and mass transfers. Thus, very hot and very cold areas will be removed. Hence, by providing the pentagonal cavity essentially heated from above with a large side opening, this study consists in using the more general fluid numerical methods such as CFD for the insightful exploration of the thermal and dynamic fields of the open classical habitat of the hot climate for energy efficiency. The state of the system over time for a large Rayleigh number will be analyzed and the impact of the opening on the transfers will be highlighted.

2. Mathematical Formulations

A schematic diagram of the system under consideration is shown in Figure 1. It consists of a side opening cavity heated from above by a uniform heat flux q. The side wall opposite the opening and the horizontal wall are, respectively, adiabatic and isothermal. Thus, we define the more real summer conditions where most of the heat introduced into the habitat comes from the roof [9] [14] [15]. The working fluid physical properties are the dynamic viscosity μ , density ρ , thermal conductivity λ and specific heat ρ c p . We assume that the flow is two-dimensional and that the physical properties of the fluid are constant unless its density in terms of gravity, which in a first approximation of Boussinesq hypothesis, varies linearly with temperature. Taking into account the above assumptions and in the stream function ψ and vorticity ω formulation, the non-dimensional equations of heat, vorticity and stream function are, respectively, as follows:

t θ + ( V θ θ ) = 0 (1)

t ω + ( V ω Pr ω ) = Pr Ra x θ (2)

( ψ ) = ω (3)

Figure 1. Schematic diagram of open cavity.

where ψ and ω are such that:

u = y ψ , v = x ψ and ω = x v y u , u and v are respectively horizontal and vertical dimensionless coordinates of the velocity V in x and y cartesian coordinates system. The reference parameters used to make the problem dimensionless are L , L 2 / α and ( q L ) / λ , which, respectively, represent the length, time and temperature gradient and the dimensionless temperature is such that: θ = λ ( T T 0 ) / ( q L ) . Thus, it appears in Equation (2) the dimensionless parameters Prandtl number, Pr and Rayleigh number, Ra defined respectively as follows:

Pr = ν / α , Ra = g β L 4 q / ( ν α λ ) ,

where ν = μ / ρ , α = λ / ( ρ c p ) , and β = ( 1 / ρ 0 ) ( ρ / T ) P

The above equations are complemented by the following initial and boundary conditions:

when t = 0

u ( x , y , 0 ) = v ( x , y , 0 ) = ψ ( x , y , 0 ) = ω ( x , y , 0 ) = θ ( x , y , 0 ) = 0 ,

when t > 0

hydrodynamic conditions on the walls

ψ = ( n ) ψ = 0 and ( n ) [ ( n ) ψ ] = ω ,

thermal conditions on the inclined walls

( n ) θ = 1 ,

thermal conditions on the right side wall

( n ) θ = 0 ,

thermal conditions on the bottom wall

θ ( x , 0 , t ) = 0 ,

where n is the external normal vector to each wall.

Conditions at the opening x = 0 and 0 y A ; where A is the dimensionless aperture (aspect ratio):

In the case of the transfers by natural convection in the open cavities, the major difficulty lies in the treatment of the magnitudes at the openings, which are unknown. To deal with these difficulties, Kettleborough [16] , Quere et al. [17] ), Penot [2] , extend the physical domain to a fictitious computational domain outside the opening. These changes, enabled them to find the minimum distances away from the opening beyond which, the imposed physical conditions do not significantly affect the natural behaviour of the flow at the openings. Thereafter, Gan [18] found that the right distances are ten times the length of the cavity for the asymmetric heating, and five times more for the symmetrical heating. Comparing the results obtained by the extended domain with the restricted domain, Chan and Tien [19] found that at high Rayleigh number, the restricted approach gives good heat transfer results and the basic flow patterns inside the cavity. The simple approach by restricted domain is used and the conditions [4] at the opening are as follows:

x v = y u = x ψ = 0 ,

ω = 0 ,

θ = 0 if u 0 and x θ = 0 if u < 0 .

The heat energy transmitted by the inclined active walls is characterized by the Nusselt number. The base or the ambient air temperature ( θ = 0 ) is used. Thus, the local and the mean Nusselt numbers of an active wall are:

N u = 1 θ P ,

N u ¯ = 1 S S N u d s

where θ P is the instantaneous local temperature of the wall and S the length of an active wall. The instantneous volume flow rate [5] is given by:

V ˙ = 1 2 0 A | u ( 0 , y ) | d y

3. Numerical Model

An implicit Euler scheme is used for time integration and the finite volume method [20] [21] is chosen for the spatial integration. Each control volume is centred on a node and a power-law scheme of Pantakar [20] [22] is used to approximate the coefficients of algebraic equations. The set of equations is solved by the iterative under-relaxation line-by-line method of Gauss-Seidel [23]. At each instant, a variation by less than 10−5 over all the grid points for each variable is considered as the convergence criterion. Thus, the computation code based on these above numerical schemes is firstly used to perform the tests with the isosceles right triangle under the work conditions of Flack [6] and Ridouane et al. [9]. Indeed, the numerical experiments that we conducted with the Rayleigh numbers less than 5 × 106 show that the transfers are mainly in the area adjacent to the active inclined walls as reported by Flack [6] and Ridouane et al. [9]. The rest of the domain remains stagnant. Thus, study in this cavity can be restricted in or under these conditions to the attic space by setting aspect ratio equal to zero as done by these authors. Secondly, we perform the tests with the partially open square cavity at Ra = 104 with aperture size AR = 0.5 and for aperture position AP = 0.5 when the heated wall is inclined φ = 90˚ from the horizontal is used as reported [5]. The comparison of our results (see Figure 2 and Figure 3(a), Figure 3(b)) with the findings of these authors when the stationary state is reached confirm the reliability of the computation code.

Spaces and time steps are chosen to satisfy the necessary conditions of Courant Friedrichs-Lewy (CFL) [24] [25]. The uniform 121 × 121 grid with the time step ∆t = 8 × 105 is then selected (see Figure 4) after the sensitivity tests performed with the mean Nusselt number calculated on the right inclined wall.

Figure 2. Temperature of the mid-plane of the triangular enclosure heated from above at Gr = 2.84 × 106: Comparison between the results by Flack [6] and Ridouane et al. [9] with the present numerical study.

(a)(b)

Figure 3. Stream functions (left) and isotherms (right) with Ra = 104, AP = AR = 0.50, φ = 90˚: comparison between the numerical results by Bilgen and Oztop [5] Figure 3(a) and the present numerical study Figure 3(b).

4. Results

The upper walls are inclined at an angle of 45˚ relative to the horizontal plane and the aspect ratio A = h/L = 0.5. The results are relative to the air whose physical properties match a Prandtl number of 0.7 and the Rayleigh number; Ra = 1 × 108. The streamlines and the isotherms are respectively incremented by Δ ψ = ( ψ max ψ min ) / n and Δ θ = ( θ max 0.02 ) / n with n, the number of increments.

Figure 4. Sensitivity of the space and time steps with the mean Nusselt number calculated on the right inclined wall.

4.1. Flow and Isotherm Patterns

The isotherms depicted in Figure 5(a) with their shapes matching the inclined walls during the first moments (t < 175), prove that the transfers are dominated by conduction even though the velocity field is non-zero. As observed in some cavities heated from above [9] [14] [15] the fluid layer in the immediate vicinity of the active walls, is driven to the confined apex angle, whereas the relatively heavy fluid near the median plane moves downward. It results two symmetric main vortices which develop and rotate in opposite directions. The rest of the fluid in thermodynamic equilibrium with the outside region of the opening remains stagnant. The temperature gradient between the outdoor and indoor is not yet sufficient to enable the aspiration of the surrounding air.

The increasing deformations with a sudden depression to the median line of the isotherms, show that the heat gradually affects the whole cavity. The transfer by convection thus develops and induces the inlet and outlet air phenomenon. The symmetry of the flow is then broken. However, for t < 225, the symmetry is preserved (see Figure 5(a), Figure 5(b)) in the attic space (triangular section). From t ≥ 175, the secondary cells form in the region near the intersection of the sloping walls just above the main vortices.

From t = 225, the diffusion depth of the incoming fresh air concomitantly increases with the appearance and development of a second counter-rotating secondary cell (dashed line) (Figure 5(c), Figure 5(d)) in the upper right corner. The main vortices (Figure 5(d), Figure 5(e)) which are thereby trapped between the secondary cells, decrease in volume and move towards the walls to yield to the rise of the stable secondary cells formed in the apex angle. The development of the trigonometric secondary cells originating from the opening are detrimental to the second counter-rotating secondary cell located in the right side of the cavity. Then, it gradually disappears and yields to the rise of the single flow of trigonometric cells with multiple cores in the rectangular part of the cavity

(a) (b) (c) (d) (e) (f) (g) (h) (i) (j) (k)

Figure 5. Evolution of streamlines (left) and isotherms (right). (a): t = 50, ψ max = 1.92, ψ min = −1.92, θ max = 0.17; (b): t = 175, ψ max = 8.96, ψ min = −8.51, θ max = 0.26; (c): t = 225, ψ max = 9.01, ψ min = −8.39, θ max = 0.33; (d): t = 275, ψ max = 8.26, ψ min = −7.96, θ max = 0.43; (e): t = 400, ψ max = 16.54, ψ min = −14.91, θ max = 0.62; (f): t = 500, ψ max = 24.55, ψ min = −22.55, θ max = 0.64; (g): t = 650, ψ max = 34.28, ψ min = −32.08, θ max = 0.55; (h): t = 900, ψ max = 66.41, ψ min = −46.68, θ max = 0.22; (i): t = 1500, ψ max = 82.68, ψ min = −65.95, θ max = 0.14; (j): t = 2000, ψ max = 60.30, ψ min = −60.83, θ max = 0.14; (k): t = 2250, ψ min = 61.53, ψ min = −45.50, θ max = 0.15.

(Figures 5(e)-(g)). This shows the renewal of the air of almost the whole rectangular part at t = 650 (Figure 5(g)). The depth of the fresh air increases with the transfers as confirmed by Gan [7]. Nevertheless, this flow is very unstable facing the very hot counter-rotating cells and the stable trigonometric cells formed in the apex angle. The latter flow down and split the trigonometric flow below (Figures 5(h)-(j)). The system is then restructured with the resorption of no stable cells by the more stable (Figure 5(k)). Also, because the fresh air increasingly invades the cavity, it opposes the heat downward migration, thus causing a tightening of the isotherms over time.

4.2. Temperature and Velocity at the Opening and Volume Flow Rate

Figures 6-8 display the temperature variations and velocity profiles over time at the opening. In time, the portion of the opening where the temperature is equal to the external ambient temperature, decreases and reaches the limit space point (x = 0, y = 0.35) that is greater than the height of the half-opening (Figure 6). The temperatures of the points above this limit position, increase monotonously with time, but from a certain date, oscillate with time due to the appearance of the more or less warm new recirculation cells at the opening as illustrated in Figure 7. Then, their decrease towards the ambient temperature shows that these points are crossed by less and less warm cells. This limit space point also corresponds to the position where the velocity (see Figure 8) still zero at high time. The increasing variation of the flow rate (Figure 9) shows that the air intake phenomenon intensifies with the transfers. The constant value is not reached in time and the oscillations indicate the deceleration of the air intake and the reorganisation of the system.

4.3. Mean Nusselt Numbers

The curves in Figure 10 show the time evolution of the mean Nusselt numbers

Figure 6. Temperature profile at the opening for different moments.

Figure 7. Time evolution of the temperature at the points (x = 0, y > 0.35).

which are substantially equal at the first moments, and their strong decrease corresponds to the pseudo-conduction mode and to the existence of the two symmetric main cells in the attic space (see Figure 5(a), Figure 5(b)). When the convection is established at high time, and in the reference to the Nusselt number ( N u r e f ¯ ) calculated on an inclined wall in the closed case, it is clear that the effect of the fresh air, whose depth extends to the right inclined wall, is seen on the variations of the Nusselt numbers and mainly on the Nusselt number ( N u r ¯ ) calculated on the right inclined wall opposite the opening. Therefore, the opening has a little effect on the heat transfer on the left active wall. Indeed, the incoming fresh fluid feeds the transfers and increases heat loss in the vicinity of the active walls. Thus, as the penetration depth of the fresh air increases, the thermal

Figure 8. Velocity profile at the opening for different moments.

Figure 9. Time evolution of the volume flow rate.

boundary layer gets thinner and the mean Nusselt numbers grow to a maximum and then drop smoothly as the thickness of the boundary layer increases.

5. Conclusion

A numerical investigation of two-dimensional unsteady laminar natural convection in a partially open pentagonal cavity was conducted with sunny boundary conditions and for a high Rayleigh number; 1 × 108. Under these conditions, the fluid flow in the cavity is unstable and enables understanding transfers by natural convection in classic semi-open habitat with a gable roof in a hot climate. When convection increases with the heat expansion at high time, phenomena of the incoming fresh air and hot air discharge begin. But, at the same time, the

Figure 10. Time evolution of the mean Nusselt numbers.

fresh incoming air opposes to the heat expansion. This perpetual competition causes the displacement of the thermal front. The volume flow rate which increases with the transfers, is a good parameter for assessing the indoor air change rate. The time when its maximum value is reached would be the right time for effective decontamination or efficient insecticide treatment. The unsteady natural convection allows predicting the adequate diameter of the fresh air intake opening or the polluted air discharge opening because the cross-section of fluids (fresh and hot fluids) at the opening is time-varying. The largest cross-section is reached before the quasi-steady state. The average Nusselt numbers show that the overall variations of the inclined walls temperature over time are at the origin of the self-organisation of the system due to the appearance and disappearance of the hot fluid cores.

Nomenclature

Nu: mean Nusselt number

V : dimensionless velocity vector

A: aspect ratio, = h/L

CP: specific isobaric heat capacity (J∙kg−1∙K−1)

Gr: Grash of number, = g β L 4 q / ( ν 2 λ )

h: side heigth (m)

hO: height of the attic space (m)

L: base (m)

Nu: local Nusselt number

Pr: Prandtl number, = ν / α

q: wall heat flux (W/m2)

Ra: Rayleigh number, Ra = Gr∙Pr

S: length of an incline wall (m)

t: dimensionless time

TO: temperature of the base wall (K)

u, v: horizontal and vertical dimensionless velocity coordinates

x, y: horizontal and vertical dimensionless coordinates

Greeks Symbols

α : thermal diffusivity (m2∙s−1)

β : coefficient of thermal expansion (K−1)

γ : angle (rad)

λ : thermal conductivity (W∙K−1∙m−1)

μ : dynamic viscosity (kg, m−1∙s−1)

ν : kinematic viscosity (m2∙s−1)

ω : dimensionless vorticity

ψ : dimensionless stream function

ψ min , ψ max : minimum and maximum value of ψ

ρ : density (kg/m3)

θ : dimensionless temperature, = λ ( T T O ) / ( q L )

Cite this paper
Kpode, K. , Nougbléga, Y. , Mbow, C. and Banna, M. (2019) Numerical Analysis of the Transient Process of Flow and Heat Transfer at High Rayleigh Number in a Partially Open Habitat Heated from Above. Open Journal of Fluid Dynamics, 9, 253-267. doi: 10.4236/ojfd.2019.93017.
References
[1]   Mokhtarzadeh-Dehghan, M., Telbany, M.E. and Reynolds, A. (1990) Transfer Rates in Single Sided Ventilation. Building and Environment, 25, 155-161.
https://doi.org/10.1016/0360-1323(90)90028-P

[2]   Penot, F. (1982) Numerical Calculation of Two-Dimensional Natural Convection in Isothermal Open Cavities. Numerical Heat Transfer, 5, 421-437.
https://doi.org/10.1080/10407788208913457

[3]   Skok, H., Ramadhyani, S. and Schoenhals, R. (1991) Natural Convection in a Side-Facing Open Cavity. International Journal of Heat and Fluid Flow, 12, 36-45.
https://doi.org/10.1016/0142-727X(91)90006-H

[4]   Mohamad, A. (1995) Natural Convection in Open Cavities and Slots. Numerical Heat Transfer, Part A: Applications, 27, 705-716.
https://doi.org/10.1080/10407789508913727

[5]   Bilgen, E. and Oztop, H. (2005) Natural Convection Heat Transfer in Partially Open Inclined Square Cavities. International Journal of Heat and Mass Transfer, 48, 1470-1479.
https://doi.org/10.1016/j.ijheatmasstransfer.2004.10.020

[6]   Flack, R. (1980) The Experimental Measurement of Natural Convection Heat Transfer in Triangular Enclosures Heated or Cooled From Below. Journal of Heat Transfer, 102, 770-772.
https://doi.org/10.1115/1.3244389

[7]   Gan, G. (2000) Effective Depth of Fresh Air Distribution in Rooms with Single-Sided Natural Ventilation. Energy and Buildings, 31, 65-73.
https://doi.org/10.1016/S0378-7788(99)00006-7

[8]   Dubovsky, V., Ziskind, G., Druckman, S., Moshka, E., Weiss, Y. and Letan, R. (2001) Natural Convection Inside Ventilated Enclosure Heated by Downward-Facing Plate: Experiments and Numerical Simulations. International Journal of Heat and Mass Transfer, 44, 3155-3168.
https://doi.org/10.1016/S0017-9310(00)00324-0

[9]   Ridouane, E.H., Campo, A. and McGarry, M. (2005) Numerical Computation of Buoyant Airflows Confined to Attic Spaces under Opposing Hot and Cold Wall Conditions. International Journal of Thermal Sciences, 44, 944-952.
https://doi.org/10.1016/j.ijthermalsci.2005.03.008

[10]   Wu, W. and Ching, C.Y. (2009) The Effect of the Top Wall Temperature on the Laminar Natural Convection in Rectangular Cavities with Different Aspect Ratios. Journal of Heat Transfer, 131, Article ID: 052501.
https://doi.org/10.1115/1.2993138

[11]   Kpode, K., Sow, M.L. and Mbow, C. (2016) Unsteady Natural Convection with Summer Boundary Conditions in a Habitat at High Rayleigh Number and at High Time. Energy and Buildings, 121, 72-77.
https://doi.org/10.1016/j.enbuild.2016.03.068

[12]   Hirunlabh, J., Kongduang, W., Namprakai, P. and Khedari, J. (1999) Study of Natural Ventilation of Houses by a Metallic Solar Wall under Tropical Climate. Renewable Energy, 18, 109-119.
https://doi.org/10.1016/S0960-1481(98)00783-6

[13]   Al-Helal, I., Waheeb, S., Ibrahim, A., Shady, M. and Abdel-Ghany, A. (2015) Modified Thermal Model to Predict the Natural Ventilation of Greenhouses. Energy and Buildings, 99, 1-8.
https://doi.org/10.1016/j.enbuild.2015.04.013

[14]   Asan, H. and Namli, L. (2000) Laminar Natural Convection in a Pitched Roof of Triangular Cross-Section: Summer Day Boundary Conditions. Energy and Buildings, 33, 69-73.
https://doi.org/10.1016/S0378-7788(00)00066-9

[15]   Kent, E.F. (2009) Numerical Analysis of Laminar Natural Convection in Isosceles Triangular Enclosures for Cold Base and Hot Inclined Walls. Mechanics Research Communications, 36, 497-508.
https://doi.org/10.1016/j.mechrescom.2008.11.002

[16]   Kettleborough, C. (1972) Transient Laminar Free Convection between Heated Vertical Plates Including Entrance Effects. International Journal of Heat and Mass Transfer, 15, 883-896.
https://doi.org/10.1016/0017-9310(72)90228-1

[17]   Quere, P.L., Humphrey, J.A. and Sherman, F.S. (1981) Numerical Calculation of Thermally Driven Two-Dimensional Unsteady Laminar Flow in Cavities of Rectangular Cross Section. Numerical Heat Transfer, 4, 249-283.
https://doi.org/10.1080/01495728108961792

[18]   Gan, G. (2010) Impact of Computational Domain on the Prediction of Buoyancy-Driven Ventilation Cooling. Building and Environment, 45, 1173-1183.
https://doi.org/10.1016/j.buildenv.2009.10.023

[19]   Chan, Y. and Tien, C. (1985) A Numerical Study of Two-Dimensional Laminar Natural Convection in Shallow Open Cavities. International Journal of Heat and Mass Transfer, 28, 603-612.
https://doi.org/10.1016/0017-9310(85)90182-6

[20]   Versteeg, H.K. and Malalasekera, W. (2007) An Introduction to Computational Fluid Dynamics: The Finite Volume Method. Pearson Education, London.

[21]   Ferziger, J.H. and Peric, M. (2012) Computational Methods for Fluid Dynamics. Springer Science & Business Media, New York.

[22]   Patankar Suhas, V. (1980) Numerical Heat Transfer and Fluid Flow. CRC Press, Boca Raton, FL.

[23]   Barrett, R., et al. (1994) Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods. Philadelphia, PA, xvii + 118.
https://doi.org/10.1137/1.9781611971538

[24]   Harlow, F.H. and Welch, J.E. (1965) Numerical Calculation of Time-Dependent Viscous Incompressible Flow of Fluid with Free Surface. The Physics of Fluids, 8, 2182-2189.
https://doi.org/10.1063/1.1761178

[25]   Srinivas, K. and Fletcher, C. (2002) Computational Techniques for Fluid Dynamics: A Solutions Manual. Springer Science & Business Media, New York.

 
 
Top