Spatio-temporal POD of turbulent channel flow

Spatio-temporal Proper Orthogonal Decomposition of turbulent channel flow


Proper Orthogonal Decomposition is applied to the wall layer of a turbulent channel flow () in a novel way, where empirical eigenfunctions are defined in both space and time. Due to the statistical symmetries of the flow, the eigenfunctions are associated with individual wavenumbers and frequencies. Self-similarity of the dominant eigenfunctions, consistent with wall-attached structures, is established. We show that the most energetic modes are characterized by time scales in the range 200-300 viscous wall units. The new decomposition provides a natural measure of the convection velocity of structures, with a characteristic value of 12 in the wall layer. Finally, we show that the energy budget can be split into specific contributions for each mode and provides a new closure formulation for the quadratic terms.

1 Introduction

Proper Orthogonal Decomposition (POD) was first introduced in turbulence by Lumley (1967). Its derivation stemmed from the Karhunen-Loève (KL) decomposition (Loève, 1977) which represents a square-integrable centered stochastic process in the time domain as an infinite linear combination of orthogonal functions. If is defined over a finite range, then


where is stochastic and are continuous functions. It is important to note that represents a stochastic variable and not a sample. Realizations of will be noted . The functions are the eigenfunctions of the covariance function , where the operator refers to expectation with respect to the measure of . In the Karhunen-Loève derivation, the variable corresponds to time, but it could indicate any other variable - such as space.

Lumley (1967) adapted the decomposition to Fluid Mechanics: the samples were constituted by flow realizations, and the ergodicity assumption was used to replace the covariance function corresponding to an ensemble average with one obtained by the time average so that the KL transform was generally applied to space. He considered the spatial autocorrelation tensor , where now and represent the deterministic variable (space) and represents the ensemble average (which is simply the time average here). therefore represents the spatial autocorrelation tensor at zero time lag. As pointed out by George (2017), it is important to realize that the spatial auto-correlation tensor is different from its sampled estimation , where is the number of samples. In general, exact eigenfunctions of cannot be computed since the can only be approximated. However, for spatially homogeneous flows such as the channel flow in horizontal directions (Aubry et al., 1988), exact solutions are known a priori since POD modes are Fourier modes in the homogeneous directions.

A key insight given by George (2017) is that if the process is stationary in time, the KL or POD modes are Fourier modes in the deterministic variable . Following George (2017), it appears that the proper way to apply Proper Orthogonal Decomposition would be to apply it in both time and space, reverting to ensemble average (as in the original definition) to evaluate the covariance function. This will allow us to identify modes that can be associated with individual spatial wavenumbers and temporal frequencies. Such decompositions could be useful to identify key features and instability mechanisms underlying the flow dynamics and eventually attempt to control them. One particularly promising approach for this is the resolvent analysis (McKeon, 2017), which is based on the Navier-Stokes equations. In contrast, the present approach is data-based, and can therefore provide a complementary viewpoint.

In the present paper, we apply fully spatio-temporal Proper Orthogonal Decomposition to the turbulent channel flow. Wall turbulence is characterized by a variety of spatio-temporal scales interacting in a highly complex fashion (Robinson, 1991; Dennis, 2015). It is well known that the flow is characterized by an alternation of high and low-speed streaks aligned with the flow, with a typical spacing in the spanwise direction and a length of wall units (Kim et al. (1971), Stanislas et al. (2008), Jiménez (2013)), which are units based on the fluid viscosity and wall-friction velocity , and will be denoted with a + throughout the paper. These streaks are associated with longitudinal vortical structures pushing low-speed fluid upwards and bringing high-speed fluid downwards, which results in a strong contribution to the Reynolds stress (Kim et al., 1971). This contribution is highly intermittent in space and time, with ‘bursts’ of turbulence production. Eulerian measures of the bursts yield time scales on the order of 100 wall units (Blackwelder & Haritodinis, 1983). It is not clear how this time scale is related to the characteristics of the coherent structures, which are associated with a self-regenerating cycle of about 300 wall units (Jiménez, 2013).

Classic spatial POD has been previously used to investigate the boundary layer (see e.g. Aubry et al. (1988), Moin & Moser (1989), Podvin & Lumley (1998)). Spatial eigenfunctions are determined from second-order statistics. The amplitudes of eigenfunctions are characterized by a mixture of time scales and can only be computed by projecting the full instantaneous field onto the spatial eigenfunctions. In contrast, the new decomposition directly provides spatio-temporal patterns.

