Spectral properties of exact random solutions to Burgers' equation for modified Thomas initial conditions

Computers ,4t Fluids Vol. 16, No. 2, pp. 147-173, 1988 Printed in Great Britain. All rights r--~-ecrved

0045-7930/88 $3.00+0.00 Copyright © 1988 Pergamon Journals Ltd

SPECTRAL PROPERTIES OF EXACT RANDOM SOLUTIONS TO BURGERS' EQUATION FOR MODIFIED THOMAS INITIAL CONDITIONSt STEVEN KELETI a n d X B REED JR Chemical Engineering Department, University of Missouri-Rolla, Rolla, MO 65401, U.S.A.

(Received 6 January 1987; in revised form 1 July 1987) Abstract--Random sawtooth Thomas [1] initial conditions can be specified on a fixed spatial mesh d or on a random one d~ (with average mesh spacing (d~)) termed modified Thomas initial conditions. Beyond the white-noise band common to both at small wave number k, the energy spectrum E(k, 0) for the modified Thomas has a smooth transition from a k -4 to a k -2 asymptote, whereas Thomas has pronounced ringing. For modified Thomas conditions, there is (i) a narrow hand (small k) white-noise spectrum E(k, t)= E o like that for Thomas conditions, (ii) a k -2 subrange of increasing bandwidth with increasing initial turbulence Reynolds' number Reo, and (iii) an exponential viscous cutoff at high wave numbers which decays more rapidly with decreasing Re0. The small time spectral transfer, T(k,t), for modified Thomas conditions does not have a spike at k -~ 2n/(d~), whereas it does for Thomas ones at k = 2~/d. The evolution of the energy and dissipation spectra and transfer and cumulative transfer spectra for modified Thomas conditions is otherwise similar in form to that for Thomas conditions for t > 1.

1. I N T R O D U C T I O N

Theoretical turbulence may be described as standing in need of the next generation of supercomputers or as suffering from the lack of random solutions to the Navier-Stokes equations. Burgers' equation, a I-D analog which is often used to test numerical methods and theoretical closures, can be solved exactly in terms of standard functions of mathematical physics by means of the Cole--Hopf [2, 3] transformation, provided: (i) The initial conditions for Burgers' equation maps under the inverse Cole-Hopf transformation into mathematically tractable initial conditions for the heat conduction equation (i.e. so that more than a purely formal solution is obtained). (ii) The boundary conditions for Burgers' equation are either to be imposed as, say, asymptotically vanishing or constant on infinite domains [4] or have counterparts which make for tractable boundary value problems for the heat conduction equation (e.g. Refs [5, 6]). Burgers devoted most of his life to the analysis of the equation, c~u/at + ucgu/c~x = va2u/cgx 2, which bears his name (see especially Refs [7-10]), including its statistical properties. Many closures have been tried on Burgers' equation prior to the use of their 3-D generalizations (Refs [I 1-13] for example). In his continuing development of the Wiener-Hermite expansion, Meecham, with a number of collaborators, has obtained significant results on Burgers' equation per se (e.g. Refs [14, 15]). In particular, he and several co-workers [16, 17] were the first to obtain an exact, (almost) statistically homogeneous, solution to Burgers' equation on an infinite domain for nonzero random piecewise constant initial conditions on a large but finite domain (evaluated at Re0 -- 200). Shih and Reed [18] observed that the Laplace solution to the heat conduction equation could also be used--certainly with more calculational difficulty--to obtain a solution to Burgers' equation (also evaluated at Reo = 200) on an infinite domain for random piecewise linear continuous initial conditions which are nonzero on a finite domain. These physically more realistic initial conditions were also used in the study of exact solutions tThis paper was presented at the 10th Symposium on Turbulence held at University of Missouri-Rolla, 22-24 September 1986. 147

148

STEVEN KELETI and X B REED JR

to Burgers' equation on a finite domain with spatially periodic [19] and with vanishing [20] boundary conditions at Re0 = 200. The above investigations provided detailed behavior of all but the fine scale structure for (u(x, t ) u ( x + r, t)) and its Fourier transform E(k, t), with ( ) denoting ensemble averaging. In Refs [18-20], histograms and the first nine (single-point) moments were obtained, with (u2(x, t ) u ( x + r, t)) also being computed in Ref. [19]. The ensembles in Refs [18-20] were large by the standards of Jeng et al. [17], consisting of 12 realizations on the IBM equipment, but small by the standards necessary to accurately compute estimators of two-point statistical properties. Moreover, the spatial mesh on which the exact solutions were numerically evaluated [18-20] was inadequate to obtain spatial correlations for extremely small separations and equivalently high wave-number spectral results. For large ensembles (50, occasionally even 100 members), Keleti and Reed [21] obtained preliminary energy, dissipation, and transfer spectra at Re0 = 400 and 1000 on an FPS-164. These and subsequent results also extended well into the high wave-number range through the use of a rapid, accurate algorithm for the evaluation of exact solutions to Burgers' equation on a very fine spatial mesh [21-23]. The initial conditions were those first employed by Thomas [1] and subsequently used by Shih and Reed [18-20]: values of velocity {u,} selected at random from a uniform ( - 1 , 1) distribution on a fixed spatial mesh {x~- xn_l} = d (d = 1) were connected linearly. The random piecewise linear continuous initial conditions generated in this way are called Thomas initial conditions. Piecewise linear continuous random initial conditions for Burgers' equation have several advantages, among which are features analogous to velocity records for homogeneous turbulent velocity fluctuations [18]. Thomas [1] initial conditions have a distribution of velocity scales but only a single length scale. This manifests itself in several ways, for instance in a spectral peak at k = 2n for small times [22]. The question that therefore arises is whether or not the statistical properties of velocity fields evolving from Thomas initial conditions are representative of a larger class of exact random solutions to Burgers' equation, this being a basic tenet of homogeneous turbulence (e.g. Ref. [24]) and more generally of statistical mechanics [25]. (For a brief discussion and a demonstration of a stronger result for a specific pair of classes of initial conditions, see Ref. [26].) An obvious way to test the hypothesis is to consider initial conditions with a distribution of length scales, as well as a distribution of velocity scales. The simplest version of these initial conditions, which we have termed modified Thomas initial conditions, may be generated by selecting the {dn) = { x , - Xo_l} from a uniform (0, 2) distribution, as well as selecting the {u,} from a uniform ( - 1, 1) distribution, with the two random variates being statistically independent and therefore uncorrelated. 2. I N I T I A L C O N D I T I O N S , V E L O C I T Y F I E L D S A N D S O M E SINGLE-POINT STATISTICAL PROPERTIES

The initial conditions are specified by the random location of {xn} at which the random velocity {un} is specified. The {un} are selected from a uniform distribution on ( - 1, l). The {x~} are selected from a uniform distribution on (0, 2). In both cases, basically, the IMSL (pseudo)random number generator algorithm was used [27, 28], except that intermediate (pseudo)random numbers ~, were used as seeds to generate the actual random numbers ~k instead of the recommended seeds, ~k-~: ~ = MOD (16807.0~k_~, 232 -- 1), ~k = MOD (16807.0~, 232 -- l). Normalization to (0, 1) was by division by 232 and the first seed was ~0 = 123457. The un on ( - l, l) were then obtained as u, = (2~n/232 - l). Similarly, dn = x. - x._ ~ on (0, 2) was obtained as d~ = 2~n/2 32 from different ~. There are pitfalls in the generation of (pseudo)random numbers [27], but the above use of a random seed has been investigated by, for instance, MacElroy and Suh [28] and shown to avoid most difficulties. The number

Spectral properties of exact random solutions to Burgers' equation

149

