CFD Online Discussion Forums - y+ and u+ values with low-Re RANS turbulence models: utility + testcase

Page 2 of 6

Show 40 post(s) from this thread on one page

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM (https://www.cfd-online.com/Forums/openfoam/)

- - y+ and u+ values with low-Re RANS turbulence models: utility + testcase (https://www.cfd-online.com/Forums/openfoam/86562-y-u-values-low-re-rans-turbulence-models-utility-testcase.html)

florian_krause

July 11, 2011 05:55

Hi grandgo,
if its a division by zero, then its definitely not in

Code:

yPlus = y.y() * uTauAvg / RASModel->nu();

or in nuEff what you are using. nu is non-zero as well as nuEff = nu + nu_t. By the way, why do you use nuEff to calculate yPlus?

So division by zero happens then in

Code:

uPlus = U / uTauAvg;

which makes sense if I look at the output vaina74 posted: avgUGradWall: 0. Clearly, the velocity gradient at the wall

Code:

U.boundaryField()[patchi].snGrad()

gives you zero.

Now, my question is, why is this? I dont know your case, but is it possible by any chance, that this is at time 0 and you initialize your internal velocity field as (0 0 0)? just a guess...

Best,
Florian

grandgo

July 11, 2011 07:14

hi florian,

Quote:

By the way, why do you use nuEff to calculate yPlus?

i'm using nuEff, because i need the tool for a LES case.

Quote:

nu is non-zero as well as nuEff = nu + nu_t

are you sure? i commented the "uplus line" out and it didnt work as well. same error message.

Quote:

Now, my question is, why is this? I dont know your case, but is it possible by any chance, that this is at time 0 and you initialize your internal velocity field as (0 0 0)? just a guess...

i changed the internal velocity field a little bit from (0 0 0) to (0.001 0 0), but it's no use. the error occurs not at time t = 0 anyway....

best regards
grandgo

florian_krause

July 11, 2011 07:40

Hi grandgo,

my understanding is, kinematic viscosity nu shouldn't be zero, since its a physical property, therefore nuEff is not zero. I also think that you should divide by nu and not nuEff even if its a LES case, see

http://www.cfd-online.com/Wiki/Dimen...e_%28y_plus%29

and the source code of the yPlusLES utility.

It really puzzles me where should be a division by zero if you comment out

Code:

uPlus = U / uTauAvg;

since you only divide by nu apart from that.

Is this the only case where you get the error message?

Could you give me a small case where you get the same error message?

Best,
Florian

grandgo

July 11, 2011 11:50

Quote:

Originally Posted by florian_krause (Post 315653)

Code:

uPlus = U / uTauAvg;

since you only divide by nu apart from that.

Is this the only case where you get the error message?

Could you give me a small case where you get the same error message?

Best,
Florian

hi florian,

i don't know why nuEff instead of nu is used in yPlusLES utility. i guess it makes no difference.

and it's not the only case, where i get an error message. i picked the pitzDaily case in incompressible/pisoFoam/les and the same problem there.

best,
grandgo

florian_krause

July 11, 2011 13:50

Thats weird! I just run the pitzDaily case a few second until t=0.017 (random) and used my own plusPostLES utility-it works!

Here are the relevant code snippets of my plusPostLES:

Code:

int main(int argc, char *argv[])

{

    timeSelector::addOptions();



    #include "setRootCase.H"

#   include "createTime.H"

    instantList timeDirs = timeSelector::select0(runTime, args);

#   include "createMesh.H"



    forAll(timeDirs, timeI)

    {

        runTime.setTime(timeDirs[timeI], timeI);

        Info<< "Time = " << runTime.timeName() << endl;

        fvMesh::readUpdateState state = mesh.readUpdate();



        wallDist y(mesh, true);

        // Wall distance

        if (timeI == 0 || state != fvMesh::UNCHANGED)

        {

            Info<< "Calculating wall distance\n" << endl;

            Info<< "Writing wall distance to field "

                << y.name() << nl << endl;

            y.write();

        }



...



#       include "createPhi.H"



        singlePhaseTransportModel laminarTransport(U, phi);



        autoPtr<incompressible::LESModel> sgsModel

        (

            incompressible::LESModel::New(U, phi, laminarTransport)

        );



        volScalarField::GeometricBoundaryField d = nearWallDist(mesh).y();

        volScalarField nuEff = sgsModel->nuEff();



        const fvPatchList& patches = mesh.boundary();

        

        dimensionedScalar utauAvg("utauAvg", dimVelocity, 0);



        scalar nPatch = 0;

                

        forAll(patches, patchi)

        {

            const fvPatch& currPatch = patches[patchi];



            if (typeid(currPatch) == typeid(wallFvPatch))

            {

                yPlus.boundaryField()[patchi] =

                    d[patchi]

                   *sqrt

                    (

                        sgsModel->nu().boundaryField()[patchi]

                       *mag(UMean.boundaryField()[patchi].snGrad())

                    )

                   /sgsModel->nu().boundaryField()[patchi];

                const scalarField& Yp = yPlus.boundaryField()[patchi];



                uTau.boundaryField()[patchi] =

                  sqrt

                  (

            sgsModel->nu()

                   *mag(UMean.boundaryField()[patchi].snGrad())

                   );



                const fvPatchScalarField& utauWall = 

                  uTau.boundaryField()[patchi];

                

                dimensionedScalar utauTmp("utauTmp", dimVelocity, average(utauWall));

                

                utauAvg += utauTmp;



                nPatch ++;



                Info<< "Patch " << patchi

                    << " named " << currPatch.name()

                    << " y+ : min: " << min(Yp) << " max: " << max(Yp)

                    << " average: " << average(Yp) 

                    << " avgUMeanGradWall: " <<  average(mag(UMean.boundaryField()[patchi].snGrad())) << nl << endl;

            }

        }



        utauAvg /= nPatch; 

        

        Info << "Avg. shear velocity is: "

             << utauAvg << "; # of walls: " << nPatch << "." << endl;

       



        yStar = y.y() * utauAvg / sgsModel->nu();

        uStar = U / utauAvg;

Maybe some more info. I am using OF-1.7.x and didn't set any additional FPE options in my bashrc file.

Sorry, at the moment there is nothing else coming to my mind what might cause the error.

Regards,
Florian

grandgo

July 11, 2011 14:22

Quote:

Originally Posted by florian_krause (Post 315723)

Thats weird! I just run the pitzDaily case a few second until t=0.017 (random) and used my own plusPostLES utility-it works!

Here are the relevant code snippets of my plusPostLES:

Code:

int main(int argc, char *argv[])

{

    timeSelector::addOptions();



    #include "setRootCase.H"

#   include "createTime.H"

    instantList timeDirs = timeSelector::select0(runTime, args);

#   include "createMesh.H"



    forAll(timeDirs, timeI)

    {

        runTime.setTime(timeDirs[timeI], timeI);

        Info<< "Time = " << runTime.timeName() << endl;

        fvMesh::readUpdateState state = mesh.readUpdate();



        wallDist y(mesh, true);

        // Wall distance

        if (timeI == 0 || state != fvMesh::UNCHANGED)

        {

            Info<< "Calculating wall distance\n" << endl;

            Info<< "Writing wall distance to field "

                << y.name() << nl << endl;

            y.write();

        }



...



#       include "createPhi.H"



        singlePhaseTransportModel laminarTransport(U, phi);



        autoPtr<incompressible::LESModel> sgsModel

        (

            incompressible::LESModel::New(U, phi, laminarTransport)

        );



        volScalarField::GeometricBoundaryField d = nearWallDist(mesh).y();

        volScalarField nuEff = sgsModel->nuEff();



        const fvPatchList& patches = mesh.boundary();

        

        dimensionedScalar utauAvg("utauAvg", dimVelocity, 0);



        scalar nPatch = 0;

                

        forAll(patches, patchi)

        {

            const fvPatch& currPatch = patches[patchi];



            if (typeid(currPatch) == typeid(wallFvPatch))

            {

                yPlus.boundaryField()[patchi] =

                    d[patchi]

                   *sqrt

                    (

                        sgsModel->nu().boundaryField()[patchi]

                       *mag(UMean.boundaryField()[patchi].snGrad())

                    )

                   /sgsModel->nu().boundaryField()[patchi];

                const scalarField& Yp = yPlus.boundaryField()[patchi];



                uTau.boundaryField()[patchi] =

                  sqrt

                  (

            sgsModel->nu()

                   *mag(UMean.boundaryField()[patchi].snGrad())

                   );



                const fvPatchScalarField& utauWall = 

                  uTau.boundaryField()[patchi];

                

                dimensionedScalar utauTmp("utauTmp", dimVelocity, average(utauWall));

                

                utauAvg += utauTmp;



                nPatch ++;



                Info<< "Patch " << patchi

                    << " named " << currPatch.name()

                    << " y+ : min: " << min(Yp) << " max: " << max(Yp)

                    << " average: " << average(Yp) 

                    << " avgUMeanGradWall: " <<  average(mag(UMean.boundaryField()[patchi].snGrad())) << nl << endl;

            }

        }



        utauAvg /= nPatch; 

        

        Info << "Avg. shear velocity is: "

             << utauAvg << "; # of walls: " << nPatch << "." << endl;

       



        yStar = y.y() * utauAvg / sgsModel->nu();

        uStar = U / utauAvg;

hi florian,

why are you using UMean now? and what's Umean?

and whats the difference between utauavg and utauavg.value() ?

florian_krause

July 11, 2011 15:41

Quote:

why are you using UMean now? and what's Umean?

sorry for the confusion with UMean in the code snippet... I used the utility for a fully developed turbulent pipe flow and time-averaged the velocity field. With that I calculated the mean (time-averaged) velocity gradient at the wall... I only wanted to show you the code of the utility I used for the pitzDaily case.

Quote:

and whats the difference between utauavg and utauavg.value() ?

I don't understand what you mean by this! Are you talking about this one

Code:

utauAvg += utauTmp;

Its just a dimensionScalar with the dimensions of velocity, see defintion of utauTmp. Its not a class with a member value().

Regards,
Florian

grandgo

July 11, 2011 15:50

Quote:

Originally Posted by florian_krause (Post 315739)

Code:

utauAvg += utauTmp;

Its just a dimensionScalar with the dimensions of velocity, see defintion of utauTmp. Its not a class with a member value().

Regards,
Florian

no, in the original code you user this code

Code:

Info << "  avg. friction velocity uTau is: "

             << uTauAvg.value() << " (averaged over " << nPatch << " wall(s))" << nl <<endl;

i'm just wondering...

i don't know what the problem could be. i'm a newbie in c++ programming

the utility would fit perfectly to what i need.
i just go on trying and testing

best
grandgo

florian_krause

July 11, 2011 16:13

Quote:

Its just a dimensionScalar with the dimensions of velocity, see defintion of utauTmp. Its not a class with a member value().

OK now I should say sorry for the non-sense... :(

It is a class, would have to check the documentation for the exact definition... but think it comes from the generic dimensioned Type class dimensioned<Type> with member function value(). That means, if you only use uTauAvg it will output name, value and dimension.

Even I am not a C++ expert...
sorry again!

grandgo

July 13, 2011 05:09

hi florian,

i'm struggling with the modification :)

the thing is: i'm not changing much at all and it still doesn't work.

here is what i changed in the createFields.H file:

Code:

        autoPtr<incompressible::LESModel> sgsModel

        (

            incompressible::LESModel::New(U, phi, laminarTransport)

        );

and in the .C file i substituted RASModel->nu() with sgsModel->nu().
that's all.

the openfoam version is 1.7.1.
i used the utility with the tutorial case: /incompressible/pisoFoam/les/pitzDaily, too. no chance.

i dont understand this.

EDIT: i think, my problem is history.

these two lines in the modified utility caused the problems in my case:

Code:

1) yPlus = y.y() * uTauAvg / nuEff;

2) uPlus = U / uTauAvg;

the problem with the first line was, that i used nuEff instead of sgsModel->nu(). and nuEff did become zero (nuEff = nu + nuSgs, with negative values of nuSgs possible). so i changed to sgsModel->nu().

the second line caused problems, because uTauAvg was zero at time = 0. as florian assumed, my 0/U internelField value is (0 0 0). so the solution is either to change the internalField value, or, as i did, change the code of the utility a little bit:

Code:

timeSelector::addOptions(false, true);

now the 0 directory is excluded in the calculation and everything is fine.

best
grandgo

jona	November 17, 2011 17:05

plusPostRANSUtitility for OF 2.0

1 Attachment(s)

Hi *,

for converting the code of Florian to OF 2.0, "RASModel->nu()" had to be explicitely treated as "volScalarField&".

The corrected code is attached.

Have fun,

Jona

fakekarma

March 8, 2012 05:28

Hi jona,

I've tried to build the utility as follow:

:~/OpenFOAM/OpenFOAM-2.0.1/applications/utilities/postProcessing/wall/plusPostRANS$ wmake libso
wmakeLnInclude: linking include files to ./lnInclude
Making dependency list for source file plusPostRANS.C
SOURCE=plusPostRANS.C ; g++ -m64 -Dlinux64 -DWM_DP -Wall -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -O3 -DNoRepository -ftemplate-depth-100 -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/meshTools/lnInclude -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/transportModels -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/turbulenceModels -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/turbulenceModels/incompressible/RAS/RASModel -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/turbulenceModels/incompressible/LES/LESModel -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/turbulenceModels/LES/LESdeltas/lnInclude -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/finiteVolume/lnInclude -IlnInclude -I. -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/OpenFOAM/lnInclude -I/home/ag-peinke/OpenFOAM/OpenFOAM-2.1.0/src/OSspecific/POSIX/lnInclude -fPIC -c $SOURCE -o Make/linux64GccDPOpt/plusPostRANS.o
'libNULL.so' is up to date.

so it seems to be no error, but I should now use it it, since it doesn't work simply typing plusPostRANS in the case folder. Should I do something more???

Thanks for your help,

Cheers,

Elia

wouter

March 8, 2012 17:29

hello,

type wmake, not wmake libso

hope this works
Wouter

fakekarma

March 9, 2012 02:47

Hallo Wouter,
with the command "wmake" works as it should for compiling.
Thanks very much.

Elia

chegdan

April 16, 2012 17:24

Hello All,

First of all, thanks to Florian for making this, I have used some other yPlus utilities floating around here and this one was particularly interesting. I was attempting to get a uPlus and yPlus for something other than a pipe flow and I have been reading this thread thinking about the problem. Have made a modification to the plusPostRANSUtility utility provided (for OF 2.1.x) and I nearly get the same answer on Florian's test case (running pisoFoam) with my utility and his. However, it is (edit was) dreadfully slow since it searches for the nearest uTau over all wall faces for each cell.

Right now it calculates a new y value (we will call it cellToFaceDistance), compares it to the current already calculated y value created using wallDist, and if there is a user specified fraction difference (i.e. (cellToFaceDistance-y)/cellToFaceDistance) than it is close enough and grabs that uTau for that cell and stops searching.

I have uploaded the code to github and it can be downloaded with

Code:

git clone git://github.com/chegdan/yPlusUplus.git

It is currently made for OF 2.1.x

EDIT: However, after thinking a bit...I'm not quite sure if this sort of u+ vs y+ is useful in cases where there is boundary layer separation around objects. For something where the streamlines are nice close to the object, it makes sense. For something with separation, it may not be best to non-dimensionalize the velocity field with respect to the friction velocity but the mean velocity. Either way, it was a good bit of code and if anyone can find it useful, enjoy.

Cheers.
Dan

The King

April 21, 2012 12:58

Hi! Did you found the error. I get the same, nothing found yet...

eysteinn

June 11, 2012 10:11

Quote:

Originally Posted by florian_krause (Post 301087)

Dear all,

This issue was already discussed in some length in various threads. Therefore I finally reviewed and attached my plusPostRANS utility.

It calculates y+ and u+ values when using one of the available low-Re RANS turbulence models. I hope the code as well as the output will be self-explaining

I also attached a case. It is a straight pipe (wedge) with periodic inlet/outlet boundary conditions. The Reynolds number based on the mean axial velocity is Re=5300. You can run the case using pisoFoam.

I would appreciate any feedback, it might be still improved at some point! :)

Best regards,
Florian

Hi Florian,

I have one question regarding your example case.
I noticed that you use zero gradient for epsilon and nut at the wall. Can you ellaborate on why that is?
(I somehow expected epsilon = 1e-9 or something small as you are solving for epsilon tilda and calculated for nut.)

Best regards,
Eysteinn

owayz

July 24, 2012 22:33

Hi Jona and Florian,
Thanks for the utility I am sure it will help a lot of people.
I am also interested in calculating the y+ for my case. But I have compressible simulation.
I want to know two things:
1- If its a transient simulation shouldn't we use the mean values of U, rho and other variables?
2- How could I extend this utility to calculate the same numbers of compressible flow?

Regards,
Ali

owayz

July 26, 2012 08:42

Compressible y+

1 Attachment(s)

Hi,
I have made some small changes in the utility, developed by Florian/Jona. I have tried to make it work with compressible simulations.
Please check it, if you see its working as it is supposed to work.
Regards,
Awais Ali
P.S: For my case it reports a very low y+ values on the wall (unexpectedly). But the value that it writes seems to be Ok.

owayz

August 13, 2012 18:56

compressibleyPlusLES

1 Attachment(s)

Hi All,
I have written this utility (again by modifying the utility previously developed) to calculate y+ and u+ values for compressible LES simulation.
Please let me know if you find this utility working. For me its working fine. I have used on my own case. The yPlus values are as they were expected to be.

Regards,
Awais

All times are GMT -4. The time now is 03:54.

Page 2 of 6

Last »

Show 40 post(s) from this thread on one page