In this paper we focus on a relatively moderate Reynolds number , based on the fluid viscosity, channel half-height and friction velocity. At these moderate Reynolds numbers our focus will be on the wall layer as large scales such as those observed by Smits et al. (2011) are not present in the flow. The rest of the paper is organized in the following manner. Section 2 gives details of the methodology for the new variant of POD. Section 3 presents the POD results. The contributions to the turbulent kinetic energy equation associated with each mode are examined in Section 4, followed by a conclusion in Section 5 which summarizes the key observations and results.

2 Description of the procedure

2.1 Numerical database

We use the finite-volume code SUNFLUIDH to simulate the incompressible flow in a channel. Further details on the configuration and numerical scheme can be found in Podvin & Fraigneau (2014). Periodic boundary conditions are enforced in the streamwise () and the spanwise direction (), and the channel lengths along these directions are and , where is the channel half-height. The total domain is discretized using a grid of size , with a uniform grid for the horizontal directions and and a hyperbolic tangent-distributed mesh for the wall-normal direction . The velocity components in the streamwise, wall-normal and spanwise direction will be denoted respectively by or equivalently by . The simulation is conducted at , based on the friction velocity and channel half-height , which corresponds to the simulation of Moser et al. (1999). The database consists of 16500 flow realizations, separated by . The time spanned by the database represents about viscous time units.

In most of the paper, we focus on the domain , but we also considered the full boundary layer where . Due to the strong decrease of the energy spectra with the streamwise wavenumber, and given the typical extent of longitudinal streaks in the wall layer (about 600-1000 wall units (Jiménez, 2013)), the streamwise extent of the domain was limited to representing 900 wall units. The full spanwise extent of the domain (about 1800 wall units) was considered.

The methodology is summarized in Figure 1. Since POD is also applied in time, the autocorrelation tensor needs to be computed from several independent realizations of the same experiment. An ideal sample would consist of the Fourier transform in space and time of the velocity field corresponding to different databases obtained at the same Reynolds number. In practice, owing to the cost of the simulation, we split a single database into several contiguous chunks of length , each of which constitutes a sample. The samples are therefore not independent realizations, but we assume that the period is sufficiently larger than the characteristic time scales of the flow, or at least sufficiently large to allow separation of the dominant time scales. Yet cannot be too large in order to allow for a reasonable number of samples to be constituted. We note that as in Moin & Moser (1989) and Podvin (2001), the number of samples is doubled by considering spanwise reflections and quadrupled by considering both lower and upper parts of the channel (including spanwise reflections).

Each sample then consists of a block of velocity fields uniformly sampled at a rate . For most of the test cases, we have fixed the value of to be 3000 wall units. Results are reported for four datasets extracted from the single database, which represent different domain sizes and sampling rates. The characteristics of the different configurations are given in Table 1. The sampling rate and number of samples were varied in in the configurations D1-D3. These datasets correspond to the wall layer, . The fourth one (D4) corresponds to the full boundary-layer height, .

0.14 100 60 29.5 2950
D2 0.14 100 120 29.5 2950
D3 0.14 500 132 5.9 2950
D4 1 500 132 5.9 2950
Table 1: Description of different test cases; denotes the number of snapshots in a given block, denotes the number of samples, is the time gap between successive snapshots, and is the total time period for a given block (expressed in viscous units).

2.2 Full spatio-temporal POD

Owing to the homogeneity of the statistics in the spatial directions and and their stationarity, POD modes are Fourier modes in the horizontal directions as well as in time. Here, for each configuration, the Fourier transform of the velocity is computed in the temporal and the homogeneous spatial directions for each sample corresponding to a set of fields. The th component of the velocity field () in the physical space is then written as


where denotes the th velocity component in the Fourier space corresponding to the streamwise wavenumber , spanwise wavenumber and frequency . Proper Orthogonal Decomposition is then applied in the wall-normal direction for each triad :


where, for any -uplet ,


Here, represents the usual ensemble average over the space of all possible flow realizations, is the Kronecker-delta function, denotes a specific POD mode for each triad, and * denotes complex conjugation. Note that is independent of as the three velocity components are stacked into a single vector. The complex, stochastic coefficients are uncorrelated and their variance is equal to .