of (pseudo)random numbers {d.} and {u.} was significantly less than the cycling that can arise, assuring negligible correlation between the (pseudo)random numbers [27, 28]. The modified Thomas initial conditions are selected with the uniformly (0, 2)-distributed {d,} being renormalized such that over each ensemble member ( d n ) = d = 1. Each member is selected on a periodic domain of length 2L = 200d, the domain size selected for the Thomas initial conditions [18, 19, 21, 22]. Using a solution over a periodic domain allows an investigation of the spectral properties of Burgers' equation without addressing issues concerning the use of the Fourier transform over a finite domain. The differences between the velocity profiles for Thomas and modified Thomas initial conditions are shown in Fig. 1 for short times. The values of the {Un} are identical in Figs l(a) and (b), being selected at random from a uniform ( - 1, 1) distribution, as explained above. The fixed spatial mesh {Xn -- Xn-l} = d = 1 for t = 0 is apparent in Fig. l(a), as is the random spatial mesh {Xn - Xn-I} = {d.}, ( d . ) = 1 in Fig. l(b) with ix. - x._,} selected as explained above. Both Figs 1(a) and (b) show the temporal evolution of velocity profiles at Re0 = (d.)uo/v = 400 with u0 = sup[u.l = 1. The effect of the distribution of spatial scales on the velocity field and thence on the spectral properties for modified Thomas initial conditions were the object of the present study. The modified Thomas initial conditions have the same initial turbulence energy as do the Thomas ones, a consequence of the {u.} being selected from the same distribution and independently of the {dn} (see Appendix A). The distribution of that energy is quite different, however, and one would therefore expect an effect of the initial conditions on the turbulence energy decay. Such is the case. Three inviscid ingredients provide a basis for the explanation (see, e.g. Ref. [29]): (i) corners u, move at a constant speed until they form or are swept up by a shock, (ii) segments with negative slopes form shocks, and (iii) ]u,(t)[ decays for shocks imbedded between positive slopes. The distribution of initial length scales of modified Thomas conditioned insures that more shocks will form much sooner than for Thomas initial conditions. Once a shock is formed, its energy decays and thus the energy of modified Thomas initial conditions initially decays faster than that of Thomas initial conditions. By t ~ 1, in contrast, most negative slopes for Thomas initial conditions have become shocks and an accelerated energy decay is observed, relative to modified Thomas initial conditions (Fig. 2). The last shock can form much later for modified Thomas initial conditions, thus delaying the approach of modified Thomas initial conditions to the asymptotic decay rate which it has in common with Thomas initial conditions. There is no appreciable effect of Reynolds' number on the energy decay rate because Re0 is sufficiently high [24]. Both initial conditions have an energy decay rate which asymptotically approaches t -2/3, consistent with the original theory of Burgers [9] and with other later research (e.g. Refs [30-32]). Other basic single-point statistical properties are tabulated in Table 1 for comparison with those for Thomas initial conditions [22]. The rms velocity u'(t) decays monotonically, and the effect of Re0 is in the third or fourth decimal place. The variance of Ux increases initially and then decreases and, in contrast to u', depends markedly on Re0, with (Ux)ma x2 being much higher for Re0 = 1000 than for 400. The Taylor microscale 2(t) may be determined from energy decay rates. The values in Table 1 were obtained from j2(t ) _ (u2) _

o

E(k, t) dk k2E(k, t) dk

o

estimated by ½N-,

/½N-I

~, E(k = 2nf /NAx, t)/ /~=ok2E(k = 27rf/NAx, t).

f=0 CAF

16,2~"

150

ST~V~N K~L~Taand X B R~D JR

o

0

ro o

II "~

.~

%%

0,,-~

~8 o

....... x4,

.r" o

b~

. . . . . .~ . . . . . ~. ............... 0 ~

~ 0

d

d

I

o

6

~

o

T

T

d

d

,~!OO'leA

o.,

oI

T

Spectral properties of exact random solutions to Burgers' equation

-0.6 -0.7

151

Energy decoy ~

\

~

.

.

-0.8 -0.9

~\

_ bd - 1 . 0

