CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Main CFD Forum

Energy Spectrum !!!

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree5Likes

Reply
 
LinkBack Thread Tools Display Modes
Old   April 6, 2009, 04:16
Exclamation Energy Spectrum !!!
  #1
New Member
 
Yilei Shi
Join Date: Apr 2009
Posts: 5
Rep Power: 7
stanking is on a distinguished road
I have a question about turbulent energy spectrum !!! Now I have the data of position and velocity at grid.

the data is like this:
X Y Vx Vy

I have already used FFT to transform Vx and Vy to V_kx and V_ky.

Then I calculate the kx = 2*pi*N / L and kx = 2*pi*N / L. (N = 1,2,....)

Here is a question. N is from 0 or from 1???

Now the data is like this:
kx ky V_kx V_ky

How can I calculate the energy spectrum??

1. How can I calculate the k? k = sqrt(kx^2+ky^2), is right?

2. V(k) = sqrt(V_kx^2 + V_ky^2), is right??

3. plot(k,V(K),k^(-5/3))

please tell me the detail of the calculation procedure!!!

thanks!!!

Last edited by stanking; April 6, 2009 at 06:21.
stanking is offline   Reply With Quote

Old   April 6, 2009, 05:48
Default
  #2
Senior Member
 
lakeat's Avatar
 
Daniel WEI (老魏)
Join Date: Mar 2009
Location: South Bend, IN, USA
Posts: 688
Blog Entries: 9
Rep Power: 10
lakeat is on a distinguished road
Send a message via Skype™ to lakeat
I remember I have asked the same question, and there are indeed some good detailed answer to them, so please search the forum and see if you can get them.

Daniel
__________________
~
Daniel WEI
-------------
NatHaz Modeling Laboratory
Department of Civil & Environmental Engineering & Earth Sciences
University of Notre Dame, USA
Email || My Personal CFD Blog
lakeat is offline   Reply With Quote

Old   April 6, 2009, 06:04
Thumbs up
  #3
New Member
 
Yilei Shi
Join Date: Apr 2009
Posts: 5
Rep Power: 7
stanking is on a distinguished road
Sorry, I have searched all of the forum, but didn't find the article of your question about turbulent energy spectrum calculation. So can you give me a link? Or please explain the last 3 steps!

Thanks!
stanking is offline   Reply With Quote

Old   April 6, 2009, 09:45
Default
  #4
Member
 
Newton KF
Join Date: Mar 2009
Posts: 36
Rep Power: 7
NewtonKF is on a distinguished road
The book from Pope (turbulent flows) have this kind of information... It is very complex, but has all energy spectrum detais... including the equations for the two point correlation...
Good luck...
NewtonKF is offline   Reply With Quote

Old   April 6, 2009, 21:56
Default
  #5
New Member
 
Yilei Shi
Join Date: Apr 2009
Posts: 5
Rep Power: 7
stanking is on a distinguished road
Someone can answer the last 3 questions????

we assume the following are right.

1. How can I calculate the k? k = sqrt(kx^2+ky^2), is right?

2. V(k) = sqrt(V_kx^2 + V_ky^2), is right??

Then we got two 32*32 matrix for k and V(k) (for example each direction has 32 points)

* + ! ... The k matrix is like this. It's symmetric along the diagonal. This means k(i,j) = k(j,i).
+ * @ ... Now I calculate the E(k). It should be E(k) = V(k(i,j)) ^2 + V(k(j,i))^2 ( Because of k(i,j) = k(j,i) ), is right??
! @ * ...

Last edited by stanking; April 7, 2009 at 04:37.
stanking is offline   Reply With Quote

Old   April 7, 2009, 16:02
Default
  #6
New Member
 
Nathan Grube
Join Date: Mar 2009
Posts: 13
Rep Power: 7
ngrube is on a distinguished road
Let's call Vx u and call Vy v. let us call the Fourier transforms U(kx,ky)=fft(u(x,y)) and V(kx,ky)=fft(v(x,y)).

First I would use a 2D FFT (from FFTW of matlab or whatever you want) to get U(kx,ky) and V(kx,ky). Then you have your spectral info in terms of wavenumber vectors (kx,ky). You probably want a scalar k=sqrt(kx**2+ky**2).

To get your info in terms of a scalar k, you loop through kx and ky and look at the energy in ecah mode and add it to the appropriate "k" bin.

So, define an array E(k)=[0,0,0,....], loop through kx and ky, compute k from kx and ky, round to an integer k, then do E(k) = E(k) + U(kx,ky)**2+V(kx,ky)**2. This turns your 2D data into the 1D E(k) plots that you usually see.

To answer your other question, usually N goes from -M/2 to (M/2)-1 or else from 0 to M where M is the umber of grid points in one direction. This varies depending on FFT package. You need to read the documentation with your particular FFT library or Matlab, etc. Then I would recommend doing an FFT of some simple sinusoids to make sure you get the results you expect and that your understanding of N is correct.

I hope that is correct. I am in a hurry and doing this from memory without proofreading. Good luck.