The eigenfunctions and eigenvalues are obtained by solving the following eigenproblem, where the autocorrelation is estimated from taking an ensemble average:


Here, represents the single velocity vector with the three components stacked in it. Upon discretization in , represents the single eigenvector containing the three components. By construction, the eigenvectors are orthogonal with respect to the Euclidean inner product, and can be interpreted as the energy content in each mode. Similarly, for each triad , the modes (i.e., different ’s) are sorted by energy. Since the left-hand-side of Equation (5) is a Hilbert-Schmidt integral operator, its eigenvalues are real and non-negative. The entire procedure is schematically represented in Figure 1. The eigenvalues associated to the triplets can then be gathered and sorted globally according to their magnitude. We will call mode number the global index associated with the sorted modes over spatial wavenumbers , frequency and quantum number , with lower mode number indicating a larger eigenvalue. In all that follows, we focus only on the most energetic modes.

Figure 1: Schematic of the methodology. Each grey slab on the left represents an instantaneous realization of the flow in the physical domain (the wall-normal direction corresponding to the thickness of the slab), while each slab on the right corresponds to the Fourier transform of a block of realizations in time as well as in the horizontal directions and . For each Fourier triad the autocorrelation tensor is computed by averaging over the blocks and an independent eigenproblem is solved.

3 Results from spatio-temporal POD

3.1 POD eigenvalues

Figure 1(a) shows the top 5000 eigenvalues for each configuration. We note that this represents a tiny fraction of the total number of eigenvalues defined in the spatio-temporal space, which is , where is the number of instantaneous fields contained in a sample, and are respectively the numbers of streamwise and spanwise wavenumbers, and is the number of grid points in the wall-normal direction. For D3 this corresponds to about degrees of freedom. The largest eigenvalue where and corresponds to that of the mean mode (discussed in the next section). In the rest of the paper we will focus on the eigenfunctions. We note that the eigenfunctions corresponding to capture about of the turbulent kinetic energy, which is in agreement with Moin & Moser (1989)’s results.

Figure 2: a) Log-log plot of the top 5000 eigenvalues for D1-D4 configurations. An asymptotic fit for the decay rate of the spectra for D3 is also shown. b) Percentage of turbulent kinetic energy (TKE) captured by the most energetic modes for D3.

Figure 2a shows that the spectra in D1 and D2, which correspond to a different number of samples, are essentially indistinguishable from each other. This indicates that the number of samples of ensemble averaging appears to be sufficient for the convergence of the dominant eigenvalues. The spectra in D1 and D2 are similar to that in D3 at low values of , but have higher energy levels at large values of , which is likely to be due to aliasing effects as the fields are sampled at higher rates there than in D3. The spectrum has an asymptotic decay of . Unlike the other configurations, D4 corresponds to the full boundary layer . As expected, energy levels are higher but the shape of the spectrum and its asymptotic decay rate are the same as in the wall region , which suggests self-similarity (see next section). Figure 2b shows the fraction of turbulent kinetic energy (TKE) captured by the first modes for each (as we are considering fluctuations, the contribution from the mean mode was set to zero). The fraction of TKE captured by the most energetic mode numbers is relatively high, but increases only slightly for mode numbers higher than , which highlights the complexity of the flow. This is confirmed by Table 2, which shows that the first eigenvalues capture about of the TKE which is increased only to when modes are included, for the D3 configuration. The trend remains the same for the D4 configuration, where eigenvalues capture around of TKE which increases only to when eigenvalues are included.

Number of eigenvalues () % of TKE (D3) % of TKE (D4)
Table 2: Percentage of TKE in the first eigenvalues for D3 and D4 configurations.
Figure 3: a) Energy spectrum () as a function of ; the data is compared with Moser et al. (1999)’s energy spectrum at (with appropriate rescaling). b) Energy spectrum () as a function of for (D1-D4).

Figure 3 shows the dependence of the integrated spectra with respect to the frequency and spanwise wavenumber for the different configurations. For the D3 configuration, the spatial wavenumbers and the frequency , expressed in wall units, are related to the mode indices as


Figure 3a shows how the sum of the eigenvalues over streamwise wavenumber and frequencies () varies as a function of . Results show a good agreement with the rescaled standard energy spectrum obtained by Moser et al. (1999) at . As expected, the spanwise spectrum over the full layer has more energy in the lower wavenumbers and less energy at higher wavenumbers compared to that in the wall layer.

