[byrong]

Hello Foamers,

I'm quite new at CFD and OpenFOAM. I am trying to simulate turbulent pipe flow using the channelFoam solver. My case was set using channel395 as a base. At first, I tried to get the turbulence (Re_tau = 170) with non perturbed initial conditions, but I only got laminar flow. Then, I used the perturbU, and perturbCylinder applications, but the oscillations were damped as the time increased. Then, I created my own perturbations by adding some noise to the U file in the 0 directory, but still got no turbulence. Would you mind helping me with a hint?

Cheers,

Byron

[alexeym]

Originally Posted by byrong

... by adding some noise to the U ...

What kind of noise did you add? It is well-known that uncorrelated noise (let's say simple Gaussian noise) in U decays very quickly and does not lead to production of eddies.

[byrong]

Thank you for your quick response Alexey,

Yes, I have been adding Gaussian noise, I also used perturbU, and perturbCylinder, but still got no turbulence. Do you have any suggestions of what should I do in order to obtain turbulent structures?

Cheers,

Byron

[alexeym]

Well, it partly depends on the nature of your case. For example you can use boxTurb utility to initialize flow (as it uses FFT, number of cells in certain directions should be power of 2).

"...turbulent pipe flow..." is rather broad term, you should be more specific.

[byrong]

Thank you Alexey,

I will try to do what you recommend. About pipe flow, I am trying to simulate turbulent flow at different Reynolds numbers. This pipe has a 30*D length, and I am using cyclic boundary conditions in the inlet and the outlet.

Regards,

Byron

[byrong]

Hello Alexey,

...For example you can use boxTurb utility to initialize flow...

I already tried to use the boxTurb application. It seems great since it uses the discrete Fourier transform. However, when I run it, I get this error.

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0

Reading field U

Reading boxTurbDict



--> FOAM FATAL ERROR:
calculated number of cells is incorrect

From function Kmesh::Kmesh(const fvMesh& mesh)
in file Kmesh/Kmesh.C at line 84.

FOAM aborting

#0 Foam::error:rintStack(Foam::Ostream&) in "/opt/openfoam211/platforms/linuxGccDPOpt/lib/libOpenFOAM.so"
#1 Foam::error::abort() in "/opt/openfoam211/platforms/linuxGccDPOpt/lib/libOpenFOAM.so"
#2 Foam::Kmesh::Kmesh(Foam::fvMesh const&) in "/opt/openfoam211/platforms/linuxGccDPOpt/lib/librandomProcesses.so"
#3
in "/opt/openfoam211/platforms/linuxGccDPOpt/bin/boxTurb"
#4 __libc_start_main in "/lib/i386-linux-gnu/libc.so.6"
#5
in "/opt/openfoam211/platforms/linuxGccDPOpt/bin/boxTurb"
Aborted (core dumped)

[alexeym]

To understand the reason of the error we should look at Kmesh.C:

    boundBox box = mesh.bounds();
    l_ = box.span();

    vector cornerCellCentre = ::Foam::max(mesh.C().internalField());
    vector cellL = 2*(box.max() - cornerCellCentre);

    vector rdeltaByL;
    label nTot = 1;

    forAll(nn_, i)
    {
        nn_[i] = label(l_[i]/cellL[i] + 0.5);
        nTot *= nn_[i];

        if (nn_[i] > 1)
        {
            l_[i] -= cellL[i];
        }

        rdeltaByL[i] = nn_[i]/(l_[i]*l_[i]);
    }

So as you can see calculated mesh size is estimated by the thickness of the last cell. If you have grading in the mesh obviously Kmesh will fail as it actually expects cells to be of the same size (per direction).

You can try to modify the code of Kmesh like this (this part is located just below forAll loop I've posted):

    if (nTot != mesh.nCells())
    {
        FatalErrorIn("Kmesh::Kmesh(const fvMesh& mesh)")
        << "calculated number of cells is incorrect.\n"
            << "Mesh size: " << mesh.nCells() << "\n"
            << "Calculated size: " << nTot
            << abort(FatalError);
    }

for the error message to be more informative:

--> FOAM FATAL ERROR: 
calculated number of cells is incorrect.
Mesh size: 60000
Calculated size: 249600

    From function Kmesh::Kmesh(const fvMesh& mesh)
    in file Kmesh/Kmesh.C at line 76.

(after modification you need to run wmake in $FOAM_SRC/randomProcesses)

So I guess, you've tried channel395 mesh and it does have grading. If you'd like to have grading in your mesh and initialize it with boxTurb, you can initialize flow on the uniform mesh and then use mapFields utility to interpolate between meshes (uniform and graded).

[byrong]

Thanks a lot Alexey.

That was really helpful. I actually have grading in my mesh. I'll try to do what you recommend.

Cheers

[morard]

Hi Byron,

an easy solution would be to increase your Re in order to get turbulence. Afterwards, just map the velocity field to another case with lower bulk velocity.

Regards

[byrong]

Thank you Dejan,

I will also try to do what you suggest.

Cheers.

[cedric_duprat]

Morning Byron,

I've ran such a case few years ago.
You could have a look there :
http://www.cfd-online.com/Forums/ope...pe-flow-2.html

#39

One idea is to use a coarse grid first, with cyclic BC to let the turbulence develop then map this results into your fine grid.

Notice that in PerturbU you can also increase the amplitude of your initial noise.

Then if I remember correctly, I've managed also add a swirl in the flow to help the turbulence develop. I switched this swirl off after few second of calculation.

Then, as Dejan suggested, increase your Re seems to be a good idea in order to get turbulence.

I hope this will help you,

Cedric

[byrong]

Thank you Cedric,

I will try to do what you suggest.

Regards,

Byron

[byrong]

Thank you Cedric, Alexey and Dejan,

I was finally able to get the turbulence. What I did was to first create turbulence in a square channel without bounded walls which fits within my pipe, and has the same length with the boxTurb utility. Then, I perturbed my wall-bounded pipe with perturbCylinder, and I mapped the velocity fields of the box within the pipe. That was my initial condition. Then I ran the simulation for 500 seconds, and I finally got some continous turbulence.

Cheers.

[byrong]

[QUOTE=cedric_duprat;486826] I've ran such a case few years ago. [Quote]

Hello Foamers,

I have been trying to obtain the averages of the velocity fields from the different time directories to obtain the U+ vs y+ curve by using the sample utility. However, I has become very laborious and I have been trying to make a script for this task but it will take me a while. Would you mind suggesting me an easier way to obtain this data?

Regards.

[aylalisa]

Dear byrong,


please have a look at this thread started by Florian Krause: http://www.cfd-online.com/Forums/ope...tml#post301087

The utility has been changed/expanded during time by others.


formulas:

http://www.cfd-online.com/Forums/ope...tml#post430703

http://www.cfd-online.com/Forums/ope...tml#post479837


I've used combination of different versions:

yPlusLES.tar

I am not clear about line 173 in file yPlusLES.C
Does this line map from the closest perpendicular boundary field face to the current face to have for computation of ?

Is it correct to map to faces of cells inside the volume to finally compute and ?

was at first computed in the cells next to the wall, therefore I understand for this case that you need according to

could be computed with gradient between current cell and neighbour cell instead of mapping the value from closest boundary field cell. But obviously must be mapped from boundary cell. Could you explain me why? I have a problem with the underlying physics.

You've mentioned that you 've used perturbU as well. I haven't read the paper ("Coherent structure dynamics in near wall turbulence") yet.
Is there a possibility to validate the perturbed field U (before using it)?


Regards,
Aylalisa

[byrong]

Hello Aylalisa,

First, thanks a lot for your references. They are really helpful.

I may not be the best person to reply your questions since I am just becoming familiar with openFoam.

If you want to validate perturbU, you should check Eugene de Villier's PhD thesis. He coded the perturbU utility. You can find his thesis in the following link http://foamcfd.org/resources/theses.html.

Cheers,

Byron

[cedric_duprat]

Dear Aylalisa,

I'm not sure about the line you are talking about but regarding this code:
http://foam.sourceforge.net/docs/cpp/a08929_source.html
I guess you are talking about these lines:

    81             yPlus.boundaryField()[patchi] =
    82                 d[patchi]
    83                *sqrt
    84                 (
    85                     nuEff.boundaryField()[patchi]
    86                    *mag(U.boundaryField()[patchi].snGrad())
    87                 )
    88                /nuLam.boundaryField()[patchi];

In this formula, is calculate like this:



The class member snGrad is the gradient of velocity normal to each face of the wall patches (you could find more information about this snGrad in the forum). and is your closest point to the wall in the normal direction (volScalarField::GeometricBoundaryField d = nearWallDist(mesh).y().

About PerturbU, one easy way to validate (or to see what it's doing) is to applied it to one time step. You can then visualize this time step.
You should/will see a nice parabolic profile and if you zoom at the wall, you'll see the sinusoidal noise added close to the wall in the spanwise and wall normal direction.

I hope this will help you,

Cedric

[aylalisa]

Dear Cedric,

Your explanation is very good for me. I've already tried to understand the code (C++ beginner) for computation of y+ and concluded that (with this formula) I can compute y+ 'only' for cells next to the walls (???). Your explanation seems to confirm this.

I've mixed up the uploaded versions for computation of yPlus, excuse me! Therefore you could not find the code line I've mentioned.

First I've started with florian krause's code: plusPostRANSUtility.tar.gz
(http://www.cfd-online.com/Forums/ope...tml#post301087)

Following the thread I stumbled over http://www.cfd-online.com/Forums/ope...tml#post354995
--> yPlusUPlus!

My question refers to line 180 in the version 'yPlusUPlus' offered by chegdan, #35. I would like to get yPlus and uPlus for the domain, to create plots yPlus/uPlus at different positions in the domain.
Therefore I use the code from chegdan with some question marks in mind(line 180, in my personal version it was line 173).

I see the sinusoidal noise close to the walls! Thanks!


Best regards,
Aylalisa

[wchen]

Originally Posted by byrong

Thank you Cedric, Alexey and Dejan,

I was finally able to get the turbulence. What I did was to first create turbulence in a square channel without bounded walls which fits within my pipe, and has the same length with the boxTurb utility. Then, I perturbed my wall-bounded pipe with perturbCylinder, and I mapped the velocity fields of the box within the pipe. That was my initial condition. Then I ran the simulation for 500 seconds, and I finally got some continous turbulence.

Cheers.

Hi Byron,

I am trying to run exactly the same case as you, turbulent pipe flow with Re_tau = 180. I am also new to OpenFoam. Could you please elaborate more on how to modify the perturbCyl code to fit into your case? I downloaded perturbCyl from the forum but am not sure how to interpret the code.

Cheers,
Winson Chen

[byrong]

Originally Posted by wchen

Hi Byron,

I am trying to run exactly the same case as you, turbulent pipe flow with Re_tau = 180. I am also new to OpenFoam. Could you please elaborate more on how to modify the perturbCyl code to fit into your case? I downloaded perturbCyl from the forum but am not sure how to interpret the code.

Cheers,
Winson Chen

Hello Winson,

I assume you already got a proper mesh. However, to obtain a proper perturbation, I kind of cheated because I tried thousands of times using the perturbCylinder, and it did not work. In fact, the sine waves always became laminar. Therefore, what I did was the following:
1. Create a square channel with only two plates at the top and at the bottom which fits within the pipe you have.
2. Generate an initial perturbation using the boxTurb utility. Check a tutorial to use boxTurb first since the mesh must be divided in powers of two in each direction.
3. Run the case with the same velocity and viscosity parameters as the ones you will use in your pipe maybe for 5 to 10 cycles within the flow domain, or until you reach constant turbulence.
4. Once you got the turbulence map the velocity field of the squared channel to the pipe geometry with developed flow as an initial condition (both geometries should have the same length). To do this you require to use the mapFields utility, and create a mapFieldsDict within the system directory (there are plenty of tutorials of this online).
5. Run your pipe case. At the begining, the initial condition will look weird, but as long as the flow develops, you will get some nice turbulence. Run the case through the pipe domain for at least 10 complete cycles. Then use your last time step as a new initial condition, and run your case again for another 5 cycles. By doing this I got really nice statistics. Indeed, my statistics were very close to some DNS statistics that my supervisor gave me.

I hope it helps.

If you still need to use perturbCyl, you might have to write to Eugene de Villers, since he coded that utility, but in my case it did not work. And if you make it work properly, please let me know.


Cheers,

Byron