
[Sponsors] 
April 21, 2014, 07:23 
Energy spectrumwrong spectrum shape

#1 
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Hi,all
Recently I got stuck with energy spectrum. I sampled the data across the spanwise direction(z) in channel flow, Uz, then I want to plot the energy spectrum, I adopt the method here Energy spectrum calculation as Robert said,here is the code Code:
clear all; [z Ux Uy Uz]=textread('Zleft_U.xy','%f%f%f%f'); N=128; Lz=6.28; kz=(2*pi*(1:128)/Lz)'; kzz=round(kz); Uzk=fft(Uz,N)/N; %Fs=2*pi; %f = kz/2*linspace(0,1,N/2+1); E=0.5*4*((2*abs(Uzk(1:N/2+1))).^2).*pi.*(kzz(1:N/2+1)).^2 plot(kzz(1:N/2+1),E) 

April 21, 2014, 08:47 

#2 
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Usually,the spectrum is averaged in time. How long should the time period be? As soon as the curve becomes smooth?


April 21, 2014, 13:25 

#3 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 1,247
Rep Power: 17 
first, you should plot in a double logarithmic scale, then you must sample the velocity field in a sufficient number of file, the sampling time depending on the turnover characteristic time. At each time compute the Fourier coefficients and then do the time average


April 21, 2014, 20:35 

#4  
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Quote:
Here, I have another question: I have read some papers about this and found that the usual method is to do the spectrum average in time and plane parallel to the wall. I'm quite confused about the average in plane(say the plane y+=10 for example). If I want to plot the spanwise spectrum(z direction in channel flow), is that means I should do average of each z position, with all the x positions(z fixed, and average the values at different x)? 

April 21, 2014, 22:51 

#5  
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Quote:


April 22, 2014, 03:35 

#6 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 1,247
Rep Power: 17 
the averaging along a plane is used only for homogeneous flows (for example channel with periodic conditions) in order to make the statistics more "robust".


April 22, 2014, 04:12 

#7 
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
so If I want to plot the spanwise spectrum(z direction in channel flow), is that means I should do average of each z position, with all the x positions(z fixed, and average the values at different x)? Just as I have asked?


April 22, 2014, 04:31 

#8 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 1,247
Rep Power: 17 
well, assume the channel has: xstreamwise (periodic), ynormal to walls, zspanwise (periodic) directions.
 Fix some y value (or y+ value)  Chose if you want the spectrum along Kx or Kz.  For example, along the streamwise direction the frequencies are Kx=nx*2*pi/Lx (nx=1,...namx), then for each z value of the y=constant plane perform the FFT and store the Fourier coefficients.  Perform the average of the FOurier coefficients over z, this way you have the coefficients depending only on Kx.  Repeat the procedure for all the temporal fields 

April 22, 2014, 06:07 

#9  
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Quote:


April 24, 2014, 06:43 