Figure 3b shows the sum of the dominant eigenvalues over all modes (i.e., ) as a function of frequency for the different configurations D1-D4. As observed previously, aliasing effects can be observed for the fields sampled in time at a lower rate (D1-D2), but the trends are very similar. A marked increase in the energy is observed in the frequency range of , which corresponds to time scales of 200-300 viscous units. This value is in good agreement with the duration of the regeneration cycle identified by Hamilton et al. (1995) and Jimenez & Moin (1991) in minimal domains, as well as the investigations of larger-scale domains reported in Jiménez (2013). A similar frequency peak is observed in the full boundary layer (D4), but with a shift towards slightly higher frequencies: the maximum is located at , while it is located at , in the region .

The predominance of the characteristic frequency at is confirmed by Table 3, which shows the top 30 eigenvalues (denoted by ) with the corresponding Fourier mode indices in space and time and quantum mode for the configuration D3. Although a large fraction (18 modes) of the first 30 modes are associated to , which corresponds to long time scales that are outside the scope of our analysis, all the other modes are characterized by frequencies between and , which correspond to the previously identified time scale of about 200-300 viscous units. The most energetic of these modes corresponds to , i.e., . We note that all the modes in the table are characterized by , which corresponds to the largest eigenvalue for a triad , and , which corresponds to the streamwise-averaged flow.

1 0 1 24.97645
2 0 1 0.008133
3 0 1 0.004205
4 0 1 0.003754
5 0 1 0.001489
6 0 1 0.001277
7 0 1 0.001222
8 0 1 0.001128
9 0 1 0.000923
10 0 1 0.000655
11 0 1 0.000591
12 0 1 0.000552
13 0 1 0.000523
14 0 1 0.000505
15 0 1 0.000457
16 0 1 0.0004568
17 0 1 0.0004559
18 0 1 0.000451
19 0 1 0.000415
20 0 1 0.000412
21 0 1 0.000396
22 0 1 0.0003806
23 0 1 0.0003643
24 0 1 0.000349
25 0 1 0.0003475
26 0 1 0.0003445
27 0 1 0.000335
28 0 1 0.000314
29 0 1 0.000295
30 0 1 0.000289
Table 3: Top 30 most energetic eigenvalues along with their corresponding .
Figure 4: Horizontally-averaged ( and ) temporal Fourier transform of the streamwise velocity component as a function of in a log-log scale at a) and b) . A vertical line is shown at to indicate the peak in both figures.

To provide a comparison with standard analysis tools, Figures 4a and b show the temporal Fourier transform of the streamwise velocity component averaged along horizontal directions for two different heights: and . Both spectra show, although not distinctly, a peak around the characteristic frequency . This illustrates the usefulness of the new POD implementation.

Figure 5: Eigenvalue distribution of the top 5000 eigenvalues for D3 configuration as a function of a) and b) .
Figure 6: a) Energy spectrum summed over frequency , , for different streamwise wavenumbers as a function of . b) Energy spectrum summed over spanwise wavenumber , , for different streamwise wavenumbers as a function of .

Figures 5a and 5b respectively show the distribution of the top 5000 eigenvalues along and for the D3 configuration. Each eigenvalue is represented by a dot. Both figures show that the peaks observed in frequency and wavenumber space are not created by an accumulation of small eigenvalues, but correspond to coherent, more energetic structures. We note that Figure 5b also shows smaller, but noticeable peaks of at around 0.009 and 0.0135, which correspond to the harmonics of the frequency .

Figure 5(a) shows that the energy spectrum integrated in frequency space slowly decreases in an apparently self-similar manner with respect to both streamwise and spanwise wavenumbers. Figure 5(b) shows the corresponding variations in the frequency space of the energy spectrum integrated in the spanwise wavenumber space for different streamwise wavenumbers . For the mode, owing to spanwise reflection invariance, the plot is symmetric in the frequency space, and the characteristic frequency and its harmonics can be clearly identified. These peaks observed at are still present in the spectra, which confirms that the characteristic frequency is not an artefact of the streamwise average. However non-zero streamwise wavenumbers are also characterized by a broad peak, which represents convection effects. In general, defining a convection velocity for the turbulent wall layer is not straightforward (Krogstad et al. (1998); Alamo & Jimenez (2009)), but the question can be more easily addressed in the present framework where each POD eigenfunction is naturally associated with a phase velocity . The spectra can be used to define a global convection velocity at each streamwise wavenumber using




