Anomaly in energy spectrum

arunsmec · May 30, 2017, 00:59

I am computing energy spectrum for computationally obtained turbulence data, which is statistically homogeneous in 2 directions (X & Z). The energy spectrum is calculated for a homogeneous plane, using Fourier transforms of the spatial data. $\vec{V}(x,z)\rightarrow\hat{V}(k_x,k_z)$ .
The calculated spectrum shows a 'bump' at higher wavenumbers (not the highest). It is prominent in a semilog plot (attached).
I am wondering what could have resulted in this anomalous behaviour. I assume it is associated with some error in calculating the spectrum but has nothing to do with the physics. However, I am unable to find any error in the calculation.
I would like to know if anybody in the community has encountered a similar erroneous spectrum, and what could cause such anomaly?

piu58 · May 30, 2017, 01:22

You may look whether the wavelength of your bump region coincides with the size of the mesh.

arunsmec · May 30, 2017, 01:49

Hi piu58,
I should have mentioned my mesh and domain earlier.

There are 512x128 points in the XZ plane which has a physical dimension $16\pi\times4\pi$ . Hence the maximum wavenumber for both $k_x\&k_z$ is 32. The computed

for various combinations of $k_x\&k_z$ are binned in nearest integer values of $k=\sqrt{k_x^2+k_z^2}$ . As you suspected, the bump corresponds to

.
I understand that under-resolved data results in oscillations at highest wavenumbers.In this case, I am perplexed by the smooth spectrum beyond

.

FMDenaro · May 30, 2017, 04:29

What problem are you solving, a plane channel flow? And what about the Re_tau of your problem? what about y+ of the computed spectrum?

What I do not understand is the fact you have the spectrum extending over the Nyquist frequency... Or am I wrong in understanding your comment?

Then, the magnitude of E(k) is very high, why?

arunsmec · May 30, 2017, 05:22

The problem under consideration is temporal mixing layer.
On the x-axis I have $k=\sqrt{k_x^2+k_z^2}$ , which goes up till $32\sqrt{2}$ . That is why the spectrum extends beyond the Nyquist limit 32. However, I am unsure about the correctness of the same.
As of now, both k and E(k) are not scaled.
My objective of plotting the spectrum is to make sure I am resolving the scales sufficiently.

FMDenaro · May 30, 2017, 05:42

Quote:

Originally Posted by arunsmec

The problem under consideration is temporal mixing layer.
On the x-axis I have $k=\sqrt{k_x^2+k_z^2}$ , which goes up till $32\sqrt{2}$ . That is why the spectrum extends beyond the Nyquist limit 32. However, I am unsure about the correctness of the same.
As of now, both k and E(k) are not scaled.
My objective of plotting the spectrum is to make sure I am resolving the scales sufficiently.

No, I strongly suggest to do the analysis along x and z direction, separately. The spectra must be computed only up to the Nyquist frequency along each direction. Furtermore, have you performed also a statistical averaging (plane and time)? Have you waited enough time to let the flow develop?

arunsmec · May 30, 2017, 06:03

Quote:

Originally Posted by FMDenaro

Furtermore, have you performed also a statistical averaging (plane and time)? Have you waited enough time to let the flow develop?

The spectrum corresponds to fluctuations from a plane averaged mean, during the self-similar evolution of the mixing layer. I have not performed a time averaging. It is observed that the shape of the curve remains the same at all instants.

Quote:

Originally Posted by FMDenaro

I strongly suggest to do the analysis along x and z direction, separately. The spectra must be computed only up to the Nyquist frequency along each direction.

Does this mean calculating one-dimensional spectrum?

To summarise, can I attribute the "bump" to the method where I am going above the Nyquist limit?

FMDenaro · May 30, 2017, 06:31

Compute the 1d spectra but perform the spatial averaging along the normal direction

arunsmec · May 31, 2017, 01:24

Quote:

Originally Posted by FMDenaro

Compute the 1d spectra but perform the spatial averaging along the normal direction

I infer this as average along Z for

, and vice versa. In one of the threads in this forum, it is suggested that a Fourier transform is performed first and the coefficients are averaged. Is that the right way to do it? Or should I take average before transformation?

FMDenaro · May 31, 2017, 04:16

Quote:

Originally Posted by arunsmec

I infer this as average along Z for

