Power spectrum of LES: understanding of the result
 

Good afternoon,

I know that it is a very common doubt, it has been spoken several times... but I couldn't get a good idea about what should I expect from a power spectra density result.

The procedure to obtain my spectrum has been:
- At a given position, I have obtained the velocity value each time step dt, so my sampling frequency is Fs=1/dt
- I have done the same with several points in the spanwise direction, I suppose that the flow is isotropic in that direction
- This information, I give it to a PSD Welch algortihm, who provides me a frequency range, and the PSD for each frequency value.
- I apply an average in the spanwise direction to smooth the PSD

When I plot it, I obtain something like the attached figure.

As far as I know, the PSD it usually used to estimate if the mesh, the numerical setup of the LES... is good enough, but I don't know how to understand this.

If somebody can help me, I would be very gratefull. Many thanks in advance,

Agustin

First of all, I suggest to do a check on your spectrum calculation to see if it works correctly. You can build your signal by sampling the sin(kx) function in [0:2pi] so that once you fix a value for k you have only a peack in your spectrum.

Second, the spectrum you plot can not be discussed if you do not give details of your flow problem and mesh resolution. The shape is very strange.

I agree with Professor Denaron comment. This spectrum does not look like a turbulent kinetic energy spectrum. I would start with a an isotrooic flow field; for example a channel flow. You can choose the center line (far from the walls). You can use the DNS data set from the jhon Hopkins turbulence data base and test your spectrum (code) with the one they provide.
As far as understanding; profesor Davison has suggested in several papers that the energy spectrum computed from the two point correlation is the right approach to validate your les computation. This is also in line with Pope's textbook.
I am looking forward to more contributions from other members and professor Denaro.

Yes. Correlations and spectra are mutually related by a mathematical transformation. This topic can be read in the book of Pope. A breif summary is here https://www.io-warnemuende.de/tl_fil...Chap4_WS08.pdf

However, understanding if a resulting spectra is mathematically well computed does not imply that it is physically meaningful. For this reason I suggest to split the study first in a numerical assesement of the procedure and then in an application on a well controlled flow problem.

Hello, sorry for not answering during the weekend.

As Filippo told me, I have tried to use the PSD Welch with different sinusoidal signals, and I have checked that I obtain a peak of frequency in the specified one. I is, the code of the algorithm is ok. Now I have to see if my input is ok or not. The mistake in my spectrum may be here.

I have my velocity signal, that I can decompose in a mean and fluctuating terms. With the fluctuating ones, I can get the turbulent kinetic energy as q=1/2 * (u'^2 + v'^2 + w'^2).

So, I get the Fourier transform's for u'u', v'v' and w'w', and I get the FT of the turbulent kinetic energy with the previous expression. Is is like this? One question that comes to my mind is if we can get the expression of q along the time and get its FFT. Would we obtain the same results?

Using the sine function you get exactly 1 peak and the other coefficients are zero at machine accuracy?? How about the value of the peak?

why do you perform the FFT of the kinetic energy?? How do you perform the periodic windowing of the time-signal??

I agree with Julio, I suggest to do a controlled case using the channel flow problem and performing the spectral analysis in streamwise and spanwise direction. Only once you can superimpose the spectra with the existing database you can be sure.

Just a final consideration, if you want to characterize a LES spectrum, you must analyze the spatial energy spectra to see the correspondance with the LES filter. The time step characterizes a different Nyquist cut-off that, usually, lies withint the Kolomogorov time-scale.

Thanks for you answer. I attach my sinusoidal signal (20 Hz), with a sampling frequency of Fs=1000 Hz, and the output of the Welch algorithm. The other frequencies, as you can see, specially the higher ones, are very low, in the order of 10E-11 and lower.

About the FFT on the TKE, I was just asking. I give as input the velocity in this algorithm.

I am totally a newbie in all these things in LES, and I want to understand these concepts properly.

Your spectrum is not correct...to see better do not use the continuous line but only symbols, you will see other contribution out of the peak at 20 Hz.
You must se only one non-zero coefficients, the others being O(10^-13).

So, your first task now is clear, focus on the right FFT procedure.

Second, if you want to get into the LES field you must start learning the basis...for example, read the intro of book of Sagaut. You must be aware of the fact the LES has the aim of reducing the degree of freedom and that is obtained by filtering the spatial small flow scales. As a consequence, you must assess your results by looking at the spatial correlation and spectra.
For this reason you have to start with well controlled flow problem.

Hello again,

for sure, when I apply the FFT to the sinus signal I got what you say, however I have some discrepancies close to the peak of 20 Hz.

If I apply directly the FFT to the fluctuating part of the velocity. I obtain complex results, and in order to obtain the term u'u', I multiply by the conjugate. After that, I get the turbulent kinetic energy in the frequency domain as usual, k=0.5*uu+vv+ww. At the end, what I obtain is something like the image I attach.

I forgot to say that my case is a natural convection boundary layer, so, I don't know if the profile I obtain is due to that fact, or there is something bad in my procedure.

Thanks for your replies!

Until you are not able to reproduce the correct sinus response in the PDS, you cannot apply your analysis to a flow problem...
You must see only one spike at the correct frequency and zero (at machine accuracy) everywhere. So, check if you passed correctly the periodic window to the FFT, it seems you get frequencies aliased.

Has your buoyancy-driven flow problem homogenous directions ?

This is a plot of the spectrum for the sin kx function at k=1

Yes. the flow ascends in the Y direction (streamwise), the normal direction to the wall is X, and the spanwise is the Z direction. It is a simple boundary layer, so it should be homogeneous in the Z direction, where I impose periodicity.

As you can see, I think that I am in the right way. In the FFT of sin(2*pi*t), I can see clearly the peak. That peak only contains the real part of the solution. If I plot the magnitude of the complex values, I get the other one. There are several orders of magnitude of difference. Could I consider that it is ok?

No, when you plot the magnitude, that is the modulus of the Fourier coefficient, you must get 0.25. The other coefficients are less than O(10^-14).
Conversely your second plot has a wrong magnitude and other coefficients are aliased.

Mmm... I am checking the values I get, and it seems that in that modulus, I should divide at a certain point by my number of samples, is it? For the moment, I was not using a window, but using the Hann one, I get really close of what you mean, only the closest points have a non negligible value.

Check the correct number of nodes you pass to the FFT and the correct value of the period... I think that is your problem

Hello again,

I have checked my code, and I obtain a peak of 0.5. I think I am missing something... What I do is

My max FFT value is 0.5 at 1Hz, as you can see in the figure. I think I am missing a factor of two in the input of my fft.

EDIT: well, the missing factor I think is the power of 2, otherwhise I am plotting a module, not the square of my signal. I should use fftSignal*conj(fftSignal)

yes, you can compare with the portion of my Matlab script I posted here

https://www.cfd-online.com/Forums/ma...tml#post565863

Buongiorno Filippo, sorry for bothering you again.

I am taking a look into your code, and this part:
In my case, I am getting a time-correlation, and in yours is a spatial one, but the theory behind is the same. I don't know what "l" means, but I suppose is an index related with position. Your sum_u, is the turbulent kinetic energy, so "k" and "l" are related with the x,y,z direction, aren't they?

At the end, you get ek1_uu dividing since you are averaging in space with dz.

The cycle is over "l", an index for the spanwise direction in a channel flow problem. This way I perform a statistical averaging of the one-dimensional spectra along the stream-wise direction at several position.
Supplementary averaging is useful to get statistically meaningful results.
For a time analysis, you should consider the same approach, using several periodical time-windows, performing spectra on each one and then averaging.