Peaks in the spectra are located at for respective streamwise wavenumbers of . This corresponds to a global convection velocity of 12, with a slight upward shift observed with increasing streamwise wavenumbers. This agrees well with the classical results (Kreplin & Eckelmann (1979), Krogstad et al. (1998), Wallace (2014), Alamo & Jimenez (2009)).

3.2 POD eigenfunctions

Figure 7: Comparison of the mean streamwise profile reconstructed using the first eigenfunction for the mode pair with that of Moser et al. (1999), and with the mean profile obtained from the DNS.
(a) streamwise component
(b) wall-normal component
(c) spanwise component
Figure 8: Comparison of the eigenfunctions corresponding to , , for different configurations.

Figure 7 shows that the reconstructed mean streamwise profile coincides with previous results (Moser et al. (1999)) and with the dominant mode of the decomposition. Figure 8 shows the absolute value of each component of the most energetic mode , for , obtained for the datasets D1-D3. As expected, there is a good agreement between the different configurations, which indicates convergence of the procedure at least for the most energetic modes. In all that follows only results for D3 and D4 will be presented.

In order to gain more insight on the structure of the modes, Figures 9a, c and e show the shape of the most energetic eigenfunction components associated with the characteristic frequency for different spanwise wavenumbers. The location of the velocity maxima moves closer to the wall as the spanwise wavenumber increases, which is in agreement with previous descriptions of wall-attached structures (Alamo et al. (2006), Podvin et al. (2010)). The values of the maxima are similar for both cross-stream components, also in agreement with previous observations (Alamo et al. (2006), Podvin & Fraigneau (2017)), but the value tends to increase with the spanwise wavenumber for the streamwise and the wall-normal component, while it slightly decreases for the spanwise component. The monotonous evolution of the shape of the modes with the spanwise wavenumber suggests self-similarity, as was proposed in Podvin & Fraigneau (2017).

Figures 9b, d and f compare the eigenfunctions obtained on the domain (for D3) and (for D4) for several spanwise wavenumbers. The eigenfunctions of D4 have been rescaled to have the same energy content as those of D3 on . We observe that the eigenfunctions nearly coincide over their common definition domain, which is not a trivial result. The persistence of the eigenfunction shape with respect to the wall-normal extension of the decomposition domain shows the coherence of the most energetic motions over the entire height of the boundary layer. It also means that the restrictions of the eigenfunctions on the larger domain are orthogonal to each other on the smaller domain (since they coincide with the eigenfunctions there), which makes it possible to recover the amplitude of the eigenfunction on the larger domain directly from the projection of the velocity field in the smaller domain onto the corresponding eigenfunction. This shows the relevance of the decomposition for estimation purposes in a context of partial information (see for instance Podvin et al. (2010)).

(a) streamwise component
(b) streamwise component
(c) wall-normal component
(d) wall-normal component
(e) spanwise component
(f) spanwise component
Figure 9: Dominant eigenfunctions corresponding to and for increasing values of . Top row: Streamwise component; Middle row: Wall-normal component; Bottom Row: Spanwise component; Left column: Comparison for (D3) and different spanwise wavenumbers ( to which corresponds to between to ); Right column: Comparison between the domains (D3) and (D4) for selected wavenumbers .

Figure 10 shows the flow generated by selected modes in figure  9 associated with for two different values of . As these dominant modes have no streamwise variation (i.e., ), the three-dimensional shapes of the modes are represented by a contour plot of the streamwise velocity component () and a velocity plot of the in-plane quantities . The modes are characterized by updrafts of low-speed fluid alternating with downdrafts of high-speed fluid, associated with vortical motions in the streamwise direction, in agreement with classic observations (Corino & Brodkey, 1969). The wall-normal extension of the structure decreases as the spanwise wavenumber increases, in agreement with Townsend’s model (Townsend, 1976) of wall-attached eddies (note the similarity of Figure 10 with Figure 9 in Jiménez (2013)). Figure 11 shows that the most energetic streamwise mode convected at a velocity of also corresponds to high and low-speed streaks alternating in the streamwise as well as in the spanwise direction and associated with vortical motions in the cross-stream plane. The wall-normal position of the vortex centers depends on the vertical extension of the streaks.