PS. (kx,ky)=(0,0) corresponds to mean flow. You may wish to ignore that mode.
mgg likes this.
ngrube is offline   Reply With Quote

Old   April 10, 2009, 03:21
Default
  #7
New Member
 
Yilei Shi
Join Date: Apr 2009
Posts: 5
Rep Power: 7
stanking is on a distinguished road
thx ngrube, i have tried it. But the result is not good!!
The curve has big oszillation.
Maybe can you check my code?

clc;
clear all;
close all;

Data = load('tg3dn');

E = zeros(1024,1);

L = 1;
N = 1:32;
kx = (2 * pi .* N / L)';
ky = (2 * pi .* N / L)';
kz = (2 * pi .* N / L)';

v = reshape(Data(:,4),32,32,32);
u = reshape(Data(:,5),32,32,32);
w = reshape(Data(:,6),32,32,32);

V = fftn(v);
U = fftn(u);
W = fftn(w);

for i = 1:32,
for j = 1:32,
for m = 1:32,
k = round(sqrt(kx(i)^2 + ky(j)^2 + kz(m)^2));
E(k) = E(k) + abs(V(i,j,m))^2 + abs(U(i,j,m))^2 + abs(W(i,j,m))^2;
end
end
end

index = find(E~=0);
E_index = E(index);
loglog(index,E_index,index,index.^(-5/3));
Attached Images
File Type: jpg result.jpg (16.9 KB, 443 views)
stanking is offline   Reply With Quote

Old   April 10, 2009, 20:54
Default check your k's
  #8
New Member
 
Nathan Grube
Join Date: Mar 2009
Posts: 13
Rep Power: 7
ngrube is on a distinguished road
I don't have Matlab (or its manual) in front of me right now, but I think your kx, ky, kz look wrong. You need to check the documentation to see what the wavenumbers are and how the FFT'd data is arranged because these vary from program to program. If you have 32 points in one direction, then the largest magnitude of kx should be only 16 (if you don't see why, google "Nyquist limit") So, your kx ky kz look like they go from 0 to 32 rather than -15 to 16 (the exact upper and lower limit can be handled different ways since kx=16 is redundant with kx=-16... ). I am not sure what MatLab uses exactly. If you can't find it in the documentation, just play around and try an FFT of sin(5*x) or something to see where k=(+/-)5 is in the resulting output.

To see what I mean, suppose N = 8. Some common output orders would be arranged with k's as follows: (-4,-3,-2,-1,0,1,2,3) or (0,1,2,3,-4,-3,-2,-1) or (-3,-2,-1,0,1,2,3,4) etc. Note that in place of -4 I could write +4 and it would be equivalent. Some programs (but probably not MatLab) might even give you (-4,-3,-2,-1,0,1,2,3,4) but then -4 and +4 will each be half of the usual kx=4 contribution. I think that MatLab will use one of my first 3 suggestions. You can just use trial and error if you like.
ngrube is offline   Reply With Quote

Old   April 11, 2009, 15:39
Default
  #9
Senior Member
 
N/A
Join Date: Mar 2009
Posts: 177
Rep Power: 7
harishg is on a distinguished road
You need to do fftshift to center the mean flow frequency. Also fft in matlab gives you frequencies from -N/2:N/2-1, scaled accordingly. If you can give your email address, I will send the matlab code that I had used.
Amjad Asad, mgg and mmmn036 like this.
harishg is offline   Reply With Quote

Old   April 12, 2009, 15:32
Default
  #10
New Member
 
Yilei Shi
Join Date: Apr 2009
Posts: 5
Rep Power: 7
stanking is on a distinguished road
thx harishg, this is my email : stanking1983@hotmail.com

I'm waiting for your code!!!
stanking is offline   Reply With Quote

Old   April 24, 2009, 23:08
Default
  #11
New Member
 
Nathan
Join Date: Apr 2009
Posts: 3
Rep Power: 7
nathanpetrelli is on a distinguished road
Yilei,

I was looking at the same problem actually. Could you please send me the code Harish sent you? I would like to try it out. Would help a lot.
Thanks!

nathanpetrelli77@yahoo.com

Nathan
nathanpetrelli is offline   Reply With Quote

Old   April 28, 2009, 11:02
Default
  #12
Member
 
vishwanath somashekar
Join Date: Apr 2009
Posts: 41
Rep Power: 7
vishwa is on a distinguished road
Hi,
i am interested in calculating the spectrum as well. My geometry however, in not regular..

I am running LES using fluent and have saved the velocity fields in the midplane. what would be the procedure to calculate the spectrum...

I have read the many posts here and in books like pope. It really isn't very clear to me. I would appreciate if you could guide me in calculating the spectrum..

The .dat file that i saved for velocity, contains the following header..

Quote:
TITLE = "title"
VARIABLES = X, Y, Z, "z-velocity", "y-velocity", "x-velocity"
ZONE T="midplane", N=63808, E=60912, ET=QUADRILATERAL, F=FEBLOCK
Regards,
vishwa
vishwa is offline   Reply With Quote

Old   April 8, 2010, 14:21
Default
  #13
New Member
 
Germán González Silva
Join Date: Mar 2010
Location: Campinas-SP
Posts: 20
Rep Power: 6
germangsilva is on a distinguished road
Send a message via MSN to germangsilva
Yilei,

I was looking at the same problem actually. Could you please send me the code Harish sent you? I would like to try it out. Would help a lot.
Thanks!

germangsilva@gmail.com
germangsilva is offline   Reply With Quote

Old   April 16, 2010, 16:31
Default Info for every one ...
  #14
New Member
 
Join Date: Apr 2010
Posts: 1
Rep Power: 0
NicoG is on a distinguished road
Yilei or anyone who has the code,

it would be very nice for every newbie to have the matlab code posted here (if it is not too long).
In fact, Pope explains things in details but example on howto implement such as thing in matlab would be helpfull.

Thanks in advance for the help,

Nico
NicoG is offline   Reply With Quote

Old   July 4, 2010, 07:50
Default
  #15
New Member
 
Sumon Saha
Join Date: Jul 2010
Location: Melbourne
Posts: 1
Rep Power: 0
uttom is on a distinguished road
Yilei,

I was looking at the same problem actually. Could you please send me the code Harish sent you? I would like to try it out. Would help a lot.
Thanks!

sumonmebuet@gmail.com
uttom is offline   Reply With Quote

Old   August 26, 2010, 05:14
Default calculating energy spectra
  #16
jet
New Member
 
Join Date: Nov 2009
Posts: 12
Rep Power: 6
jet is on a distinguished road
hi,
i want to calculate an energy spectrum out of the velocity (or better out of the velocity fluctuations in one direction) measured at one point in a turbulent pipe flow calculated with OpenFoam (LES).

i get a velocity value each timestep (0.00004 s).
and here is my way to calculate the energy spectra:

1) calculating the autocorrelation of the velocity vector
2) fftshift of the autocorrelation
3) calculating abs(fft(fftshift(signal))) --> the energy spectrum values