, and vice versa. In one of the threads in this forum, it is suggested that a Fourier transform is performed first and the coefficients are averaged. Is that the right way to do it? Or should I take average before transformation?

No, you compute first the FFT and then you perform the average of the coefficients

arunsmec · June 1, 2017, 04:07

The calculated one-dimensional spectrum is attached.

A smooth curve is not obtained even after averaging FFT coefficients in the Z direction. Is it required to perform an average in time to obtain a smooth spectrum?

FMDenaro · June 1, 2017, 04:31

That looks quite correct ... have a look at sec.5.2 here
https://www.researchgate.net/publica...-uniform_grids

and Ref.42 in the reference

arunsmec · June 3, 2017, 07:18

One dimensional spectra, $E(k_x)\ \& \ E(k_z)$ , calculated from the spatial data corresponding to center plane, are attached. Here I have binned them into integer values of $k_x\& k_z$ .

Quote:

Originally Posted by FMDenaro

That looks quite correct ... have a look at sec.5.2 here

In this reference, the spectrum is calculated as
$E(k)=\frac{1}{2L}\int^{+L}_{-L}\left|\hat{u}(k,y)\right|^2dy$ ,
which seems to be an averaging in the non-homogeneous direction. If so, what difference this makes to the spectrum?

FMDenaro · June 3, 2017, 07:36

Your results seems physically acceptable.

You can have two choices:
1) the averaging along the non-homogeneous domain give the total energy content of the mixing layer for each wavenumber component.
2) computing and analysing spectra at several stations along y+. This is for example the standard procedure for a channel flow.

Both are correct and the choice depends on what you want to analyse. I used the averaging along y to compare the solutions with the paper of Lesieur.

arunsmec · June 3, 2017, 07:57

Quote:

Originally Posted by FMDenaro

Your results seems physically acceptable.

Dear FMDenaro,

Thank you for your valuable suggestions.

Can you provide some clarification regarding the 'binning' procedure?

Quote:

Originally Posted by arunsmec

The calculated one-dimensional spectrum is attached.
A smooth curve is not obtained even after averaging FFT coefficients in the Z direction. Is it required to perform an average in time to obtain a smooth spectrum?

In this plot, I used values of k, scaled according to the length of the domain and therefore are rational numbers.

Quote:

Originally Posted by arunsmec

One dimensional spectra, $E(k_x)\ \& \ E(k_z)$ , calculated from the spatial data corresponding to center plane, are attached.

In the latter plot, the spectrum is binned to integer values of k.
$E(k)=\sum_{k'} E(k'),\ \forall k'\in [k-0.5, k+0.5)$

FMDenaro · June 3, 2017, 08:12

The plot is versus the frequency or the integer wavenumber, the relation being k=n*2*pi/L (n=0,1,2..,Nmax). The Nyquist frequency is kmax=pi/h. This is for each direction Lx and Lz.
What do you mean exactly for "binning procedure"?

arunsmec · June 3, 2017, 08:27

The length of the domain in x direction is $16\pi$ and has 512 points. A Fourier transform of this spatial data gives amplitudes of 256 frequencies (n=0,1...256). Corresponding wavenumbers are k=0,1/8,2/8.....,32 and they are not all integers.

In the latter plot, I have used only integer values of wavenumber on the x-axis. The energy of non-integer wavenumbers is rounded off to nearest integer k. For example, E(1) is the sum of all u(k)^2 for k=4/8,5/8,.....,10/8,11/8.
This procedure was described in a reference, but I couldn't trace it now.

FMDenaro · June 3, 2017, 08:38

After you perform the FFT, you have a vector of the complex Fourier coefficients for each integer n=0,1,... Thus, you can plot (in a log scale) the modulus versus these non-dimensional numbers or you simply convert the integers in the dimensional spatial frequency k = n*2*pi/L.
No "binning" in such procedure...

May 30, 2017, 01:22		#2
piu58 Senior Member Uwe Pilz Join Date: Feb 2017 Location: Leipzig, Germany Posts: 744 Rep Power: 15	You may look whether the wavelength of your bump region coincides with the size of the mesh. __________________ Uwe Pilz -- Die der Hauptbewegung überlagerte Schwankungsbewegung ist in ihren Einzelheiten so hoffnungslos kompliziert, daß ihre theoretische Berechnung aussichtslos erscheint. (Hermann Schlichting, 1950)