Figure 10: Contour plot of streamwise velocity component reconstructed from eigenfunction in plane along with a vector plot for velocity components for the non-trivial dominant eigenmode corresponding to a) , , , and b) , , , .
Figure 11: Reconstruction of the most energetic mode corresponding to the peak in the frequency spectrum in figure  5(b) for , , , . Isosurfaces (in arbitrary units) of streamwise component in are shown along with two slice views at and showing a contour plot of the streamwise component and a vector plot of the spanwise and wall-normal components.

The contribution of the most energetic modes to the fluctuating stresses can be evaluated using


where the subscript ’’ is used to distinguish it as a reconstructed quantity using most energetic modes. The bar indicates average over space and time and denotes the usual ensemble operator. Note that the mean modes corresponding to have been excluded from the reconstruction as the stresses are associated with fluctuating components. Figures 12a, b and c show the three components of the root-mean square velocity (rms) and Figure 12d shows the Reynolds stress as a function of reconstructed using 200 and then 3000 modes. The plots are compared with the present simulation results and the DNS results of Moser et al. (1999). A non-negligible fraction of the turbulent intensity is recovered with the most energetic 200 modes: 60% for the streamwise velocity, 30% for the cross-stream components, and about 25% of the Reynolds stresses. The gain obtained using 15 times as many modes (3000) is relatively small (40% increase for the spanwise component and 30% for the wall-normal component), which highlights the complexity of the flow.

Figure 12: Reconstruction of the turbulent intensities and Reynolds stress using 200 and 3000 modes. Comparison with Moser et al. (1999)’s and our DNS data is also shown.

4 Energetics of modes

With the new variant of Proper Orthogonal Decomposition, projection of the Navier-Stokes equations onto the basis of eigenfunctions no longer yields a dynamical system for the temporal amplitudes of the spatial eigenfunctions, but makes it possible to evaluate the contributions made by the different terms of the momentum equations for each mode .

For consistency with standard analysis, we decompose the th component of the velocity field into its mean part which, as we have seen above, can be represented by the mode and a fluctuating part .

We rewrite the equations as


where (mean streamwise velocity) and if and only if and .

To obtain a budget for each mode we project Equation (10) onto the corresponding mode, i.e., we take the Fourier transform in the space, take the inner product in the wall-normal direction, and apply an ensemble average. This gives


where denotes complex conjugation.

  • The first term yields a contribution . We note that this contribution is purely imaginary.

  • The second and third terms represent the interaction of the mode with the mean velocity profile. The second term is associated with the classic production term and is expected to represent a source of energy for the fluctuations:


    The third term represents the convection effect of the mean field:

  • The last term on the left-hand-side of the equation corresponds to the viscous diffusion term defined as:


    Integrating the last term by parts, one has


    Except for the last term corresponding to a boundary effect, all contributions to are real and negative, which corresponds to an energy loss for the mode, as expected.

    Note that and only require information about the mode (one can think of them as purely diagonal operators) and can be directly evaluated from the Proper Orthogonal Decomposition, while the other two contributions require additional information about the coefficients .

  • As in the classic derivation (Aubry et al. (1988)) , the contribution from the pressure term represents the influence of the pressure at the upper boundary of the wall layer and can be expressed as

    where is the Fourier transform (along , and ) of pressure at height , which needs to be evaluated from the DNS. This term is proportional to the wall-normal intensity of the structure at the top of the layer and depends on the velocity-pressure correlation at the top of the layer. It represents an external forcing term which corresponds to the interaction of the wall layer with the outer region (Aubry et al., 1988).

  • Finally, the quadratic interactions can be evaluated as


    requires information about the triple correlations of the coefficients , and involves first derivatives of the eigenfunctions. This term characterizes the energy transfer from the different modes to the mode . It can be seen as a non-isotropic extension of the energy transfer function defined in isotropic turbulence (Zhou, 1993), and involves triads of modes in the space. It also corresponds to the forcing term of the resolvent analysis (McKeon, 2017).

    If all other terms are known, the quadratic interaction terms can be determined from the budget equation:


    The real part of the equation represents a balance between the production and dissipation term , which depend only on the characteristics of the mode , and the interaction terms due to pressure and convection by velocity fluctuations , which characterize how the mode interacts with the full flow. The imaginary part of the equation can be seen as a phase dispersion relation linking the frequency of the mode () with the different physical mechanisms.

