[hiter]

Hi everyone, I am focusing on the generation of HIT by using the method proposed by Rogallo.
In fourier space, the mass conservation i.e. kx*ux + ky*uy + kz*uz = 0 reaches the mechine accuraty (to the power -14).
However, after inverse Fourier tranfer, DIV*U approximaly O(1), which is not divergence free at all.
Also, I found that after the inverse Fourier of the U(k), the velocity field U(X) has image part.
I comfirmed that using both part of U in real space U(real, image ), by FFT, in Fourier Space the equation indeed valid, kx*ux + ky*uy + kz*uz to the O(-12).

I am quite sure that the generated random velocity in wave number space obeys the symmerty rules for both (real part and the image part).

Does anyone have experience about this?

Thanks for helping!

[cfdnewbie]

There is a sign error in Rogallo's original paper, that might be the problem. You should recalculate the k1 k2 k3 summands, I believe it was the k3 one...
maybe that causes the divergence issues.
Otherwise, check your FFT for a simple for -and backward problem.

[hiter]

Originally Posted by cfdnewbie

There is a sign error in Rogallo's original paper, that might be the problem. You should recalculate the k1 k2 k3 summands, I believe it was the k3 one...
maybe that causes the divergence issues.
Otherwise, check your FFT for a simple for -and backward problem.

hi cfdnewbie,
Thanks. I also noticed that error,. Yes, as you mentiond, it is K3. Actually, I added the minus sign.
Currently, I performedd the 3D FFT/IFFT by MKL. The results are in accord with the samples of other references.

[cfdnewbie]

hm, how do you calculate the divergence in physical space, i.e. how do you build the derivatives?

[hiter]

Originally Posted by cfdnewbie

hm, how do you calculate the divergence in physical space, i.e. how do you build the derivatives?

Thanks, cfdnewbie
I both used the second order and the fourth order. i.e.,
- U(i-1,j,k) + U(i+1,j,k)
- V(i,j-1,k) + V(i,j+1,k)
- W(i,j,k-1) + W(i,j,k+1)
This part is no problem, I think.

Currently, I really guess it is related to IFFT.
As your previous suggested, I start to check the IFFT from 1D case.
I generate 1D array, following the realtionship C^*(-k)= C(k).
After IFFT, the array still possess image part. it is really strange.
Beside the above metioned symmerty constrain, maybe i need to add other condions, to get a real array, after IFFT.

[hiter]

YES, When I add the constrains for the random generated image part Imag(0)=Imag(N1/2) = 0.0.
I can obtain the desired resluts.

[cfdnewbie]

Originally Posted by hiter

Thanks, cfdnewbie
I both used the second order and the fourth order. i.e.,
- U(i-1,j,k) + U(i+1,j,k)
- V(i,j-1,k) + V(i,j+1,k)
- W(i,j,k-1) + W(i,j,k+1)
This part is no problem, I think.

Just keep in mind that while you compute the divergence globally (spectrally) in wave space, your local derivatives in physical space are a lot less accurate. So while you might be divergence free in wave space (e-15), I wouldn't expect better than something like e-4 ... e-6 with second and 4th order in physical space.

Currently, I really guess it is related to IFFT.
As your previous suggested, I start to check the IFFT from 1D case.
I generate 1D array, following the realtionship C^*(-k)= C(k).
After IFFT, the array still possess image part. it is really strange.
Beside the above metioned symmerty constrain, maybe i need to add other condions, to get a real array, after IFFT.

Yes, then maybe you didn't fill the wave number array for the ifft correctly? did you put the Nyquist and mean mode correctly? Could you use Matlab for a quick check?


Here is an easy way to check the full 3D algorithm: Use the Taylor Green vortex initial condition (3D). The velocity field is analytically divergence free, you can calculate it by hand. Then transfer to wave space, calculate the divergence there. That should also be exactly 0. Now return to wave space, using your algorithm, and reevaluate the divergence!
That should help you identify the error!