for plotting the spectra i need a frequency-axis and i don't know how to calculate it!
my sampling freq should be 1/0.00004 = 25000 [Hz] ?

what should be my first value on the freq axis? what is my delta freq? and what is the last axis value?

Thanks!
jet is offline   Reply With Quote

Old   August 26, 2010, 10:02
Default
  #17
Senior Member
 
lakeat's Avatar
 
Daniel WEI (老魏)
Join Date: Mar 2009
Location: South Bend, IN, USA
Posts: 688
Blog Entries: 9
Rep Power: 10
lakeat is on a distinguished road
Send a message via Skype™ to lakeat
I guess, when we based on Taylor frozen hypothesis, the fluid particle will move a certain distance in the time step = 0.00004s, and this distance is equivalent to spacial flow structure length in its homogeneous direction (w/o Taylor frozen hypothesis), and hence it is related to the wave number.

Just some thoughts.
__________________
~
Daniel WEI
-------------
NatHaz Modeling Laboratory
Department of Civil & Environmental Engineering & Earth Sciences
University of Notre Dame, USA
Email || My Personal CFD Blog
lakeat is offline   Reply With Quote

Old   January 16, 2011, 11:57
Default code request
  #18
New Member
 
Join Date: Jan 2011
Posts: 1
Rep Power: 0
jitahz2 is on a distinguished road
could you please also send me the Matlab code?

jitahz@yahoo.com


tnx!
jitahz2 is offline   Reply With Quote

Old   March 13, 2012, 14:46
Default Doubts!!
  #19
New Member
 
Clara
Join Date: Sep 2010
Location: Bruxelles
Posts: 8
Rep Power: 5
risku9 is on a distinguished road
Hi everyone!!

I've been visiting this web for a long time trying to understand how to plot the turbulence spectra without good results!!

Can anybody send me the code to see what am I doing wrong? gsclara@gmail.com

Thank you!
risku9 is offline   Reply With Quote

Old   May 4, 2012, 12:55
Default Calculating Energy Spectra
  #20
New Member
 
Join Date: Jun 2009
Posts: 16
Rep Power: 7
ghoshi1983 is on a distinguished road
Hi Guys,
I am also having problems with plotting energy spectrum. Could somebody plz send me the code (zietapk@hotmail.com) also so that i can compare it with mine and figure out what i am doing wrong. Thanks in advance!!
ghoshi1983 is offline   Reply With Quote

Reply

Tags
energy spectrum, fft

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
turbulent energy spectrum problem cfd.newbie Main CFD Forum 1 June 20, 2008 12:00
turbulent energy spectrum cfd.newbie CD-adapco 1 June 19, 2008 23:48
turbulent energy spectrum cfd.newbie FLUENT 0 June 18, 2008 18:34
LES correlation and turbulent energy spectrum Fabian Main CFD Forum 4 October 18, 2005 02:04
Energy Spectrum Emad Khalifa Main CFD Forum 3 June 30, 2003 16:03


All times are GMT -4. The time now is 23:09.