Figure 13 represents the contributions of the different terms to the equations for the largest modes (the contributions of the modes corresponding to the mean flow were set to zero). These modes are all characterized by a zero streamwise wavenumber () and a quantum number . All terms are evaluated directly, except the quadratic term which is evaluated using Equation (17). We have checked (not shown here though) that the quadratic terms could not be evaluated correctly by direct computation limited to the first 200 modes, as higher-order contributions were significant, which is typical of the closure problem.

For the sake of clarity, the left plots show mode number from and the right plots show from . The top row (Figures 13a and b) shows that the real part of the production term is essentially balanced by the sum of the real part of the dissipation and the quadratic terms. For the less energetic modes (Figure  13b) we note that the quadratic terms are mostly negative, which is consistent with the idea of a positive energy transfer from the large scales (most energetic modes) to the small scales. Figure 13 shows that the pressure interaction term is very small compared to the other terms, in agreement with Aubry et al. (1988)’s derivation where it is modelled as a stochastic term of small amplitude. If we neglect the influence of the pressure term we have

(a) Real part, modes 1-30
(b) Real part, modes 30-200
(c) Imaginary part, modes 1-30
(d) Imaginary part, modes 30-200
Figure 13: Energy budget for the first 200 most energetic spatially fluctuating modes.

Figures 13c and d show that for the most energetic modes, as , the frequency of the mode is directly related to the imaginary part (denoted by ) of the nonlinear contributions to the mode :


This corresponds to the following closure approximation:


As observed above, we can see that triple correlations between the different POD modes are essential to characterize the energy transfer between the scales. Computing these correlations is difficult, however if can be determined from Equation (17), the decomposition offers a new way to reconstruct the quadratic terms using


Substitution of the velocity decomposition into the left-hand side of the equation requires performing a cumbersome convolution on all coefficients , but the right-hand-side provides a straightforward expression of the quadratic terms in the stochastic POD basis of coefficients , where coordinates are solely determined from second-order statistics. At each scale level , the effect of the quadratic terms is to stretch and to rotate the corresponding velocity mode by a factor . Direct modelling of the distribution of could therefore lead to new formulations of turbulence models in the decomposition framework.

5 Conclusion

Spatio-temporal Proper Orthogonal Decomposition has been applied to the wall layer of a turbulent channel flow. The decomposition represents an efficient data reduction technique which is adapted to large simulation databases. It brings to light typical features of wall turbulence in a straighforward manner, but also provides a fresh viewpoint on the flow organization. Due to symmetry properties, the decomposition singles out empirical eigenfunctions for each frequency and horizontal spatial wavenumber. Besides time scales superior to 3000 wall units, which our limited implementation of POD did not allow us to characterize fully, we have shown that the most energetic modes were characterized by a time scale on the order of 250 wall units. About 30% of the turbulent kinetic energy was captured by the 200 most energetic modes. The most energetic modes were found to have a self-similar shape that appeared largely independent from the wall-normal extent of the decomposition domain, which shows the coherence of the motions over the height of the boundary layer. Convection velocities could be directly defined from the POD spectrum. A global convection velocity on the order of 12 was identified in the wall layer, in good agreement with previous approaches.

Finally, substitution of the decomposition into the Navier-Stokes equation and numerical computation of the different contributions of the modes to the turbulent kinetic energy budget highlighted the key role played by quadratic interactions and allowed us to propose a new closure formulation to model the contribution of these interactions. Such relationships, which need to be further explored in a careful manner, could be useful to derive new turbulence models. We hope that this work will pave the way for comprehensive investigations of wall turbulent flows using spatio-temporal Proper Orthogonal Decomposition.