Modified Thomas ~ ( R e o " 400,1000)

Thomas NNN\

o

(Re(:} "400,1000)

~N\ \

O _1

X\\~t-2/3

-1.2 -1.3 -1.4

-1.5

i

i

t

I i , i i I I

i

I

I I I I I [ 10

[

I

Ti me Fig. 2. Comparison of energy decay curves for Thomas initial conditions and for modified Thomas

initial conditions.

With 2 initially decreasing and then increasing, consonant with the behavior of (u2), Re~ decreases sharply at first--the decrease being the more precipitous, the higher the Re0--and then gradually increases. An increase in Re0 by a factor of 2.5 produces an asymptotic increase in Re~ by a factor of 1.6. Burgers (eqn (94) of Ref. [9]) derived 2 2 = 3Vt,

which together with the (u 2) asymptote of t -2/3 gives Re~ = ( u 2 ) U 2 2 / v ~ 3u2tl/6/v 1/2, from which follows Rea (Re0 = 1000) = (1000/400)1/2 = 1.6. Re~ (Re0 = 400) Although the evolution of Re~ is available in the literature, for instance, from modified 0-4th cumulant expansions [33] (t °'°~°'~3) and from the finite difference computations of Mizushima and Saito [34] (t°J°6), the results presented here (t °'~°) are the first computations accurate enough to confirm the Re~ ~ t ~/6/v ~/2 asymptote for Burgers' equation. Burgers (eqn (92b) of Ref. [9]) also found ( u ~ ) ,,. t 5/3. These asymptotes may be compared with those calculated from the last four points and presented in Table 2. Typical values of Re~ for wind tunnel and locally isotropic turbulence for comparison with values in Table 1 are

Table 1. Modified Thomas statistics u' =

t 0 0.1 0.2 0.3 0.5 I 2 3 5 10 15 20 25

(u~) '/2


~. = (
Re~ = u'~.lv

Reo=400 Reo=1000 Reo=400 Reo=1000 Reo=400 Reo:1000 Reo:400 0.46904 0.46818 0.46715 0.46587 0.46258 0.45020 0.41326 0.38009 0.33293 0.27110 0.23861 0.21774 0.20244

0.46904 0.46856 0.46780 0.46674 0.46376 0.45180 0.41469 0.38129 0.33382 0.27166 0.23903 0.21810 0.20274

2.36421 1.71851 2.15868 2.60907 3.48001 5.69081 6.07744 4.46540 2.52408 0.92464 0.47594 0.30608 0.21366

2.36421 2.85349 4.28623 5.61764 8.07615 14.14681 15.34581 11.25774 6.37163 2.32827 1.19714 0.76941 0.53736

0.30505 0.35714 0.31795 0.28842 0.24797 0.18872 0.16763 0.17987 0.20956 0.28193 0.34586 0.39357 0.43796

0.30505 0.27738 0.22595 0.19692 0.16319 0.12012 0.10586 0.11364 0.13225 0.17804 0.21846 0.24864 0.27658

57.23219 66.88275 59.41207 53.74629 45.88266 33.98509 27.71062 27.34645 27.90740 30.57163 33.00973 34.27918 35.46457

Reo=1000 143.08047 129.97195 105.70135 91.91029 75.68051 54.27158 43.89930 43.32853 44.14654 48.36505 52.21941 54.22796 56.07402

STEVENKELETIand X B Rr~DJR

152

Table 2. Asymptoticevolutionof modifiedThomasstatistics Reo=,)oo Reo=looo ( u 2) (u~) 2 Re~

0.3187 t o,6367 36.1461 t 1.5945 0.09391 t +0,4789 21.2070 t +0 ~606

0.3210 t -°'6379 91.3100 t.1.5960 0.05929 t +0.4790 33.5920 t +0~6o~

50 by Van Atta and Chen [35], 70 by Uberoi [36] and 500 by Kistler and Vrebalovich [37], compared to about 2000 by Grant, Stewart and Moilliet [38]. The above statistical measures may be compared with those for Thomas initial conditions [22] as follows (see also Fig. 2). The variance (Ux2) for Thomas initial conditions [22] is initially much smaller but grows to a much higher maximum and then decays to lower values than do those for the modified Thomas. The initial decay of 2(Q to a minimum is more rapid for Thomas, but then its gradual growth recovers to that for modified Thomas. Consequently, Rex for Thomas initial conditions begins roughly a factor of 2 higher but ends (at t = 25) slightly lower than Re~ for modified Thomas.

3.

COMPUTATIONAL

METHODS

Here and in Ref. [22] u(x, t) has been computed with high accuracy (greater than 10 decimal places) on a fine spatial mesh (200/216) for a quite large ensemble (50 members, which gave results that were imperceptibly distinct from preliminary studies with 100 members), allowing spectral resolution far into the viscous cutoff. Numerical evaluation of exact solutions to Burgers equation are often avoided due to their complexity [39]. In addition, exact solutions can have magnitudes of intermediate computations which vary over an extremely wide range, a range which increases exponentially with decreasing viscosity v. Because the range of magnitudes grows beyond the range of most computers even at small Reynolds' numbers and because the arguments of some terms can become imaginary, the integrals in exact solutions are usually computed numerically instead of in closed form (see, for example Refs [16-18, 39]). Here and in Refs [21-23], in contrast, a closed-form Laplace solution has been used. As an algorithm for efficiently handling the problems associated with this solution is not available in the open literature, it is presented in Appendix B. The consistent discrete Fourier transform pair which was used in the present computations was: N-I

a(k = 2rcf /NAx, t) = Ax ~ u(x = nAx, t ) e x p [ - i 2 n f n / N ] , n=0

f = { 0 , 1. . . . . N - l } , l N-l U(X = nAx, t) = N A x r~=o~(k = 2nf /NAx, t)exp[i2nfn/N], n = { 0 , 1. . . . . N - l } . Computations were implemented on the FPS-164 using FFT techniques [40] such as have been used in turbulence by Van Atta and Chen [41]. Energy spectra were computed from:

E(k = 2nf / N Ax, t) = :~*(k, t )f~(k, t ) /N Ax, which is the Fourier transform of the autocorrelation function l

N-I

Q(r = nAx, t) - N A x f=~o E(k = 2 n f / N A x , t)exp[i2nfn/N], n={O, 1. . . . . N - l } ,

Q(r, t) = (u(x, t ) u ( x + r, t)).

Spectral properties of exact random solutions to Burgers' equation

153

The dissipation spectrum D(k, t) can be computed by first computing u~(x, t) and then employing the Fourier transform (in this case, FFT) or by recognizing that D(k, t ) = k2E(k, t). The latter is obviously simpler and computationally less intensive. That it is also more accurate at a given spatial mesh may be seen by recognizing that a spatial mesh just fine enough to replicate the steep velocity profile through a shock is insufficiently fine to provide a good representation of the velocity gradient itself. A demonstration of the higher accuracy in the high wavenumber spectrum has been given in Appendix C. The spectral transfer T(k, t) = _1 ik[a*(k, t)u2(k, t) - u2*(k, t)£t (k, t)]/NAx

is the Fourier transform of --ld[(u(x, t)u2(x + r, t ) ) -- (uS(x, t ) u ( x + r, t))]/dr. The cumulative transfer spectra is S ( k = 2nfo/NAx, t) = -

T(k', t) dk" ff-Axe~o T(k = 2 n f / N A x , t),

fo < ½N (k < n/Ax). These terms appear in the spectral form of the von K~trm~m-Howarth equation for Burgers' equation, dE(k, t)/dt = T(k, t) - 2vk:E(k, t).

4. E N E R G Y

SPECTRA

The energy spectrum of Burgers' equation for modified Thomas initial conditions is characterized by (i) a low wave-number white-noise energy-containing subrange, 0 < k < k_2 (Fig. 3), (ii) a k -2 subrange, k_2 < k < kexp (Fig. 3), and (iii) an exponential viscous cutoff, kexp< k (Fig. 4). The small white-noise domain is characteristic of energy spectra evolving from both Thomas [21, 22] and modified Thomas initial conditions; although the modified Thomas has more total energy, the spectra here and in Keleti and Reed [22] have the same low wave-number level. The location of the transition wave number k_2 at which the white noise gives way to the k -2 subrange decreases with increasing time at a given Re0, a qualitative observation which is quantified at 400 and 1000 in Figs 3(a) and (b). Burgers was the first to obtain the k-5 spectrum (eqn (85b), Ref. [9]). The k-5 subrange is usually discussed by comparing and contrasting it with the k -5/3 subrange for three-dimensional isotropic turbulence [24, 42-45]. We elect instead to discuss it in straightforward terms of properties of Fourier transforms. The Fourier transform of Burgers' equation is:

d~(k, t)/dt+

ik[.I~_oo ~(m,t)fl(k--m,t)dml/4n = -vk'~(k,t),

with ~(k, t) the Fourier transform of u(x, t). The asymptotic behavior of ~(k, t) may be seen by considering the limits k ~ 0 and k ~ oo in this equation. The former yields: lim dfi(k, t )/dt = O,

k~0

STEVEN KELETI and X B REED JR

154

o.1-__. (a)

3

-5

5 / 10 15 25

-10

w

o o ._1

/ 1

/ 3 5

-15

IO

/

15 / 2O 2G

-20

-25

I

I LIIIII

I

I IIIIIl

0,1

I

I

I llllll

10

i

l

I I IIIJl

10

I

l

1 illlJI

100

1000

k

o,

(b)

o

-5 101

"' "°~"~ - l O

~°5

/

o

20

25

-15

I l I Illj

01

I

I

I IIIIII

I

[

I I I1111

1.0

10

I

I

i

Illlj

I

100

I

k

I Ijllll

1000

k

Fig. 3. Evolution of energy spectra E(k, t) for t = 0, I, 3, 5, 10, 15, 20 and 25. Log-log plot emphasizing white noise and k-2 subranges: (a) at Reo = 400; (b) at Reo = 1000.

155

Spectral properties of exact random solutions to Burgers' equation 0,

(e)

-5

-10 v IJJ o o~

o

-15 25

-20

-25 11 0

100

I

t

I

I

I

I

I

I

I

I

200

500

400

500

600

700

800

900

t000

1100

I 700

I 800

I 900

I 1000

k 0 ¸ (b)

-5

o _1 -1 0

25 -15 0

I 100

I 200

I 300

I 400

I 500

I 600

I 1100

k Fig. 4. Evolution o f energy spectra for t = 0, I, 3, 5, 10, 15, 20 and 25. Log-,~,mi plot emphasizing exponential viscous cutoff: (a) at Re0 = 400; (b) at R% = |000.

156

STEVENKELETIand X B REEDJR

under the mild restriction that the energy density of the field be bounded, i.e. lim (1/2L) L ~ ,~c

u 2 (x, t)dx < C. L

Although this implies nothing about the wave-number dependence of E(k, t ) = ~*(k, t)~ (k, t) for small k, it does indicate that the large scale structures are slowly varying and that the mean velocity is time invariant, i.e. ~(0, t) = ( u ) = const. The latter yields a 1-D version of the weak turbulence limit: lim ~(k, t) = fi(k, O)exp(-vk2t), k~o~

provided the inertial term (the convolution term with prefactor k) grows less rapidly than the dissipation term (the energy term with prefactor k2). This is a variant of plausible arguments sometimes invoked for spectral behavior at high wave numbers [24, 44, 45]. Unfortunately for such arguments, there is the further implication that lim E(k, t) = E(k, 0)exp(-2vk2t), k

~oo

in clear disagreement with the results in Fig. 4. The exponential viscous cutoff makes clear that convolution terms cannot be ignored in the high wave-number limit and must be handled with care in the theoretical estimation of spectral behavior. Consequently, equations which have solutions which behave like solutions to Burgers' equation-equations which are unforced or are weakly forced on scales larger than typical intershock distances (macroscales)--lead to fallacious high wave-number spectra if the convolution term is neglected on the basis of physical arguments, dimensional analysis, and ordering arguments predicated only upon local scales. Finally, the mathematical form of the high wave-number energy spectra evolving from Thomas [21, 22] and from modified Thomas (see especially Fig. 4) initial conditions agrees with that for a single Burgers' "eddy" [46, 47]; Mizushima [30] also obtained an exponential viscous cutoff as a consequence of similarity and invariance of the large scale motions. Convolution produces broadening and smoothing effects (pp. 143-149 of Ref. [48]). For a function y(x), Y(k)=

y(x)exp[-ikx]dx= -

oo

y'(x)exp[-ikx]dx -

ik

oo

follows from integration by parts and the periodicity o f y ( x ) (or from (y(x) ~ 0 as Ix[ ---,oo) and from y'(x) being non-impulsive (or from y'(x) being absolutely integrable). Thus, k Y ( k ) - ~ O as k --* oo and k Y ( k ) ,~ 0([k[ -I) or Y(k) ,,, 0(tk[-2). Repeating the procedure, if the mth derivative of a function becomes impulsive, then for large k its transform behaves as [kl -m and its energy spectrum as k-2m. As our solution to Burgers' equation for v > 0, t > 0 is everywhere smooth (all derivatives exist), its energy spectra must fall off faster than any polynomial at large k. This the above exponential viscous cutoff does. The energy spectrum of a function with corners (discontinuities in the function derivative) has an asymptote of k-4 [as does the Thomas initial condition E(k, 0)], whereas one with jumps (discontinuities in the function) has an asymptote of k-2 (which is the asymptote to which the energy spectra E(k, t) quickly evolve in the absence of viscosity). Modified Thomas initial conditions have mostly corners, but with {dn} selected on (0, 2), jumps are not a priori excluded (as they are for Thomas initial conditions), and one sees a transition from a k -4 asymptote to a k -2 asymptote in Fig. 3 at t = 0. The k -2 subrange passes smoothly into the exponential viscous cutoff at a transition wave number k,xp, which decreases with time. The bandwidth of the k 2 subrange remains virtually constant (Fig. 3), with the migration of the center frequency to the lower wave numbers thus accounting for the energy decay.

Spectral properties of exact random solutions to Burgers' equation

157

An exponential viscous cutoff for 3-D isotropic turbulence has been proposed by Dugstead (p. 279 of Ref. [43]) and derived by Tatsumi and Kida [31] for the modified 0-4th cumulant expansion. Both Corrsin [49] and Pao [50] have obtained exponential functions for high wave numbers with, however, powers of k in the exponential; Tatsumi, Kida and Mizushima [51] obtained exp ( - f l k t.5) by multiple-scale cumulant expansion. The so-called Gaussian falloff has also been derived by Saffman [52], and is widely used in studies of Burgers' equation (e.g. Refs [15, 33, 53]). The exponential viscous cutoff for Burgers' equation dates from one of his early analyses (eqn (85a) [9]) in the form A/sinh2flk; Benton [46], Saffman [47] and Mizushima [30] subsequently obtained the same form. Tokunaga's [54] two-equation 1-D compressible turbulence investigation led to an exponential visocus cutoff, as have several closure hypotheses (e.g. Ref. [33]). Love [55], in the context of subgridscale modelling of Burgers' equation, shows without comment a log-semi plot with a high wave-number band which is nearly linear. Other than Love, there have been no computations into the high wave-number domain until Keleti and Reed's [21, 22] studies of the Thomas initial condition, although a belief in the exponential cutoff was widespread. The situation may be likened to that for the Navier-Stokes equations, except that the requisite delicate measurement with sufficiently high spatial and temporal resolution of fluctuating velocity fields at sufficiently high Re~ had been lacking. Indeed, until 1962, the k-5/3 spectrum had not been comfortably confirmed [38]. The difficulty cannot be overcome through direct numerical simulation of the Navier-Stokes equations because of the hardware requirements of high Re~ computations. 5. D I S S I P A T I O N S P E C T R A

With D(k, t) given by k2E(k, t) and with E(k, t) already described, little needs to be said about the dissipation spectrum. It is characterized by: (i) a quadratic low wave-number domain (Fig. 5), (ii) a white-noise subrange (Fig. 5), and (iii) a viscous cutoff of k 2 exp(-~tk) (Fig. 6). Dissipation spectra for Burgers' equation have customarily been discussed in terms of a uniform distribution of dissipation described in (ii). Equipartition of dissipation obviously requires two caveats, viz. (i) and (iii), with (i) conceivably restricted to Thomas and modified Thomas initial conditions. The dissipation spectrum D(k, t) is of special interest in the high wave-number range, for which computations for small ensembles and relatively coarse meshes [16-20] were admittedly lacking until recently [21, 22]. 6. S P E C T R A L T R A N S F E R

AND CUMULATIVE

SPECTRAL TRANSFER

Spectral transfer along the energy spectrum is at the heart of an understanding of the mechanics of homogeneous turbulence [56, 24, 43, 44]. Two-point closures which focus on the structure of homogeneous turbulence require not only empirical results for the energy spectrum, they require them also for the transfer spectrum. The temporal evolution of the two-point spatial correlation of third order, (u2(x, t)u (x + r, t)) has been computed using the Lagrangian history direct interaction approximation of Kraichnan [11] by Walton [57] and has been computed from exact solutions [19] on what, by comparison to the present spatial mesh, may be described as coarse; and as there was an ensemble of only 12 realizations, no attempt was then made to compute transfer spectra. Subsequently, T(k, t) were computed for single-point Gaussian initial conditions [58], but those results were not adequate for the present purposes because they were also for a small ensemble (12 realizations), as well as being limited in k-space and in accuracy. The spectral transfer function T(k, t) is shown in Figs 7 and 8. The former emphasizes how T(k, t) evolves in the k-2 energy subrange, whereas the latter emphasizes spectral

158

STEVENKELETIand X B REEDJg

O'

(a)

-5

i" 3'

20 25

0 -10

0

3

-15'

-20

I

II

I

~KI

i

L

II1|

I

i

I IIItJ

I

l

i

10

10

0.1

I IIII

I

I

I

I

I lillJ

100

1000

k

(b)

o

)

, I015 2O

-10

I

I

II

II[I 0.1

I

[

I

illll

i

I

I

I I1 I l l

1.0

10

i

I

I I lJlll

I 1 O0

t

I

lillll 1000

k

Fig. 5. Evolution of dissipation spectra D(k, t) ffi kiE(k, t) for t ffi 0, 1, 3, 5, 10, 15, 20, 25. Log-log plot emphasizing ¢quipartition subrang¢ for dissipation: (a) at Re0 = 400; (b) at R% ffi 1000.

Spectral properties of exact random solutions to Burgers' equation

159

-5

¢3 o

-10

0

-15

-20.

l

l

I

I

l

l

i

L

i

J

1

J

0

100

200

300

400

500

600

700

800

900

1000

1100

k

O- (b)

0

-5

,0,,;0//

-

25

-100

I 100

I 200

1 300

I 400

1 500

I 600

I T00

I 800

[ 900

I 1000

I t100

k

Fig. 6. Evolution of dissipation spectra D(k, t)•k2E(k, t) for t •0, I, 3, 5, 10, 15, 20, 25. Log-semi plot emphasizing effect of viscosity on dissipation at high wave numbers: (a) at Reo ffi 400; Co) at Re o -- 1000.

160

SamwN o

KELETIand X B REEDJR

(a)

-5

1

/

3

/ 5

o

1 10

-10

15

/

25

-15

i

-20

l i

iil

I

I

I I tllli

0.1

1

t

[

111111

I .0

I

I

I I 1111t

i 0

I

I

I [lllll

I00

I 000

k

o. ( b )

-5

I-o _I

/

lO

- i O.

15

/ 2o 25

-45

I I

F i g . 7. E v o l u t i o n

ikL 0.1

]

of spectral

i

1111111 1.0

transfer Re0

t

i i illlll

i I 0

t

i tlllkl 100

I

i

i ii~ll[ 1000

T(k, t) for t = 0, 1, 3, 5, 10, 15, 20, 25. L o g - l o g plot: (a)

at

= 400; (b) a t R e 0 = 1000.

transfer in the viscous cutoff range. The effect of increasing Re0 is to extend the far wave-number range into which energy is transferred from the low wave modes. The higher the value of Re0, the more T(k, t) is dominated by inviscid mechanics typified by E(k) ~ k -2 [cf. Figs 7(b) and (a)]. The more rapid falloffof T(k, t) with increasing viscosity in the viscous cutoff range is clear by comparing Figs 8(a) and (b). Initial spectral transfer

161

Spectral properties of exact r a n d o m solutions to Burgers' e q u a t i o n

o]

(al

-5-

r--io-

o cn 0

0 /

._I

i 3

5 / I0

15

25

-15-

]

[

I

I

i

I

[

I

[

I

I

100

200

300

400

500

600

700

800

900

1000

1100

-20 i l 0

k

o

(b)

-5

o

-10

2o

/

25 -I 5 0

[

I

I

i

I

I

I

I

[

I

I

100

200

500

400

500

600

?'00

800

900

1000

1'100

k

Fig. 8. Evolution of spectral transfer

T(k,t) for t = 0 , 1, 3, 5, 10, 15, 20, 25. (a) at Reo = 4 0 0 ; (b) at Re0 = 1000.

functions T(k, 0) are quite noisy, but their general location (Figs 7 and 8) provides a ready reference for comparing T(k) at different Re0, as well as different t, much as E(k, 0) does for E(k, t; Re0). Because T(k, t) and S(k, t) take on negative values, log-log and log-semi coordinates can be used only if the absolute values are plotted or if the negative values

162

STEVENKELETIand X B REEDJR

are omitted. The latter option was chosen here. The log-log plots (Fig. 7) of T(k, t) do not show any power law subranges but do reflect the more extensive viscous cutoff at Re0 = 400 and the more extensive inviscid subrange at Re0 = 1000. The log-semi plots (Fig. 8), in contrast, suggest exponential cutoffs for T(k, t), just as E(k, t) has an exponential viscous cutoff. However, the log-semi plots of dissipation spectra (Fig. 6) also give the appearance of being exponential but are known to be of the form k2E(k) ~ k 2 exp(-ilk). Upon closer scrutiny, several of the T(k, t) lines (Fig. 8) may be seen to have a slight curvature, comparable to that of the D(k, t) plots (Fig. 6). Consequently, one infers that T(k, t) ~ Pm(k)exp(-flk) for high wave numbers, with Pro(k) a polynomial in k. For a single "eddy", the asymptotic curve according to Benton [46] would have the form (k 2 - k ) e x p ( - f l k

).

(a)

0 1 3 5 10

-5

20" 25

o

-I 0

0 .-I

-15

-20

~ illili

i

I

i l llilJ

0.1

t

I

I llllll

1.0

I

i

I llllll

10

ill

100

1000

k

o

-10

{b)

I

J111ll

i

0.1

I

J lilllJ

I

J

J IllllI

1.0

I

10

l

l I ilili

I

I O0

l

l l

llJ

1000

k

Fig. 9. Evolution of cumulative transfer spectra S(k, t) -SkoT(k ', t)dk' for t = 0, 1, 3, 5, 10, 15, 20 and 25: (a) at Re 0 = 400; (b) at Re 0 = 1000. =

Spectral properties of exact random solutions to Burgers' equation

163

(a)

-5.

o -10. O

_1 10' 15 20' 25

-1 5.

-20

I1

1

I

I

I

I

I

I

I

I

I

I

0

100

200

300

400

500

600

700

800

900

1000

1100

k

o

(b)

25 -15 0

I

I

I

1

I

I

I

I

I

I

I

100

200

300

400

500

600

700

800

900

t000

1100

k

Fig. 10. Evolution ofcumulative transfer spectra S(k, t) = -S~T(k', t)dk" for t = O, 1, 3, 5, 10, 15, 20 and 25: (a) at Re o --400; (b) at Re 0 = 1000.

-0.11

-0.1 0

-0.09

-0.08

-0.07

-0.06

-0.05

--0.04

•

--0.0:3 •

--0.02

-0.01

0.00

0.01 .

0.02.

0.05

0.04

005-

0.06.

0.07'

0.08-

0.09-

01 0

0.1 1

01

T(k)

\

0.5 0.3 0.2 0. I 0 \

10

S{k)

k

~

10

0 I

100

2 3

1100

2.

•

-0.021

-0.018

--0.015 •

-0.01

-0.009

-0.006-

-0.003.

0.000

0.003-

0.006-

0.009"

(b)

0.1

T(k)

1 .0

S(k)

k

"10 "15 "20 "25

I 0

15

100

//2025

10

Fig. II. Evolution o f transfer a n d c u m u l a t i v e transfer spectra. S e m i - l o g plot at Re 0 = 400: (a) t = 0, 0.1, 0.2, 0.3, 0.5, 1, 2, 3 a n d 5; (b) t = 10, 15, 20 a n d 25.

(a)

0.01 2.

o.o15.

1000

X

,,..,

Z

F:

0.1

1.0

-0.11

-0.t0.

•

-0.09 •

-0.08 -

k

100

1 2

3 5

1000

- 0.021

- 0 . 0 1 8

-0.015

-0.012

-0.009

-0006

-0.003

0.000

0.003.

0006

0.009

,

Ol

, , , , , i

(b)

,

,

,

, , , ,

1%

,

,

,

k

, , , i ,

20 25

15

, ~, dl,~l"" J ' ~ " - ~

1'o

,

.......

,

,

, , , , ,

1~o

,

,

i015

.

Fig. 12. Evolution of transfer and cumulative transfer spectra. Semi-log plot at Re0= 1000: (a) t =0, 0.1, 0.2, 0.3, 0.5, 1, 2, 3 and 5; (b) t = 10, 15, 20 and 25.

10

-0.06

-0.07

'1 '3 "5

T(k)

0\\

0.2\ 0.1\

6.3\

0.5

S(k)

'0

] (a)

-0.05

- 0 . 0 4

-0.03

-0.02

-0.01

0.00

0,01

0.02

0.05

0.04

0.05

0.06

0.07

0.08

0.09.

0.I0

0.11

0.012.

o.o15.

. . . .1. .0 0 0

o

.8

o

=..

B

o

o

E

166

STEVENKELETIand X B REEDJR

The cumulative spectra are shown in Figs 9 and 10, with similar comments holding. Semi-log plots of both T(k, t) and S(k, t) are displayed in Figs I I (Re0 = 400) and 12 (Re0 = 1000), for short times in Figs ll(a) and 12(a) and longer times in Figs l l(b) and 12(b). The forms, except for the "noisy" nature of T(k, t), are quite like those for 3-D transfer spectra [35, 36, 59, 60, 51]. The most spectacular observation has little physics but reflects the well known smoothing role played by integration. More stringent windowing-or french curve smoothing--would yield T(k, t) vs k resembling even more the published less noisy 3-D results. The only window we have used for T(k, t) is that used throughout, as in Ref. [22], which does not affect the visual character of the curves. Aliasing effects are seen at the highest wave numbers and a discussion of the windowing and aliasing effects is presented in Appendix C. For the modified Thomas initial conditions, the spectrum of initial eddy sizes assures that for small k, its spectral transfer is greater than that for Thomas initial conditions (due to the greater energy distributed in the larger "eddies"). This "signature" of the modified Thomas initial conditions persists through all of the times for which the spectra were evaluated. The other feature distinguishing T(k, t) for modified Thomas from Thomas initial conditions is the lack of spikes for small times.

7. DISCUSSION The modified Thomas initial conditions introduced here provide for a distribution of initial spatial scales, whereas a period of evolution is necessary for Thomas initial conditions to evolve into a state having a distribution of spatial scales. This has consequences for the energy decay rate (Fig. 2), for the evolution of the Taylor microscale 2(t), and for the corresponding turbulence Reynolds' number Re~(t)= u'2/v (Tables l and 2). The implications of the distribution of initial spatial scales--which may be likened to initial eddy sizes--on the initial turbulence energy spectrum are profound. For all but short times, nevertheless, the energy spectra for the two classes of initial conditions are qualitatively quite similar. This suggests that they are two-point statistical properties which, when appropriately normalized, are characteristic of random solutions to Burgers' equation ("turbulence") and not characteristic of initial conditions ("specific grid geometry") [24]. These spectral results (Figs 3 and 4), together with those for Thomas initial conditions [22], thus invite formulation of a conceptual framework which would provide a unifying mathematical description of the results for Burgers' equation which would be (virtually) independent of initial conditions. The dissipation spectra for the two classes of initial conditions, being simply and directly related to the energy spectra, have corresponding similarities and dissimilarities. The most noteworthy feature of D(k, t) for modified Thomas initial conditions is the dissipative equipartition on k 2 < k < ke,p, which it has in common with that for Thomas initial conditions [22]. These observations themselves are hardly new. What is novel is the discovered limit kexp(t,Re0) and the quantitative demonstration of the asymptote k2exp(-ot(t)k) as a function of Re0. Spectral transfer results for Burgers' equation have seldom been presented (but see, e.g. Ref. [12]), other than very limited ones [58], before our study of Thomas initial conditions [21, 22]. Except for short times, spectral transfer for Thomas [22] and for modified Thomas are hardly dissimilar. Initial transients, as expected, are strongly dependent on initial conditions. In conclusion, through the use of numerical procedures to compute the exact solution to Burgers' equation with random initial conditions, spectral solutions over a range of 20 decades--up to 16 decades in the viscous cutoff at Re0 = 400--have been computed for the modified Thomas initial conditions. These results reveal unequivocally an exponential viscous cutoff for Burgers' equation.

Spectral properties of exact random solutions to Burgers' equation

167

Acknowledgements--The first author was supported during the course of this research by a Shell Fellowship, for which we are both indebted. It is also a pleasure for us to acknowledge the encouragement of J. W. Johnson, and the cooperation of D. W. Dearth, both of whom helped in different ways to hasten completion of the research. Finally, we are indebted to H. Haddad who carried out extensive computations on earlier versions of the velocity and spectral codes [21].

REFERENCES I. J. H. Thomas, Numerical experiments on a model system for magnetohydrodynamic turbulence. Phys. Fluids 11, 1245-1250 0968). 2. J. D. Cole, On a quasi-linear parabolic equation occurring in aerodynamics. Q. appl. Math. 9, 225-236 (1951); [see also P. A. Lagerstrom, J. D. Cole and L. Trilling, Problems in the theory of viscous compressible fluids. Cal. Inst. Tech., Pasadena, CA (1949)]. 3. E. Hopf, The partial differential equation ut + uux = vu,~. Commun. pure appl. Math. 3, 201-230 (1950). 4. E. R. Benton and G. W. Platzman, A table of solutions of the one-dimensional Burgers equation. Q. appl. Math. 30, 195-212 (1972). 5. E. Y. Rodin, On some approximate and exact solutions of boundary value problems for Burgers' equation. J. Math. Anal. Appl. 30, 401-414 (1970). 6. X B Reed Jr and Y.-C. Shih, Application of the Cole-Hopf transformation to boundary value problems for Burgers' equation. To be submitted. 7. J. M. Burgers, Mathematical examples illustrating relations occurring in the theory of turbulent motion. Verh. Kon. Nederl. Akad. v. Wetenschaffen, Amsterdam, Afd. Natuurk. 17, 1-53 (1939). 8. J. M. Burgers, A mathematical model illustrating the theory of turbulence. Advances in Applied Mechanics, Vol. 1, (Edited by R. yon Mises and T. yon Karman), pp. 171-199 Academic, New York (1948). 9. J. M. Burgers, Correlation problems in a one-dimensional model of turbulence I-IV. Proc. Acad. Sci. Amsterdam 53, 247-260, 393-406, 718-742 (1950). 10. J. M. Burgers, The Nonlinear Diffusion Equation. Reidel, Dordrecht (1974). 11. R. H. Kraichnan, Lagrangian-history statistical theory for Burgers' equation. Phys. Fluids 11, 265-277 (1968). 12. W. H. Reid, On the transfer of energy in Burgers' model of turbulence. Appl. Sci. Res. A6, 85-91 (1957). 13. H. W. Wyld Jr, Formulation of the theory of turbulence in an incompressible fluid. Annls Phys. 14, 143-165 (1961). 14. W. C. Meecham, P. Iyer and W. C. Clever, Burgers' model with a renormalized Wiener-Hermite representation. Phys. Fluids 18, 1610-1616 (1975). 15. W. C. Meecham and A. Siegel, Wiener-Hermite expansion in model turbulence at large Reynolds numbers. Phys. Fluids 7, 1178-1180 (1964). 16. D. T. Jeng, A Study of Burgers' model equation of turbulence. Magistral Thesis, University of Minnesota (1965). 17. D. T. Jeng, R. Foerster, S. Haaland and W. C. Meecham, Statistical initial-value problem for Burgers' model equation of turbulence. Phys. Fluids 9, 2114-2120 (1966). 18. Y.-C. Shih and X B Reed Jr, Solution to the piecewise linear continuous random initial value problem for Burgers' equation. Phys. Fluids 28, 2088-2099 (1985). 19. Y.-C. Shih and X B Reed Jr, Exact solution to Burgers' equation for spatially periodic boundary conditions. Physicochem. Hydrodyn. 6, 853-861 (1985). 20. Y.-C. Shih and X B Reed Jr, Solution to Burgers' equation on a finite domain with random initial and zero boundary conditions. Symposium of Fluid Dynamics, (Edited by A. L. Addy, W. L. Chow, N. S. Vlachos, and R. A. White), pp. 301-338. Department of Mechanical and Industrial Engineering, University of Illinois (1984). 21. S. Keleti and X B Reed Jr, Behavior of energy and dissipation spectra and of spectral transfer for exact random solutions of Burgers' equation for initial turbulence Reynolds' number of 400. 5th Symposium on Turbulent Shear Flows, Cornell University (1985). (Abstract available from authors.) 22. S. Keleti and X B Reed Jr, Evaluation of spectra for large ensembles of exact solutions to Burgers' equation for Thomas initial conditions. To be submitted. 23. S. Keleti and X B Reed Jr, Exact solutions to Burgers' equation at high Reynolds' numbers: II: Efficient accurate numerical evaluation. Devs Mech. 14a, 63~58 (1987). 24. G. K. Batchelor, The Theory of Homogeneous Turbulence. Cambridge University Press, Cambridge (1953). 25. R. C. Tolman, The Principles of Statistical Mechanics, 1st edn. Oxford University Press, London (1959). 26. X B Reed Jr and Y.-C. Shih, A serendipitous discovery. To be submitted. 27. D. E. Knuth, The Art of Computer Programming, Vol. 2 Seminumerical Algorithms, pp. 9-37 Addison-Wesley, Reading, MA (1981). 28. J. M. D. MacElroy and S. H. Suh, Concentration and noncontinuum effects in size-exclusion partitioning. AIChE Syrup. Ser. No. 248, 82, 133-148 (1986). 29. P.-L. Tseng and X B Reed Jr, Exact solutions to Burgers' equation at high Reynolds' numbers: I. Inviscid solutions. Devs Mech. 14a, 57-62 (1987). 30. J. Mizushima, Similarity law and renormalization for Burgers' turbulence. Phys. Fluids 21, 512-514 (1978). 31. T. Tatsumi and S. Kida, The modified cumulant expansion for isotropic turbulence at large Reynolds numbers. J.Phys. Soc. Jap. 49, 2014-2025 (1980). 32. S. Kida, Asymptotic properties of Burgers turbulence. J. Fluid Mech. 93, 337-377 (1979). 33. J. Mizushima and T. Tatsumi, The modified zero-fourth cumulant approximation for turbulence. J. Phys. Soc. Jap. 50, 1765-1773 (1981).

168

STEVEIq KELETI and X B REED JR

34. J. Mizushima and Y. Saito, Numerical study of the initial condition dependence of Burgers' turbulence. Phys. Fluids 28, 1294-1298 (1985). 35. C. W. Van Atta and W. Y. Chen, Measurements of spectral energy transfer in grid turbulence. J. Fluid Mech. 38, 743-763 (1969). 36. M. S. Uberoi, Energy transfer in isotropic turbulence. Phys. Fluids 6, 1048-1056 (1963). 37. A. L. Kistler and T. Vrebalovich, Grid turbulence at large Reynolds' numbers. J. Fluid Mech. 26, 37-47 (1966). 38. H. L. Grant, R. W. Stewart and A. Moilliet, Turbulence spectra from a tidal channel. J. Fluid Mech. 12, 241-263 (1962). 39. C. Basdevant, M. Deville, P. Haldenwang, J. M. Lacroix, J. Ouazzani, P. Orlandi and A. T. Patera, Spectral and finite difference solutions of the Burgers equation. Comput. Fluids 14, 23-41 (1986). 40. M. Martell, Computing very large FFTs efficiently on the FPS-164. Checkpoint 2, pp. 5-18. Floating Point Systems Inc. (Oct. 1970). 41. C. W. Van Atta and W. Y. Chen, Correlation measurements in grid turbulence using digital harmonic analysis. J. Fluid Mech. 34, 497-515 (1968). 42. A. N. Kolmogorov, The local structure of turbulence in incompressible viscous fluid for very large Reynolds' numbers. C.r. l'Acad. Sci. I'URSS 30, 301-305 (1941). [See pp. 151-161 of Turbulence: Classic Papers on Statistical Theory (Edited by S. K. Friedlander and L. Topper). Interscience, New York (1961)]. 43. J. O. Hinze, Turbulence, 2nd edn. McGraw-Hill, New York (1975). 44. A. S. Monin and A. M. Yaglom. Statistical Fluid Mechanics, Vol. 2, (Edited by J. L. Lumley). MIT Press, Cambridge, MA (1975). 45. H. Tennekes and J. L. Lumley, A First Course in Turbulence. MIT Press, Cambridge, MA (1972). 46. E. R. Benton, Solutions illustrating the decay of dissipation layers in Burgers' nonlinear diffusion equation. Phys. Fluids 10, 2113-2119 (1967). 47. P. G. Saffman, Lectures on homogeneous turbulence. In: Topics in Nonlinear Physics (Edited by N. J. Zabusky), pp. 485-614. Springer, New York (1968). 48. R. Bracewell, The Fourier Transform and Its Applications. McGraw-Hill, New York (1965). 49. S. Corrsin, Further generalization of Onsager's cascade model for turbulent spectra. Phys. Fluids 7, 1156-1159 (1964). 50. Y. H. Pao, Structure of turbulent velocity and scalar fields at large wavenumbers. Phys. Fluids 8, 1063-1075 (1965). 51. T. Tatsumi, S. Kida and J. Mizushima, The multiple-scale cumulant expansion for isotropic turbulence. J. Fluid Mech. 85, 97-142 (1978). 52. P. G. Saffman, On the fine-scale structure of vector fields convected by a turbulent fluid. J. Fluid Mech. 16, 545-572 (1963). 53. K. Yamamoto and I. Hosokawa, Energy decay of Burgers' model of turbulence. Phys. Fluids 19, 1423-1424 (1976). 54. H. Tokunaga, The statistical theory of one-dimensional turbulence in a compressible fluid. J. Phys. Soc. Jap. 41, 328-337 (1976). 55. M. D. Love, Subgrid modelling studies with Burgers' equation. J. Fluid Mech. 100, 87-1 I0 (1980). 56. W. Heisenberg, Zur statistischen Theorie der Turbulenz. Z. Phys. 124, 628-657 (1948). 57. J. J. Walton, Integration of the Lagrangian-history approximation to Burgers' equation. Phys. Fluids 13, 1934-1935 (1970). 58. Y.-C. Shih and X B Reed Jr, Can exact random solutions to Burgers' equation be normally distributed? Int. Conf. of Physicochem. Hydrodyn., Oxford University (1987). 59. I. Proudman and W. H. Reid, On the decay of a normally distributed and homogeneous turbulent velocity field. Philos. Trans. R. Soc. A247, 163-189 (1954). 60. K. N. Helland, C. W. Van Atta and G. R. Stegen, Spectral energy transfer in high Reynolds number turbulence. J. Fluid Mech. 79, 337-359 (1977). 61. A. V. Oppenheim and R. W. Schafer, Digital Signal Processing. Prentice-Hall, Englewood Cliffs, NJ (1975).

APPENDIX A Derivation o f Initial Turbulence Energy The expectation for the initial turbulence energy for the Thomas and modified Thomas initial conditions can be determined from the expectation of the energy for an individual interval, as the velocity endpoints are statistically independent. Taking the linear velocity defined between the endpoints u,_ ~ and un over the segment d,, one has at t = 0 an energy density of: En=

[(Un - - U n

i)x/d.+u.12dx,

=dn[u~ +u.u. , + u~_d/3. From statistical independence and the assumption that the endpoints are selected with zero mean, one obtains an energy of: E(t = O) = (En)/(dn) = (2/3)(u:,) = 2/9, because (u~,)= I/3 for the uniform distribution of the Thomas and modified Thomas initial conditions.

Spectral properties of exact random solutions to Burgers' equation

169

APPENDIX B Computational Algorithm The difficulties with numerical evaluation of exact solutions to Burgers' equation arise primarily from the fact that the velocity is finite but the numerator and denominator o f the Cole-Hopf transformation may individually underflow or overflow. Consequently, the computation of the velocity at x = xpo~t~onj using the solution appearing in Ref. [18] requires an accurate estimate of the magnitude o f the largest terms ~n the summation, for not all terms are numerically significant. Indeed, the most important terms are those corresponding to the domain surrounding the solution position Xpo~t~on~. Thus, an estimate of the largest terms to within sufficient orders of magnitude to ensure that factorization oiz numerator and denominator produces computations falling within the dynamic range o f the computer can be obtained from these terms. The initial conditions are repeated beyond + L so as to allow the computation of a periodic solution. The algorithm requires subroutines for the accurate computation of the following functions:

erfc(x)=(2~/(~)) f f e x p ( - t 2 ) d t serfc(x) = (2/~/(g)) exp(x 2)

F(x) = e x p ( - x 2)

e x p ( - t 2) dt

exp(t 2) dt

exp,(x) = e x p ( x ) - 1 (computed accurately for small x) ilog2(x ) = integer [log2(x)] = base 2 exponent field of floating point number x 2x = floating point 1.0 with x in exponent field. The arguments of the routine are: x n = left bound of linear segment n of N segments;

xl'= - N / 2 ; xn+~'= xn + dn; a n = slope o f linear segment n of N segments; an'= (Un+~ - un)/d.;

bn = intercept of linear segment n of N segments;

bn'= un - andn; cn = integration constant for segment n o f N segments; c I := O;

cn'= cn-I + x.(x.(a. -- a n_ 0/(4v) + (b. - bn_ 1)/(2v)); (or, if segments are piecewise continuous);

cn,= cn_ i + x.(bn - bn_ 0/(4v); xvosition.j= j t h solution position, x~ < Xposition,j<: XN; Jmax = number of desired solution positions;

uj,= velocity at xposit~on,j ; u~,j,= velocity gradient at Xposition,j; and the algorithm is as follows: begin compute-velocity (an, bn, xn, xposition.j,uj, u~,j, t, v, N, Jm.0; '*constants; nb~ts= 51; /* number of unsigned mantissa bits of machine; d~man~t= - I 0 ~ ;

/* most negative machine number;

e2 = log~(2); ed~.a~= 0.8; din. =

-- drain;

7rI = ~/(n)/2;

e~. = 2.0;

edmin= 1.4; e 3 = (?lbits +

drain = --2nbit~e2; 1)e2;

?test'= 1;

:*invariant expressions;

s:=~/(vt);

s_t,=l/(2sO;

v_l,=l/(4v);

s2,=2vt;

s_2,=l/s2;

s3:=2v;

/*The range of effect Ax is estimated from the movement of a viscous shock front such that x > (t + ve3)suplu I(x, t = 0)1; Ax.'= 1.5(t + re3);

170

STEVEN KELETI and X B REED JR

/*loop to evaluate Jm~x velocity positions; nstart: =

I;

for j : = 1 to Jm~ by 1 do; Xst=XX;

X : = Xposition.j;

Xstartt=X - - A x ;

Xend: = X + A X ;

/*find largest n such that x. < Xstart; n : = nstart~

if x. ~ x~t~. then do; until (x. > X~ta.) n ,= n + l;

n,= n - 1; end;

else

until (x. ~< x~t.a) n:=n - I;

/*find largest n such that x. < x; until (x. > x ) n , = n + 1;

nposition t = n - - l ;

ds,= dsmallest;

/*estimate scale of terms; for

n,=

nposition=

nes t

to

nposition+ n~ t by 1;

g~=ant + 1; fl'=X/(lg

g_l,= l/g;

h , = b . t - x;

df:=G-(a.xs

.1);

+ b.(x -h))g

kl:=(gXn + h)fls 1;

k~,=(gx.+t + h )~ls_a;

ktd=klk6

kzc"= k2k2;

i f g I> 0 then do;

iv i;

e~,= min(kls, k2s); if ktk 2 <~ - e , . m then G,=0; if e s > emin then es,=e ~+ efilog2(G); dc,= df - e~ + e2 ilog2~); end;

else do;

G,= max(k~s, k2~); if e~ < edm~xthen ed= G + e2 ilog2(G); if e~ > edmin then ed=e ~-- e2 ilog2(e~); dd= df + G + e2 ilog2(fO; end;

@ = max(d., de);

end for n; dT,=0;

dr'=0;

dg:=0;

n,=nstart;

/*loop to perform summation of terms; until x n > X~nd do; g'=ant + 1;

g_ ~,= 1/g;

h,=bnt - x;

f3'= [g_l[;

f.,=./@;

.~:=slg 1;

A , = s2f, A; k,,=(gx.+h)fls

df:=c

,;

kw=ktk 6

n - - ( a n X s + bn(x - h ) ) g _ l v

k2:=(gx~+t + h ) f t s

I -- ds;

l;

k2s'= k2k2;

/* estimate magnitude of term; if g/> 0 then do;

es,= min(kl~, k2s); if k~k 2 ~< -emi n then es=0; if es > emi. then es,= e. + e2ilog2(e~); dc,= df - es + e2ilog2(f0; end;

else do;

es,= max(k~s, k2s); if es < eama~ then e~,= es + e2ilogz(G); if e~ > edmin then e~.= G -- e2ilog2(e~); rid= dr + es + e2ilog2(f0; end;

/* evaluate term only if its magnitude is greater than dram; if d~ >/d~i . then do; /* if magnitude is larger than d~a~. rescale terms; if d~ > dm~ then do;

4,=a~ +a~;

4,=a.-4;

@ = e x p ( - d~);

d~,=d#o;

4 = 44;

Gt..,=n;

Spectral properties of exact random solutions to Burgers' equation

171

end; /* evaluate of error function (Dawson integral if ant + 1 < 0) term; if g >/0 then do; kls '= -kl~;

k2s ,= - k z , ;

/* compute: din, = ( ftn ,) exp (dr) (erffk2) - erf(k 1)), where erf(x) = (2/x/(n)) ( f exp(- t 2) dt; /* minimize error by ensuring largest magnitude "erf" term is positive; p:=kl;

p2 ,= k2;

if IPl[ ~<]P21 then do; if p~ ~<0 then do;

p:= -k2;

p2,= - k t ;

end;

Pc:= - P 6

P : = -P2;

p2:=p¢;

end; q~,= O;

wl,= erfc(pO;

end;

q := - P t P 6

w : = serfc(pt);

end;

if P2 ~<0 then do; else do;

q2'= 0;

w2,= erfc(P2);

end;

q:'= -P2P2;

w2'= serfc(P2);

end;

if q~/> q2 then do;

er,= df "t- ql;

q3 '= q2 -- qt;

fp:= -- W2;

q3'=ql - q2;

fp== wt;

if p~ ~<0 then do; else do;

end; else do;

ef:= dr + q2;

end; din:= (flnt)exp(er)(wl -- w2 +fp expt(q3)); end; else do; /* compute: tim , = ( - f O e x p ( d f ) ( V , ( k 2 ) - rl(kl)),

where Fl(x ) =

exp(t 2) dt; ql,=k2k2;

w?=F(kz);

q2".=ktkt;

w2,=F(kt);

if q~/> q2 then do; else do;

f p , = - w 2 ; end;

el:= df -+- ql;

q3,=q2-qt;

ef:= df q--q2;

q3'=ql - q2; fp:= Wl; end;

dml,=(--fOexp(ef)(wl -- w 2 +fp expl(q3));

/* compute:

din2= (f2)exp(dr)(exp(kt~) - exp(k2~));

/* compute:

din3 = (f3) exp(df) (klexp(kls) - k2exp(k2~));

if kls/> k2s then do; else do;

e p - df + kh;

q~:= k2s - kts;

fc'-- --f2;

fp:= - k 2 ; end;

ef:= df + k2s;

q3,= kls -- k2s;

fc:=f:;

fp,= kl; end;

rib'= expl(q3);

fd=f~db;

db,= exp(eO;

dm2:=fcdb;

.fp,=A(k, - k 2 +fpdb); dm3,= f pdb;

so,= --hg_l;

Sx'=-Sc-- X;

ds,=ds + sxdml +din2;

d9:=d 9 -k de3 q- 2Sxdm2+ (sc(sx -- X) "6 s2g_ 1 "k- X s -- s2)dml ;

d7,=d7 + dml ;

/* reset scale; p,= ilogz(dT);

@ = 2-P;

@ = d~ + ezp;

end; n : = n + 1;

end until; (OJO),= dss_ 2/dT;

end for; end;

(Oxx/O),=dgs_:_JdT;

u:=-s3(Ox/O);

u~.:= s3((Ox/O ) (OJO) - (OJO ));

172

STEVEN KELETI and X B REED JR

The computation o f the velocity field was accelerated by noting that two positive slope domains o f a solution for piecewise linear initial condition are only colinear if either the condition contains colinear segments or the domain overlap. Although it is not excluded, the initial conditions selected do not have any colinear portions. The velocity field was thus computed as follows. The domain was subdivided into several sections. After the solution was computed at these points, the endpoints were checked for being on a section of the solution with a positive slope and having tangents which were colinear (i.e. having relative slopes and intercepts within the accuracy of the computations "~ 10-t°). If the tangents were colinear, then all the points desired within that section were computed using a cubic spline defined by the velocities and velocity gradients o f the endpoints. The sections with non-tangent endpoints were bisected and the procedure is repeated until all the desired points were computed. The accuracy o f the solution was estimated by solving a random flow field u~(x, 0) and the flow field u2(x, O) = - u ~ ( - x , 0) at time t. All the transcendental routines were written to be accurate to the full precision o f the machine. Thus, any error represents the loss of accuracy due to numerical cancellations between the variables which appear in the solution. As u2(x, t) should equal - u ~ ( - x , t), the maximum error was found to be lu2(x, t) + u l ( - x , t)l/suplul(x, t)l < 10 -I°.

APPENDIX

C

Spectral Considerations In order to reduce the numerous spectral data (2 ~6x 8 per figure) to a manageable quantity, the original spectral results were windowed and decimated before being plotted. A normalized Hamming window was arbitrarily selected as a weighting function, such that a windowed spectra Yw(k) was formed from the spectra Y(k ) via M/2

Yw(k = 2nn/NAx) =

~

w(m)Y(k = 2n(n + m ) / N A x ) ,

m = - M/2

where the normalize,d windowing function w(m) is defined by:

w(m) =

0.54 - 0.46 cos(2nm/M)

for M = 2 , 4 , 8

0.54 M + 0.08 . . . . . 2p,and

Iml=0,1 ..... g/2 (for M = 1, Yw(k) = Y(k)). For different domains of the wave number (from n~ to ni+t) window sizes M were selected so as to produce only gentle smoothing o f the spectra. From the windowed spectra, only a manageable subset of the wave numbers were plotted by selecting every qth wave mode from the beginning o f each range r~ corresponding to a wave number k~--0.01nn~. The ranges o f the windowing and decimation used in the figures are summarized in Table 3. As a periodic domain was used in the simulation, the spectral results are free from effects due to discontinuities at the solution boundary ("picket-fence" effects) and from energy losses out of the boundary. The results are, however, subject to aliasing effects (see, for example, pp. 26-30 of Ref. [61]). At each instant in time, the exact solution, a continuous spatial velocity field Uc(X), is discretized on a mesh Ax such that u(x = nAx) = uc(x). The Fourier transforms o f u(nAx) and uc(x) are ~(k) and ~c(k), respectively, such that

u(nAx ) = ~1 ~I ,n/Ax a x ~(k )exp[iknAx] dk,

(IC)

and

•(k)=

~

ac(k +2nr/Ax).

(2C)

Aliasing effects are caused by the contribution of terms other than that for r = 0 in eqn (2C). Because I~(k)l 2 decreases rapidly with increasing k, especially in the exponential viscous cutoff, the two most significant terms in the summation in eqn (2C) are r = - 1 , 0 and a(k) = ac(k) + a¢(2n/Ax - k).

(3C)

Because E(k ) = Ia(k )I2/NAx, lim E(k ) ~-12Re(ao(n/Ax))12/NAx. So, 0 ~< E(k = n/Ax)< 41'~c(lr/Ax)12/NAx = 4Ec(n/Ax). Thus, the aliased energy spectra of the highest wave number, k = n/Ax, could conceivably attain values anywhere from 0 up to 4 times the actual energy spectra there. The effect of aliasing on the spectral transfer function can be precisely determined. Simpifying T(k) and taking Table 3. Windows and decimation used i ni ki M q

1 0 0 1 1

2 300 9.425 2 3

3 400 12.57 4 7

4 1200 37.70 4 13

5 4097 128.7 4 31

6 16097 505.7 4 37

7 32768 1029.0 ---

Spectral properties of exact random solutions to Burgers' equation

173

the high wave-number limit: r (k) = - ½~ [,~* (k ~ ( k ) - u 2, (k)~ (k)]/N

= kIm[~*(k)u2(k)]/NAx, k~m T ( k ) f N ~ I m [ r

and

+2~s)/Ax)].

~_ ~*c((~t+21tr)/Ax)s ~ffi~ u~('"Ot

Because uc(k ) and u~(x ) are real, t~c(k)= ~*(-k), u~(k ) = u~*(-k ), and lira T(k) = ~ k k,x/Ax

" ~ + 2;ts)/Ax) ) ] = O. Im [ 2Re ( E® ,~*((~t+ 2;tr)/Ax) ) 2R e( E® U2(( rffi

--

sffi

--

Thus, the high wave number aliased spectral transfer function always attains the value of zero. The effect of aliasing on the dissipation spectra D(k) can be estimated from the asymptotic viscous cutoff, exp(-~tk). Because E(k)~A[exp(-o~k)+exp(-ct(2;t/Ax-k))], computing the dissipation spectra from the aliased energy spectra results in Dl(k)= k:E(k)~ Ak2[exp(-0tk)+ exp(-~t(2;t/Ax- k))], whereas directly computing the dissipation spectra results in D2(k ) = ~*~(k)~(k)/NAx ~ A[k2exp(-0dc) + (27t/Ax-k)~exp(-ot(27t/Ax-k))]. At the high wave-number end of the spectra, D2(k)-Dl(k)= 4(n/Ax)(Tt/Ax - k)exp(-~t(27t/Ax - k))/> 0. Therefore, estimating the dissipation spectra using k2E(k) results in less error than using ,~*~(k):4x(k). The effect of aliasing on T(k) results in an error in the cumulative spectral transfer function

S(k ) = - ti~ T(k') dk'. As the aliased T(k) is smaller than the unaliased T(k) at high wave numbers, S(k) will be correspondingly greater, as seen in the upturn for high k in Fig. 10. It should be noted, however, that as the integral is estimated in practice by a summation, summation errors could become larger than the aliasing errors if the range of S(k) values exceeds the machine precision. In the results presented here, the number of velocity samples N was large enough (2 j6) to practically eliminate aliasing error. As there is a decade of exponential viscous cutoff even at Reo = 1000, the contamination due to aliasing is minimal, being seen only in the very highest wave numbers. If N were reduced by a factor of two, however, almost all of the viscous cutoffwould be removed for Reo = 1000, but one could still see that the viscous cutoff for Re 0 = 400 and the shape of the other spectral curves would be unaffected.