[hiter]

Originally Posted by cfdnewbie

Just keep in mind that while you compute the divergence globally (spectrally) in wave space, your local derivatives in physical space are a lot less accurate. So while you might be divergence free in wave space (e-15), I wouldn't expect better than something like e-4 ... e-6 with second and 4th order in physical space.


Yes, then maybe you didn't fill the wave number array for the ifft correctly? did you put the Nyquist and mean mode correctly? Could you use Matlab for a quick check?


Here is an easy way to check the full 3D algorithm: Use the Taylor Green vortex initial condition (3D). The velocity field is analytically divergence free, you can calculate it by hand. Then transfer to wave space, calculate the divergence there. That should also be exactly 0. Now return to wave space, using your algorithm, and reevaluate the divergence!
That should help you identify the error!

Thanks

I know that after IFFT, in real space, using the CD schme the error could not be O(-14). But I my case, it is really large, sometimes approximaly equals zero. So the deviation is obvious.

I am sorry now, I am only familar with fortran program. But I can trust the MKL library issued by Intel. If i made mistake, it should to be related with what you mentioned the Nyquist and mean mode.

Did you mean the oddball wavenumber? I seems to neglect this point.
This is the key point, Thanks.

Could you recommend some refereces for help? I need to extend to three dimentionals.

[cfdnewbie]

check out this link for the Taylor Green vortex :http://dept.ku.edu/~cfdku/hiocfd/case_c3.5.pdf
initial conditions are given on the first page.
check out also this dissertation, http://users.ugent.be/~dfauconn/diet...onnier_PhD.pdf and read the Taylor green vortex chapter.

To get you started on the FFT, definitely read this:
http://code.google.com/p/p3dfft/

It is a Fortran 3D FFT, it has real to complex and complex to real transform, examples etc. Use this, and that you help you a great deal.

[hiter]

Originally Posted by cfdnewbie

check out this link for the Taylor Green vortex :http://dept.ku.edu/~cfdku/hiocfd/case_c3.5.pdf
initial conditions are given on the first page.
check out also this dissertation, http://users.ugent.be/~dfauconn/diet...onnier_PhD.pdf and read the Taylor green vortex chapter.

To get you started on the FFT, definitely read this:
http://code.google.com/p/p3dfft/

It is a Fortran 3D FFT, it has real to complex and complex to real transform, examples etc. Use this, and that you help you a great deal.

I will study now.
Thanks. Really approiate.

[anzillo]

Hello Hiter

I am also working on the same problem. I found that MATLAB has the ifftn with the option of 'symmetric', which is used if the variable that is being inverse-transformed does not obey the rule of 'conjugate symmetry'.

However, that options also leads to highly divergent velocity field (10^-1). I am trying to write a code for converting the wavenumber space velocity field matrices into conjugate symmetric ones, following which I will use the normal ifftn function.

I will let you know how it goes. In the meanwhile, I would be thankful if you could inform me in case you are able to find a solution.

Success !

[hiter]

Sorry the answer so later. Acutally, I change the form of the power spectrum. So that there is little energy in the high wave number region.
In this way, the mass consvation up to 10^-4.
At first, you need to guaratee that the random field meet the conjugate symmetry.
For three dimensional, I think it is difficult to avoid the odd ball effects.
Hope still help you.

Originally Posted by anzillo

Hello Hiter

I am also working on the same problem. I found that MATLAB has the ifftn with the option of 'symmetric', which is used if the variable that is being inverse-transformed does not obey the rule of 'conjugate symmetry'.

However, that options also leads to highly divergent velocity field (10^-1). I am trying to write a code for converting the wavenumber space velocity field matrices into conjugate symmetric ones, following which I will use the normal ifftn function.

I will let you know how it goes. In the meanwhile, I would be thankful if you could inform me in case you are able to find a solution.

Success !