This work was supported by the Center of Data Science from the Paris-Saclay University. Computations were carried out at IDRIS-GENCI (project 02262).


  1. Alamo, J. C. Del & Jimenez, J. 2009 Estimation of turbulent convection velocities and corrections to taylor’s approximation. Journal of Fluid Mechanics 640.
  2. Alamo, J. C. Del, Jimenez, J., Zandonade, P. & Moser, R. D. 2006 Self-similar vortex clusters in the turbulent logarithmic region. Journal of Fluid Mechanics 561.
  3. Aubry, N., Holmes, P., Lumley, J. L. & Stone, E. 1988 The dynamics of coherent structures in the wall region of the wall boundary layer. Journal of Fluid Mechanics 192.
  4. Blackwelder, R. F. & Haritodinis, J. H. 1983 Scaling of the bursting frequency in turbulent boundary layers. Journal of Fluid Mechanics 132.
  5. Corino, E. R. & Brodkey, R. S. 1969 A visual investigation of the wall region in turbulent flow. Journal of Fluid Mechanics 37.
  6. Dennis, David J.C. 2015 Coherent structures in wall-bounded turbulence. Anais da Academia Brasileira de Ciencias 87, 1161 – 1193.
  7. George, W. K. 2017 A 50-year retrospective and the future. In Whither Turbulence and Big Data in the 21st Century?. Springer.
  8. Hamilton, J. M., Kim, J. & Waleffe, F. 1995 Regeneration mechanisms of near-wall turbulence structures. Journal of Fluid Mechanics 287.
  9. Jiménez, J. 2013 Near-wall turbulence. Physics of Fluids 25 (10).
  10. Jimenez, J. & Moin, P. 1991 The minimal flow unit in near-wall turbulence. Journal of Fluid Mechanics 225.
  11. Kim, H. T., Kline, S. J. & Reynolds, W. C. 1971 The production of turbulence near a smooth wall in a turbulent boundary layer. Journal of Fluid Mechanics 50 (1).
  12. Kreplin, H. P. & Eckelmann, H. 1979 Propagation of perturbations in the viscous sublayer and adjacent wall region. Journal of Fluid Mechanics 95.
  13. Krogstad, P. A., Kaspersen, J. H. & Rinestead, S. 1998 Convection velocities in a turbulent boundary layer. Physics of Fluids 10.
  14. Loève, M. 1977 Probability Theory. Berlin, Germany: Springer.
  15. Lumley, J. L. 1967 The structure of inhomogeneous turbulent flows. In Atmospheric Turbulence and Radio Wave Propagation. Nauka, Moscow.
  16. McKeon, B. J. 2017 The engine behind (wall) turbulence: perspectives on scale interactions. Journal of Fluid Mechanics 817.
  17. Moin, P. & Moser, R. D. 1989 Characteristic-eddy decomposition of turbulence in a channel. Journal of Fluid Mechanics 200.
  18. Moser, R. D., Kim, J. & Mansour, N. N. 1999 Direct numerical simulation of turbulent channel flow up to = 590. Physics of Fluids 11 (4).
  19. Podvin, B. 2001 On the adequacy of the ten-dimensional model for the wall layer. Physics of Fluids 13 (1).
  20. Podvin, B. & Fraigneau, Y. 2014 POD-based wall boundary conditions for the numerical simulation of turbulent channel flows. Journal of Turbulence 15 (3).
  21. Podvin, B. & Fraigneau, Y. 2017 A few thoughts on proper orthogonal decomposition in turbulence. Physics of Fluids 29.
  22. Podvin, B., Fraigneau, Y., Jouanguy, J. & Laval, J. P. 2010 On self-similarity in the inner wall layer of a turbulent channel flow. Journal of Fluids Engineering 132 (4).
  23. Podvin, B. & Lumley, J. L. 1998 A low-dimensional approach for the minimal flow unit. Journal of Fluid Mechanics 362.
  24. Robinson, S. K. 1991 Coherent motions in the turbulent boundary layer. Annual Review of Fluid Mechanics 23.
  25. Smits, A. J., McKeon, B. J. & Marusic, I. 2011 High–reynolds number wall turbulence. Annual Review of Fluid Mechanics 43 (1).
  26. Stanislas, M., Perret, L. & Foucaut, J. M. 2008 Vortical structures in the turbulent boundary layer: a possible route to a universal representation. Journal of Fluid Mechanics 602.
  27. Townsend, A. A. 1976 The Structure of turbulent shear flow. Cambridge University Press.
  28. Wallace, J. M. 2014 Space-time correlations in turbulent flow: A review. Theoretical and Applied Mechanics Letters 4.
  29. Zhou, J. 1993 Interacting scales and energy transfer in isotropic turbulence. Tech. Rep. CR-191477. NASA.
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description