Abstract
A finite Fourier transform is used to perform both linear and nonlinear stability analyses of a Darcy-Lapwood system of convective rolls. The method shows how many modes are unstable, the wave number instability band within each mode, the maximum growth rate (most critical) wave numbers on each mode, and the nonlinear growth rates for each amplitude as a function of the porous Rayleigh number. Single amplitude controls the nonlinear growth rates and thereby the physical flow rate and fluid velocity, on each mode. They are called the flak amplitudes. A discrete Fourier transform is used for numerical simulations and here frequency combinations appear that the traditional cut-off infinite transforms do not have. The discrete show a stationary solution in the weak instability phase, but when carried past 2 unstable modes they show fluctuating motion where all amplitudes except the flak may be zero on the average. This leads to a flak amplitude scaling process of the heat conduction, producing an eddy heat conduction coefficient where a Nu-R
1. Introduction
Convection in porous media is intensively studied because of its many applications in science and industry. Free convection in porous media (i.e., convection without the forcing of horizontal temperature gradients) has been studied from 1948, [1–3] made the first experiments, and [4] treated the stability of 2-dimensional cellular flow (rolls). Others are investigating stability problems with increasing complexity, for example, [5]. Modal flow and the spectral method were used to simulate unsteady motion by [6].
The spectral method is used in [7] and it gives an excellent review of approximate solutions, sometimes based on earlier findings [8, 9]. In the review paper [7] the most general case, the DLFB (Darcy-Lapwood-Forchheimer-Brinkman) equation, is studied, sometimes with Coriolis force and magnetic forces being included.
In thick natural aquifers the viscous dissipation, not already included in Darcy-Lapwood (DL) systems, is unimportant in most cases and so is turbulent dissipation except in local irregularities in the porous matrix. Porous media flows on low Reynolds numbers are therefore mostly treated without the terms of Forchheimer and Brinkman; they serve as a bridge between the turbulent and viscous and Darcy regimes. This paper is devoted to the case of simple DL systems and it is shown that turbulent-like phenomena can be encountered in such flow even though dissipation is strictly Darcy-laminar. The effect of this turbulence is on the heat flow, amplitude growth rates, and the wave length spectrum and this is studied further.
The heat flow through a porous convective layer is controlled by Nusselt’s number Nu. Finding it is therefore very important in industrial applications. In [10] flow in a box with open surface is investigated numerically for Rayleigh numbers up to 300 and it is found that there may be more than one flow cell. The author finds Nusselt’s numbers in the range of – for this flow and in some cases for eddy flow. He also finds single-cell flow for but two cells for . Multicellular flow and modal flow are also encountered in other papers on the subject [11, 12].
Stability analyses of the no-flow situation and the onset of convection give important information on the structure of the resulting flow through the composition of the eigenfunctions. They depend upon the boundary conditions; therefore a great variety of eigenfunctions have been found for the rolls of the DL equation system covering different sets of boundary conditions [12–15].
The eigenfunctions of the DL system in an infinite horizontal aquifer of constant thickness and constant temperature difference are particularly simple as they are the ordinary trigonometric functions. This makes the Fourier transform of the spectral method be an expansion of the temperature function in a series of the eigenfunctions of the Sturm-Liouville problem underlying the stability analysis [6]. This makes the identification of the different modes and their stability limits particularly easy and in this paper it is shown to give rise to a new kind of non-linear stability analysis with the starting point in an arbitrary stationary solution. From this process the Fourier amplitudes emerge that make up the average vertical temperature gradient, here called the flak amplitudes. They alone control the exponential growth rate of all amplitudes and make it possible to construct by scaling an eddy coefficient of heat conduction, very similar to eddy viscosity in ordinary turbulence. This eddy coefficient of heat conduction allows the approximated solutions described by [7] to produce realistic relationship between the porous Rayleigh number and Nusselt’s number for DL systems and thus give better estimates for the heat flow through the porous layer. As the flak amplitudes appear also in DLFB systems, it may be suggested to use the scaling procedure for them also, but it is much more complicated and is not attempted in this paper. Simulation grids up to . .
2. The Darcy-Lapwood System
Two-dimensional nonlinear thermal convection in a homogeneous horizontal aquifer, uniformly heated from below and cooled from above, is described by a differential system originally presented in [16]. In [15] it is presented in cylindrical coordinates. Here we use the notations in [10]: , dimensionless temperature. is the difference of the constant temperatures of top and bottom. is dimensionless stream function. and are dimensionless coordinates in horizontal and vertical direction. Lapwood (porous) Rayleigh number, where is heat capacity of the porous layer , is coefficient of permeability (m/s), is fluid density difference (kg/m3) corresponding to , is thickness of aquifer (m), and is coefficient of heat conduction for the porous layer.
3. The Spectral Method Using Infinite Spectrum
3.1. Equations
The boundary conditions are at , , and at ; then the temperature distribution in a cellular flow is given by the following Fourier series: is the basic horizontal wave number, is horizontal wave number, and is called the mode. Individual amplitudes must satisfy the following equations: Equation (4) indicates exponential growth in the first term; (5) shows the growth rate. The last term in (4) is a sum, that is, a quadratic form in all the temperature amplitudes in (3). This simple form emerges because the trigonometric functions in (3) are the eigenfunctions of the Sturm-Liouville problem encountered in stability analysis of the Lapwood system [6]. The quadratic form in (4) may be compared to the sums presented in [7] developed for the DLB (Darcy-Lapwood-Brinkman) system. They include the same amplitude combination in the quadratic form and the later introduced flak amplitudes may be found in [7, Equation (86)], even though their system is very different.
3.2. Linear Stability Analysis
In linear stability analysis we assume no flow to be present; this means linear temperature gradient or that all amplitudes in (3) are zero. Using the regular perturbation theory, we then assume a perturbation introduced in the form of infinitesimal amplitudes. The quadratic form in (4) will now drop out and we are left with the exponential growth terms, that is, (4) without the quadratic form. Stability requires all growth rates to be negative or zero: is the critical Rayleigh number for the mode (the th eigenfunction). The classical value . If and is a whole number defined by modes are unstable and all unstable wave numbers belonging to these modes will grow exponentially in a flow started from rest. The unstable wave numbers are within the wave number instability band, between and in the following equation, they are the zeros of (5) Maximum growth rate wave number, sometimes called the most critical wave number, is For (ca. 40) we find and it does not become until so, in a Fourier transform with , the will be the fastest growing wave number in flows between the first and second critical Rayleigh numbers. This sets the width of the flow cell.
3.3. Nonlinear Stability Analysis
The elements in the coefficient tensor in (4) can be calculated from the following algorithm: ; otherwise .
We now assume that a stationary solution to the system equations (4)–(10) exists. It has to fulfil the system equations with the term . No matter what the value of individual amplitudes is, this system must be stable against all perturbations also if a perturbation is placed on only one amplitude but the others are unperturbed. This leads to the following equation for the perturbation of the stable amplitude : This is the same equation as the linear counterpart, only with a different growth rate, which from (10) is found to be This is the non-linear growth rate instead of the linear one . The only difference between the nonlinear and the linear growth rate is that we must insert in (5) a new :
As in the linear case stability requires where the = sign means neutral stability. It must be noted that this is a necessary condition for any solution to be stable, but it is not sufficient as all amplitudes are not perturbed.
If a stable solution turns unstable (e.g., from increasing ) the disturbance of one amplitude may eventually spread to all the others. Taking themselves they have this equation Here there are two opposing forces, the strong curbing effect of the negative linear growth rate and the negative sum that carries in it the total spectral energy on the th mode. There is infinity of other quadratic terms, not shown, but these have a good chance of cancelling each other out in a time average, as they would contain the correlation coefficient between the respective amplitudes. If any amplitude on any mode gets very big the will grow strongly negative with it and then pull it down as the nonlinear growth rate turns negative. Such action may clearly be seen in the following simulations and therefore we call the ’s flak amplitudes. The amplitude combinations in the quadratic form of (4) reflect the fundamental property of the trigonometric functions that . They are the same for all two-dimensional systems having the nonlinearities in quadratic terms only.
The flak amplitudes govern the flow pattern. They grow and diminish with the energy on the corresponding unstable mode by controlling the nonlinear growth rates of all flow amplitudes. Later we see that they make up the horizontal average temperature distribution and with it the Nu. They are always negative and must not be above a certain stability value in any stationary flow that might exist. There is a single infinity of critical Rayleigh numbers; each time one is passed while the fluid is heating up, the corresponding flak amplitude has to grow to a significant value.
It may be shown that (1) and (2) have a symmetrical solution that is even; that is, all amplitudes in (4) where is an uneven number are zero. In [4, 7] this is used; the consequence is that, when only one stable solution exists for each , it is even. The physical difference of uneven and even solutions is that uneven solutions have uneven flak amplitudes with being an uneven number. Then the average heat flow at top and bottom is not the same and the porous layer may be heating up or cooling down as a whole. An even Fourier transform is more stable than the one where uneven amplitudes are allowed. But it seems inevitable that uneven amplitudes can participate in unsteady flow. Stability analysis of the even solution [6] indicate that this solution may be stable in the number range or between the two first critical numbers.
Equation (14) does not have to hold for fluctuating motion when many modes are unstable. It is however difficult to see how fluctuating motion (Figure 4) can be maintained unless the flak amplitudes do fluctuate around (14) value.
3.4. The Spectral Method Using Discrete Fourier Transform
In numerical calculations (4) cannot be truncated without introducing systematic errors. Finite Fourier transforms will be used instead and then new amplitude combinations appear in the quadratic form in (4). The new quadratic form becomes where and are the horizontal and vertical maximum wave numbers; is Kronecker’s delta; and are sign control functions that are to be taken from Table 1 as plus or minus one.
As may be seen by comparing (10) and (17), truncation of (3) may introduce a significant error. The discrete quadratic equation (16) is evaluated for discrete values of and only, so it contains amplitude combinations that are not at all present in the infinite quadratic equations (4) and (10), as explained in [6]. If the truncation is done as a cutoff at constant these amplitude combinations are totally absent.
When an infinite transform is truncated, all the amplitudes above the cut-off frequency are considered to be zero, but discrete forms may be cut off anywhere, discrete forms may be cut off anywhere. A discrete Fourier transform is a process where any number of points is transformed into an equal number of Fourier coefficients, so the transform creates a function that goes through the entire original points exactly. Therefore, for all situations, stable or unstable, there exists a finite Fourier transform for any number of and . It may be noted that FFT algorithm is a discrete transformation so turbulence simulations using this technique have the extra frequency combinations included.
4. Simulations Using the Discrete Fourier Transform
4.1. General Remarks
To simulate, as well as and has to be selected. When fluctuating motion like this is simulated using (4) with a constant value of , one is actually simulating convection in a box. Simulations using (16) and (17) may be fired up by putting a small constant on all amplitudes (spatially distributed random noise). For low Rayleigh numbers this leads to a stable solution where all the uneven amplitudes disappear. In order to give a realistic picture of the flow, simulations need to include sufficiently many modes, so all unstable modes are controlled. must therefore be not less than two times the number of unstable modes, at least in theory. Intuitively, one would expect too low simulations to be unstable for high Rayleigh numbers, but on the contrary they are more stable. The physical explanation of this is that using a discrete transform with low means that we only have many values for the and and each of these values must represent the average in the corresponding rectangle. But averaging the equations means that a new coefficient of heat conduction, similar to the eddy viscosity, appears on the scene (21). This is explained in the section on scaling.
In the very weak instability phase ( little higher than 40), only one mode is unstable and gives the maximum growth rate of the instabilities, so in (14) the period will grow fastest when we use . This is therefore used in all the simulations; it makes the results better comparable.
4.2. Simulations with NI × NJ Grid
Possible number of amplitudes is 6, , making 3 even amplitudes but one of them, , turns out to be zero. The two remaining are the and its flak amplitude .
This approximation has been studied by [4], who use a series expansion to find it, and [6] that uses (16) and (17). This rather crude approximation can take all ’s and is always stable.
The formula for this solution is Figure 1 shows the simulation for Rayleigh numbers 50–250. All simulations with produce the same stable result as (18) and the flak is the same as (14) with the equal sign. The maximum value of that (18) can give is 3; this low value of shows better than the value of the amplitudes how crude the approximation (18) is. There are in reality 5 unstable modes, but only two are in the simulation and of their amplitudes only and are nonzero. Increasing to 4 but keeping changes very little. More unstable modes are needed, Figure 2. But a stable solution is obtained for , that is, 5 unstable modes.
4.3. Simulations with NI × NJ Grid
Figure 3 shows and simulation. This grid is slightly more accurate than the . It contains 2 flak amplitudes and should therefore be able to give a realistic picture of stable stationary flows that might exist up to the third critical Rayleigh number. The possible number of amplitudes is 32 making 16 even amplitudes. Of them 9 end up with significant values to the second digit in Table 2.
In Table 2 uneven amplitudes are in lighter shade. Runs with uneven amplitudes in the spectrum show that all uneven amplitudes diminish with time and drop out. Simulated is −0.31 and ; this includes some aliasing effects from the finite transform. Equations (14) and (21) give that is −0.25 and . The results compare favourably and show how the flak amplitudes turn the unstable modes into stationary motion. Figure 3 shows the time history of the simulation.
The slight increase in A from 180 to 190 changes this completely. Now we have a stable periodically fluctuating solution, Figure 4. The grid approximation does not render stable solutions for higher Rayleigh numbers.
4.4. Simulations with NI × NJ Grid
Here we look for how the solutions in the grid look in a slightly more accurate grid that contains 3 flak amplitudes and examine the fluctuating motion more closely. The possible number of amplitudes is 60 making 30 even amplitudes.
Simulations in a grid are a little more unstable than in the grid. Fluctuating periodic flow is reached at and can be maintained until ; see Figure 5. The average amplitude array of is shown in Table 3, and the much larger standard deviations are in Table 4. Here we have the interesting result that the average flak amplitude is while all others are practically zero; this gives . Equation (14) stability limit for is −0.2411 below the simulation value. This is to be expected as this solution is not stationary. Noting the standard deviation in Table 4 (0.08), it corresponds to a fluctuation of about 0.1 around the mean value. The value of this flak amplitude will therefore be between ca. −0.1 and −0.3, with these values being on the unstable and stable side, respectively, of the stability limit (14). Figures 4 and 5 show this flak amplitude fluctuation clearly. Equation (21) gives and we see that flak amplitudes and are diminished from the stationary level but still up and active.
The simulations in the and grids show that stable solutions can be maintained only up to the second critical where periodic motion begins. Similar behaviour is observed in fluid turbulence; when the Reynolds number is increased, fluctuating motion sets in. Many vortex flow fields of this type are very well known in fluid mechanics (Reynolds, Taylor, Karman, and Kelvin-Helmholtz). With increasing instability these regular flows disappear and chaotic turbulence appears instead. This seems to happen here between the second and the third critical Reynolds numbers. This is in accordance with the findings in [4, 6]; stable motion is concluded to cease at the second critical Rayleigh number. Chaotic motion is found in [7] for that system at high fluid numbers.
5. Scaling
The amplitudes in Tables 2 and 3 have average amplitudes and high standard deviation. The effect of this is to create a flow with “eddy heat conduction,” a phenomenon similar to eddy momentum transport in fluid turbulence.
It has been shown [17] that, when a running average is taken in the convection box (Figure 6) covering a grid area of , the nonlinear terms in (1) will produce a net transport of heat into the element when the second derivative of the average temperature distribution does not vanish. This net heat flow may be represented as the divergence of an eddy heat flow vector equal to grad. This is the heat transported by fluctuations in excess of the molecular heat conduction and the convection by means of the average flow velocity. To model the flow in a coarse grid (Figure 6(a)) it is therefore necessary to include in the simulations a subgrid model that represents this heat flow, just as it is necessary to include a subgrid model to take care of the Reynolds stress tensor in macroscale models of turbulent flow.
(a)
(b)
The widely popular grid solution is one flow cell in one single block, so the underlying assumption is = constant. This is a crude simplification, similar to Prandtl’s mixing length theory. It is not quite correct, as would most likely take the highest values in the zone shown in Figure 6(b). We must conclude that in this solution there is active eddy heat conduction coefficient acting in excess of the normal heat conduction coefficient so the which we must use in (16) is scaled down in the following manner: In finding we have a closure problem just as in turbulence; principally it must be found in experiments. From the simulations we learn that in stable and oscillatory solutions the flak amplitudes always have significant values. By averaging the temperature gradient given by (3) at either top or bottom and calculating the heat flow and using the equal sign in (14) for the flak we get This result is an approximation valid for numbers larger than the second critical Rayleigh number, If is used instead of then the grid solution (18) approximates the flow field in this region. One has to remember that the has to be used instead of when calculating the heat flow from in (18); otherwise this approximation will render too low heat flow. The numbers can be estimated directly by (21); for high , is a good approximation. It should be noted that (21) means that stationary average flow can be approximated up to .
Bifurcations are known to occur in porous media [17]. The first bifurcation happens in DL systems at the first critical number. According to (14), the physical process in fluctuating flow is that the flak amplitudes maintain significant values to keep the nonlinear growth rates down and hinder the amplitudes they control from increasing without limit. Then there must be a new bifurcation each time passes a critical number. Equation (21) is the result of repeated bifurcation each time the Rayleigh number passes a critical value and one more mode is rendered unstable in the free convection. The result may be seen in Figure 7.
6. Laboratory Scales
To make laboratory tests of convection in porous media one has to scale natural convection (e.g., geothermal fields) down to laboratory dimensions; in doing this, dynamic similarity may be a problem. One dimensionless parameter, for example, the Rayleigh number, can be kept constant in the model and the prototype but not more than that. Aquifer thickness, , is typically 1000 m in geothermal reservoirs but ~1 m in the laboratory. To get a flow going the aquifer rock matrix has to be replaced by glass or plastic pearls with up to 1000-fold permeability and this brings the test from the DL regime of the aquifer into the DLFB regime. Acknowledging this fact [7] consequently uses , not . This brings up the question of dispersion. This seems to be overlooked, but dispersion must play a role in laboratory experiments with a porous matrix of relatively large glass pearls. The scale effects in the dispersion term may be estimated from the dispersion equation: This assumes the axis in the direction of the flow. is the dispersion tensor. denotes the concentration of a solute, in our case the mixing ration of hot water into cold, or the dimensionless temperature . The term is the molecular diffusion, in our case the coefficient of heat conduction as before. In natural low permeability aquifers, is very small and the dispersion term is negligible compared to the term and generally left out. Estimating the dispersion/conduction ratio DCR, with respect to dynamic similitude, we see that the DCR will not be the same dimensionless number in model and prototype regardless if DL scaling ( constant), or DLFB scaling (Ra constant) is used, the following equation shows this, Here is used for the velocity scale in the dispersion term, is the porosity parameter [7], is called the Darcy number by some, and Ra is the fluid Rayleigh number . The parameter a is called the dispersivity and is usually taken as proportional to length scale; [20] suggests the order of magnitude 0.01–0.1. To use this value, both the aquifers in nature and the laboratory would have to be glass pearls which of course is not the case, but if that is overlooked, scaling ( is the same in model and prototype) would give dynamic similarity. When scaling is applied, the ratio can easily be in the range 103–106; this makes the values much higher in the model than in the prototype as . Then heat dispersion is more important than heat conduction resulting in higher Nusselt’s number in the laboratory tests than in natural aquifers of large dimensions. DLFB systems can accordingly not be scaled to DL systems when the FB terms have a significant effect.
6.1. Nusselt Number Compared To Experiments
Trying to rescale laboratory results and compare computed and measured Nusselt’s numbers is impossible, except possibly for the two-amplitude (one block) solution, where dispersion effects may be included in the . The grid procedure can only produce numbers up to 3 in the Lapwood system. [7] bring it up to approximately 4 in their system using a cutoff frequency of , which corresponds to the grid solution. Nusselt’s numbers in natural aquifers of large dimensions tend to be higher.
Figure 7 shows Nusselt’s number based on the assumption that repeated bifurcation brings into the picture new flak amplitudes for each new mode. The formula (blue line, Figure 7) is a good approximation for . values from laboratory tests are higher, but effects of eventual scaling or forcing (see next chapter) are unknown so the triangle for the results location is only approximate for pure DL systems.
6.2. Effect of Forcing
Due to (2), all horizontal temperature gradients from outside force the DL system. Forcing creates some flow, no matter how low the number is. Forced convection is an interesting topic with many applications; in [21] there is a review of the results in this field. When heating and cooling are between two vertical walls, stability does not have to be a problem and elegant solutions can be found both analytically and numerically [22]. The work in [23] is very interesting, showing an analytical solution of a forced system with variable permeability.
In [24] there is a treatment of the heat transfer through a box with heated sides using the DLFB system with added turbulence by a subgrid model and a dispersion term included. The turbulence increases significantly but the effect of changing from 107 to 108 is even greater. They use but if their Figure 4 in [24] is rescaled to , the lower σ2 curve (higher Darcy number, Da in their notation) produces higher numbers for the same , which illustrates the scaling difficulty mentioned above and could point to the influence of dispersion. [25] studies the same problem and the conclusion is that dispersion has a significant effect.
Forcing thus makes it possible to analyze numerically and analytically situations with high heat transfers counted in numbers, situations where numerical simulations of free convection with similar numbers become unstable because the flow fluctuates. The effect of basic processes such as dispersion must generally have the same effect on the magnitude of heat flow in forced and free convection problems.
7. Discussion
The physical process of free convection flow in porous media suggests itself to be the following. For Rayleigh convection sets in as a result of a bifurcation process when small disturbances of the spatial wave number become unstable and start to grow exponentially with time. For wave numbers below the second critical wave number only one mode (the first) is unstable. The entire energy spectrum of unstable wave numbers in that mode participates in making the first flak amplitude highly negative. Then the nonlinear growth rates of all first mode wave numbers turn negative and all the flow amplitudes on the first mode start decaying. This process repeats itself until stationary cellular flow of convective rolls is achieved (Figure 2).
If the Rayleigh number is increased slightly above the second critical Rayleigh number, a fluctuating flow sets in, but the grid approximation does not show fluctuations at any Rayleigh number. In contrast hereto, the and grid discrete Fourier transforms (16) and (17) show stable periodic motion for numbers just above the , but it changes to fluctuating motion that resembles eddy motion in turbulent flow when the is slightly increased. In this fluctuating flow, the flow amplitudes shift from positive to negative indicating reversal of the flow direction in a periodic manner, but the flak amplitude stays high and fluctuates around the stability value (14). In this flow the average amplitudes are less important than the standard deviations and the correlation between amplitudes. The flak amplitudes are the only exception.
For higher Rayleigh numbers the flow becomes chaotic. Forcing will delay this considerably, so in laboratory experiments it is very important to exclude all forcing.
This flow resembles fluid turbulence in many ways. In spatially averaged equations the average amplitudes of fluctuations smaller than the grid size drop out leaving an eddy heat conduction effect. This leads to the scaling rule that makes it possible to use the grid solution and still get realistic numbers.
Turbulent fluid flow is governed by the eddy momentum transport due to the fluctuation that results from quadratic forms of the velocities in the Reynolds stress tensor. In a Darcy-Lapwood system there is eddy heat flow due to the fluctuation that results from quadratic forms of the velocities and the temperature. The flow is very slow, the nondimensional time interval from 0 to 1 can mean 30.000 years, [6] and heat conduction is faster than heat dispersion.
Systems governed by the fluid Rayleigh number do not in principle scale to Lapwood systems, and there are strong indications that systems (DLFB, DLB, and DL systems with added effects of turbulence or dispersion) scaled down to laboratory model size with run on higher heat flow due to dispersion in the pore matrix than the prototype.
As fluctuations dominate the flow at high number flows and contribute to the heat flow, stable nonfluctuating solutions, analytical as the approximation or numerical, do need something like the eddy heat flow coefficient to render correct heat flow when the convection is free.
This investigation covers convective rolls. When strong fluctuations set in, the rolls become unstable and the flow becomes three-dimensional. The rolls are still there as a background motion, but when passes 3000–4000 they have probably disappeared. In [26] a scaling approximation with a boundary layer at top and bottom ~ is used for ; it gives ~, but there is fair agreement with the results here, up to that point. To investigate three-dimensional flow, it would be interesting to perform an analysis with the finite Fourier transform of a two- or three-roll system intersecting each other. But this would still produce a quadratic form similar to what we have, so it does not have to change much. Here it is judged unlikely that it would suddenly change the ~ to ~, ([26], Figure 2).
8. Conclusions
In the Darcy-Lapwood system there are unstable modes ( nearest lower integer).
Associated with all amplitudes, on a fixed mode (constant ), there is a flak amplitude. It is the amplitude on the zero frequency of the mode with the double mode number being the same for all horizontal wave numbers . The flak amplitude is always negative, independent of the sign of the amplitudes.
The flak amplitudes have a neutral stability value as shown by (14). If a stationary solution exists, it is symmetrical with all flak amplitudes not higher (less negative) than this value.
The Fourier transforms of the spectral method do have to include active modes up to double the number of unstable modes. Using discrete Fourier transform includes several frequency combinations generated by the nonlinear terms of the system (included in the Jacobian of (1)) that the truncated infinite form does not have.
Using the discrete Fourier transform makes it possible to use very coarse grids in numerical simulations, such as the approximation, but then the flows on unstable modes not represented in the grid are aliased on the existing modes. For high coarse grids give too low Nusselt’s numbers even though the fluid flow picture seems realistic.
Simulations above the second critical result in fluctuating flow. Even though the average fluctuations may have zero average amplitude and thus do not participate in the average fluid motion, they are active in transporting heat.
The assumption that all the flak amplitudes fluctuate around the stability value results in the scaling rule (17) that makes it possible to use the solution to obtain realistic heat flow by defining an eddy heat conduction coefficient that makes up for the missing effect of the fluctuations not present in the approximate solution.
Assuming the flak amplitudes to stay on the neutral stability value on the average results in the approximate formula for . Then the approximation renders realistic , when is used. Without this scaling, is maximum in the approximation. It gives a stable flow for all .
With 2 unstable modes fluctuating motion sets in with the fluctuations more dominating in the finer grid transforms, as finer grids bring higher growth rates into the simulation.
Today’s research focuses on DLFB systems as this equation system is more accurate for high permeability small scale flow systems. Studies of these systems do not automatically include the mathematically simpler DL system because of scaling problems associated with the great difference in the porosity value that follows the higher permeability and smaller linear scale.
Strong forcing stabilizes flows that would be otherwise unstable. But flow on the first mode does not stabilize the flow in DL systems as flak amplitudes have no effect outside their own mode. A DL system will pick up disturbances from outside that fall into the wave number instability band window, when the corresponding flak is in a downswing. Natural disturbances from outside, for example, tidal and barometric effects, will therefore start new fluctuations and the system will by time become independent of the initial condition thus forgetting its past and becoming chaotic.
Conflict of Interests
The author declares that there is no conflict of interests regarding the publication of this paper.