June 3, 2017, 08:38		#18
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,896 Rep Power: 73	After you perform the FFT, you have a vector of the complex Fourier coefficients for each integer n=0,1,... Thus, you can plot (in a log scale) the modulus versus these non-dimensional numbers or you simply convert the integers in the dimensional spatial frequency k = n2pi/L. No "binning" in such procedure... arunsmec likes this.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
How can I map the spectrum energy with its associated velocity field?	juliom	Main CFD Forum	0	March 23, 2016 09:16
ATTENTION! Reliability problems in CFX 5.7	Joseph	CFX	14	April 20, 2010 16:45
turbulent energy spectrum	cfd.newbie	Siemens	1	June 20, 2008 00:48
LES correlation and turbulent energy spectrum	Fabian	Main CFD Forum	4	October 18, 2005 03:04
Energy Spectrum	Emad Khalifa	Main CFD Forum	3	June 30, 2003 17:03

May 30, 2017, 01:49		#3
arunsmec Member Aru Join Date: Feb 2012 Location: Chennai, India Posts: 40 Rep Power: 14	Hi piu58, I should have mentioned my mesh and domain earlier. There are 512x128 points in the XZ plane which has a physical dimension $16\pi\times4\pi$ . Hence the maximum wavenumber for both $k_x\&k_z$ is 32. The computed for various combinations of $k_x\&k_z$ are binned in nearest integer values of $k=\sqrt{k_x^2+k_z^2}$ . As you suspected, the bump corresponds to . I understand that under-resolved data results in oscillations at highest wavenumbers.In this case, I am perplexed by the smooth spectrum beyond .

May 30, 2017, 04:29		#4
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,896 Rep Power: 73	What problem are you solving, a plane channel flow? And what about the Re_tau of your problem? what about y+ of the computed spectrum? What I do not understand is the fact you have the spectrum extending over the Nyquist frequency... Or am I wrong in understanding your comment? Then, the magnitude of E(k) is very high, why?

May 30, 2017, 05:22		#5
arunsmec Member Aru Join Date: Feb 2012 Location: Chennai, India Posts: 40 Rep Power: 14	The problem under consideration is temporal mixing layer. On the x-axis I have $k=\sqrt{k_x^2+k_z^2}$ , which goes up till $32\sqrt{2}$ . That is why the spectrum extends beyond the Nyquist limit 32. However, I am unsure about the correctness of the same. As of now, both k and E(k) are not scaled. My objective of plotting the spectrum is to make sure I am resolving the scales sufficiently.

May 30, 2017, 06:31		#8
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,896 Rep Power: 73	Compute the 1d spectra but perform the spatial averaging along the normal direction

June 1, 2017, 04:31		#12
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,896 Rep Power: 73	That looks quite correct ... have a look at sec.5.2 here https://www.researchgate.net/publica...-uniform_grids and Ref.42 in the reference

June 3, 2017, 07:36		#14
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,896 Rep Power: 73	Your results seems physically acceptable. You can have two choices: 1) the averaging along the non-homogeneous domain give the total energy content of the mixing layer for each wavenumber component. 2) computing and analysing spectra at several stations along y+. This is for example the standard procedure for a channel flow. Both are correct and the choice depends on what you want to analyse. I used the averaging along y to compare the solutions with the paper of Lesieur.

June 3, 2017, 08:12		#16
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,896 Rep Power: 73	The plot is versus the frequency or the integer wavenumber, the relation being k=n2pi/L (n=0,1,2..,Nmax). The Nyquist frequency is kmax=pi/h. This is for each direction Lx and Lz. What do you mean exactly for "binning procedure"?

June 3, 2017, 08:27		#17
arunsmec Member Aru Join Date: Feb 2012 Location: Chennai, India Posts: 40 Rep Power: 14	The length of the domain in x direction is $16\pi$ and has 512 points. A Fourier transform of this spatial data gives amplitudes of 256 frequencies (n=0,1...256). Corresponding wavenumbers are k=0,1/8,2/8.....,32 and they are not all integers. In the latter plot, I have used only integer values of wavenumber on the x-axis. The energy of non-integer wavenumbers is rounded off to nearest integer k. For example, E(1) is the sum of all u(k)^2 for k=4/8,5/8,.....,10/8,11/8. This procedure was described in a reference, but I couldn't trace it now.