#10  
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Quote:
Sorry for my late. I tried similar as you say to calculate the twopoint correlation(just to simply see if it's right), however, the correlation exceeds 1 Code:
clear all; [x y z Ux Uy Uz]=textread('U_constantPlane.txt','%f%f%f%f%f%f'); [xx yy zz Uxm Uym Uzm]=textread('UMean_constantPlane.txt','%f%f%f%f%f%f'); [xorder]=textread('x.txt','%f'); Uxmatrix=reshape(Ux,128,256); Uxmmatrix=reshape(Uxm,128,256); ux=UxmatrixUxmmatrix; uxplane=mean(ux); for i=1:256, corrx(i)=uxplane(1)*uxplane(i)/(uxplane(1)*uxplane(1)); end plot(xorder,corrx) 

April 24, 2014, 09:11 

#11 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 1,247
Rep Power: 17 
have you normalized correctly? why don't you try the spectrum?


April 24, 2014, 12:21 

#12 
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
I think It's due to the data structure. My plane is 128*128 grid. While the output contains twice the points value(32768 values). I'll try to figure out a better way to output. I have just tried the spectrum(one temporal data), the result seems incorrect too.
Q1:Is it need to do like this: for i=1:N, fftcoefficient/N as someone suggested ? Q2:Also, I have another question. What's the mean velocity's definition? If we do the average on the plane, is that means the average should be the fff(U)dxdzdt? In the present reply , I use the temporal velocity U and timeaveraged UMean, is this right or just as I said before, should do further average? 

April 26, 2014, 09:08 

#13 
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
After numerous trying, I think I get the right(seems) spectrum.
Code:
clear all; [x y z Ux Uy Uz]=textread('Uzyplus35.txt','%f%f%f%f%f%f'); [xx yy zz Uxm Uym Uzm]=textread('UzMeanyplus35.txt','%f%f%f%f%f%f'); %[xorder]=textread('z.txt','%f');,just ignore this line N=128; Lz=6.28;%length in the z direction kz=2*pi*(N/2:1:N/21)/Lz;%frequencies, or wavenumber Uzre=reshape(Uz,N,N);%change the column of Ux into a matrix,each column corresponding to one z Uzmre=reshape(Uzm,N,N); Uzmreave=mean(Uzmre(:));%field average Uzmfield=zeros(N,N); for i=1:N, for j=1:N, Uzmfield(i,j)=Uzmreave; end end w=UzreUzmfield;%get the velocity fluctuation=instannous velocity  timeaveraged velocity wk=fft(w);%each column corresponds to a z value wkk=fftshift(wk); E=zeros(N,N); for i=1:N, for j=1:N, E(i,j)=(wkk(i,j)*conj(wkk(i,j)));%equal to abs(wkmean(i))^2 end end Emean=mean(E,2); loglog(kz,Emean); The latter one is from a paper with DNS data in channel flow,Re=2666, while my simulation has Re=2900, the deviation is somehow large, do you think it's reasonable or my method is wrong? 

April 26, 2014, 13:15 

#14 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 1,247
Rep Power: 17 
the shape seems correct, however I see one order of magnitude in the lowest wavenumber... check it carefully, I suspect you are doing the streamwise component of the spectrum while the DNS reports the spanwise...


April 27, 2014, 05:30 

#15  
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Quote:
Now, the result is correct 

April 30, 2014, 00:30 

#16 
New Member
amjad
Join Date: May 2012
Posts: 12
Rep Power: 4 
Hello,
could anyone explain the steps of this code please thanks 

June 8, 2014, 20:37 

#17  
New Member
ehk
Join Date: Sep 2012
Posts: 21
Rep Power: 4 
What is the dimension of Energy Spectra calculated in this code? Is it (m^4/s^3)?
As FFT(m/s)=m^2/s. Quote:


June 8, 2014, 23:29 

#18 
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 

June 8, 2014, 23:54 

#19 
New Member
ehk
Join Date: Sep 2012
Posts: 21
Rep Power: 4 
As I mentioned FFT[m/s]=m^2/s. If you apply FFT for nondimensional variable, the result would have the dimension of [m], as there is an integration over space to derive FFT coefficients. If you also check the literature for plots of Energy Spectra, it is seen the dimension is like [m^3/s^2]=FFT(Kinetic Energy)=FFT([m^2/s^2]). The way to nondimensionalize Energy Spectra is dividing by (\eta u_*) where \eta is the kolmogrov length scale.


June 9, 2014, 02:29 

#20  
Senior Member
Huang Xianbei
Join Date: Sep 2013
Location: CAU,China
Posts: 211
Rep Power: 4 
Quote:


Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
ATTENTION! Reliability problems in CFX 5.7  Joseph  CFX  14  April 20, 2010 15:45 
energy in sonicFoam  joern  OpenFOAM Running, Solving & CFD  0  June 16, 2009 15:53 
udf error  srihari  FLUENT  0  February 9, 2009 09:00 
Energy Spectrum in LES Simulation  Anindya  FLUENT  0  March 17, 2006 14:20 
Why FVM for highRe flows?  Zhong Lei  Main CFD Forum  23  May 14, 1999 13:22 