PANS Turbulence Model Implementation Guglielmo Minelli Chalmers - - PowerPoint PPT Presentation

▶
pans turbulence model implementation
SMART_READER_LITE
LIVE PREVIEW

PANS Turbulence Model Implementation Guglielmo Minelli Chalmers - - PowerPoint PPT Presentation

Feb 16, 2023 201 likes •547 views

Outline Intro Turbulence models The test case PANS k implementation k f Conclusion PANS Turbulence Model Implementation Guglielmo Minelli Chalmers University of Technology Applied Mechanics Vehicle Aerodynamic

slide-1
SLIDE 1

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

PANS Turbulence Model Implementation

Guglielmo Minelli

Chalmers University of Technology Applied Mechanics Vehicle Aerodynamic Laboratory

December 01, 2014

1 / 33

slide-2
SLIDE 2

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Intro Solving methods Turbulence models The implemented models The test case PANS k − ε implementation .C modification .H modification k − ε − ζ − f .C modification .H modification Wall functions Conclusion

2 / 33

slide-3
SLIDE 3

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Solving Methods

  • RANS, LES, DNS

Figure : Energy cascade

3 / 33

slide-4
SLIDE 4

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Hybrid Methods

  • Zonal approaches, switch the solving equations during the

simulation or inside the domain.

  • DES and PANS
  • DES is a hybrid LES/RANS approach.
  • Fixed zonal approach, the RANS zones are defined a priori

(close to wall)

  • PANS is a hybrid DNS/RANS approach.
  • Varying zonal approach, the RANS equations are employed
  • nly where the grid can not effort a DNS resolution.

4 / 33

slide-5
SLIDE 5

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

URANS k − ε closure equation system

∂k ∂t + Uj ∂k ∂xj = Pk − ε + ∂ ∂xj

  • ν + νt

σk ∂k ∂xj

  • ∂ε

∂t + Uj ∂ε ∂xj = C1ε ε k (Pk) − C2ε ε2 k + ∂ ∂xj

  • ν + νt

σε ∂ε ∂xj

  • C1ε = 1.44, C2ε = 1.92, σε = 1.3 and σk = 1
  • νt is updated every time step as

νt = Cµ k2 ε ; Cµ = 0.09

  • The method can catch only large oscillation such as vortex

shedding

5 / 33

slide-6
SLIDE 6

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

PANS k − ε closure equation system

∂ku ∂t + Uj ∂ku ∂xj = Pu − εu + ∂ ∂xj νu σku ∂ku ∂xj

  • ∂εu

∂t + Uj ∂εu ∂xj = C1εu εu ku Pu − C ⋆

2ε

ε2

u

ku + ∂ ∂xj νu σεu ∂ε ∂xj

  • νu = Cµ

k2

u

εu C ⋆

2ε = C1ε + fk(C2ε − C1ε)

fk = 1 √cµ ∆ Λ 2/3 ; Λ = k3/2 ε ; ∆ = (∆x∆y∆z)1/3

  • fk, depends from the Taylor scale and the cell size.
  • fk can be also fixed.

6 / 33

slide-7
SLIDE 7

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

PANS k − ε − ζ − f closure equation system

∂ku ∂t + Uj ∂ku ∂xj = Pu − εu + ∂ ∂xj νu σku ∂ku ∂xj

  • ∂εu

∂t + Uj ∂εu ∂xj = C1εu εu ku Pu − C ⋆

2ε

ε2

u

ku + ∂ ∂xj νu σεu ∂ε ∂xj

  • ∂ζu

∂t + Uj ∂ζu ∂xj = fu − ζu ku Pu + ζ ku εu(1 − fk) + ∂ ∂xj νu σζu ∂ζ ∂xj

  • L2

u∇2fu − fu = 1

Tu

  • c1 + c2

Pu εu ζu − 2 3

  • νu = Cµζ k2

u

εu c1 = 0.4; c2 = 0.65

  • Two more equation for the low Re approach
  • Viscosity depends on ζ

7 / 33

slide-8
SLIDE 8

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

The test case

Figure : Backstep test geometry dimension.

  • 2D geometry
  • U = 10m/s; ν = 1 × 10−5
  • folder you can download and copy in your run:

tutorial_backStep_PANS_fixed_fk

8 / 33

slide-9
SLIDE 9

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

PANS k − ε implementation

  • mkdir src/turbulenceModels/incompressible/RAS/ in

your own folder

  • Download kEpsilonPANS_fixedFK.gz and extract the model

in your own folder.

  • cp -r kEpsilonPANS_fixedFK kEpsilonPANS_varFK
  • change the .C and .H file name inside the

kEpsilonPANS_varFK folder just created with kEpsilonPANS_varFK.C and kEpsilonPANS_varFK.H

  • sed -i s/kEpsilonPANS_fixedFK/kEpsilonPANS_varFK/g

kEpsilonPANS_varFK.*

  • Add the new .C file to the Make/files file.
  • Ready to modify the .C and .H files.

9 / 33

slide-10
SLIDE 10

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsilonPANS varFK.C

  • We need to add the coefficient Cµ=Cmi_

Cmi_ ( dimensioned<scalar>::lookupOrAddToDict ( "Cmi", coeffDict_, 0.22 ) ), 10 / 33

slide-11
SLIDE 11

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsilonPANS varFK.C

  • We need to change the constant fk =fK_ and C ⋆

2ε =C2U to

scalarField and add the cellVolume_ scalarField

cellVolume_ ( IOobject ( "cellVolume", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar("zero",dimVolume,0.0) ), fK_ ( IOobject ( "fK", runTime_.timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ ), 11 / 33

slide-12
SLIDE 12

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsilonPANS varFK.C

C2U ( IOobject ( "C2U", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar("zero",0.0) ), 12 / 33

slide-13
SLIDE 13

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsilonPANS varFK.C

  • The inserted variables should be initialized before the

calculation as

{ cellVolume_.internalField() = mesh_.V(); C2U = C1_ + (fK_/fEpsilon_)*(C2_-C1_); .... .... }

  • Delete the old definition of C2U
  • Sobstitute bound(kU\_, fK\_*kMin\_); with

bound(kU\_, min(fK\_)*kMin\_); after the kU_ equation

13 / 33

slide-14
SLIDE 14

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsilonPANS varFK.C

  • Remember to delete the line

fK\_.readIfPresent(coeffDict());

  • At the end of each time step the factor fk and C ⋆

2ε are

recalculated as

// Re calculation of fK and C2U for PANS fK is calculated cell by // cell and upper and lower limited to avoid unfisical solution for (int i=0; i<fK_.size(); i++) { fK_[i]=(1/sqrt(Cmi_.value()))*(pow(pow(cellVolume_[i],1.0/3.0)/ (pow(k_[i],3.0/2.0)/epsilon_[i]),2.0/3.0)); if (fK_[i]>1) //upper limited fK_[i]=1; if (fK_[i]<0.1) //lower limited fK_[i]=0.1; } C2U =C1_+ (fK_/fEpsilon_)*(C2_-C1_); 14 / 33

slide-15
SLIDE 15

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsilonPANS varFK.H

  • We need to add the new foelds and coefficients in the .H file
  • Remember to delete the old declaration of fk

// Model coefficients dimensionedScalar Cmi_; //new dimensionedScalar Cmu_; dimensionedScalar C1_; dimensionedScalar C2_; dimensionedScalar sigmaEps_; dimensionedScalar fEpsilon_; // Fields volScalarField cellVolume_; //new volScalarField fK_; //new volScalarField C2U; //new volScalarField k_; volScalarField epsilon_; volScalarField kU_; volScalarField epsilonU_; volScalarField nut_; 15 / 33

slide-16
SLIDE 16

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Case folder modification

  • Download and extract the case folder backstep
  • Add the file fK to the 0 folder
  • Change the RASProperties file to

RASModel kEpsilonPANS_varFK;

  • Add the following to the controlDict file to average quantities

functions { fieldAve { type fieldAverage; functionObjectLibs ("libfieldFunctionObjects.so"); enabled true;

  • utputControl
  • utputTime;

cleanRestart False; timeStart 30; timeEnd 100; fields ( U { mean

  • n;

prime2Mean on; base time; } ); } } 16 / 33

slide-17
SLIDE 17

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Case folder modification

  • Add the following to the controlDict file to read your library

libs ("libmyIncompressibleRASModels.so");

  • Start simulation creating a log file

pimpleFoam > log&

  • Check the simulation with the gnuplot script ”Residual”

typing gnuplot Residual

  • Change the RASProperties to compare different solvers

RASModel kEpsilonPANS{\_}fixedFK; RASModel kEpsilon;

17 / 33

slide-18
SLIDE 18

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Starting from the last implemented model we can implement

the four equation model.

  • Change the value of C1 and C2 to 0.4 and 0.65 respectively
  • Insert three more constants after fEpsilon_

CL_ ( dimensioned<scalar>::lookupOrAddToDict ( "CL", coeffDict_, 0.36 ) ), 18 / 33

slide-19
SLIDE 19

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

Cmi_ ( dimensione CTau_ ( dimensioned<scalar>::lookupOrAddToDict ( "CTau", coeffDict_, 6 ) ), Ceta_ ( dimensioned<scalar>::lookupOrAddToDict ( "Ceta", coeffDict_, 85 ) ), 19 / 33

slide-20
SLIDE 20

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Ceps1_ is now a scalar field that change at each timestep

CEps1_ ( IOobject ( "CEps1", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar("zero",0.0) ), 20 / 33

slide-21
SLIDE 21

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • ζ and f can be introduced after Ceps1_

zitaU_ //parameter of zeta transport equation ( IOobject ( "zitaU", runTime_.timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ ), fU_ //parameter of poisson equation ( IOobject ( "fU", runTime_.timeName(), mesh_, IOobject::MUST_READ, IOobject::AUTO_WRITE ), mesh_ ), 21 / 33

slide-22
SLIDE 22

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Length and time scale are now introduced

lU_ //Length scale parameter ( IOobject ( "lU", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar("zero",dimLength,0.0) ), tU_ //Time scale parameter ( IOobject ( "tU", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar("zero",dimTime,0.0) ), 22 / 33

slide-23
SLIDE 23

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Bounding for f is now introduced

BfU_ //needed to bound fU ( IOobject ( "BfU", runTime_.timeName(), mesh_, IOobject::NO_READ, IOobject::AUTO_WRITE ), mesh_, dimensionedScalar("zero",dimTime/(sqr(dimTime)),0.0001) ), 23 / 33

slide-24
SLIDE 24

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Initialization of the introduced quantities before the member

function section

CEps1_=1.4*(1.0+0.045/(sqrt(zitaU_))); cellVolume_.internalField() = mesh_.V(); C2U = CEps1_ + (fK_/fEpsilon_)*(CEps2_-CEps1_); tU_=max(kU_/epsilon_, CTau_*sqrt(nu())/sqrt(epsilon_));//time scale lU_=CL_*max(pow(kU_,3.0/2.0)/epsilon_, Ceta_*pow(nu()*nu()*nu()/epsilon_, 1.0/4.0)); //length scale 24 / 33

slide-25
SLIDE 25

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Introduce modification to ε equation

tmp<fvScalarMatrix> epsUEqn ( fvm::ddt(epsilonU_) + fvm::div(phi_, epsilonU_)

  • fvm::laplacian(DepsilonUEff(), epsilonU_)

== CEps1_*G*epsilonU_/kU_

  • fvm::Sp(C2U*epsilonU_/kU_, epsilonU_)

); 25 / 33

slide-26
SLIDE 26

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Introduce equation for ζ and f

//Unresolved zita equation tmp<fvScalarMatrix> zitaUEqn ( fvm::ddt(zitaU_) + fvm::div(phi_, zitaU_)

  • fvm::laplacian(DzitaUEff(), zitaU_)

==

  • fvm::Sp(G/kU_, zitaU_)

+ fU_ + fvm::Sp(epsilonU_*(1.0-fK_)/kU_, zitaU_) ); zitaUEqn().relax(); zitaUEqn().boundaryManipulate(zitaU_.boundaryField()); solve(zitaUEqn); bound(zitaU_, min(zitaU_)+0.0001); //bounding of zeta (dimensionless) Info << "zita done" << endl;//to check where the simulation stop //Unresolved fU equation tmp<fvScalarMatrix> fUEqn ( (fvm::laplacian(fU_))

  • fvm::Sp(1/sqr(lU_), fU_)
  • 1.0/(tU_*sqr(lU_))*(C1_+C2_*G/epsilonU_)*(zitaU_-2.0/3.0)

); fUEqn().relax(); fUEqn().boundaryManipulate(fU_.boundaryField()); solve(fUEqn); bound(fU_, min(BfU_));//bounding of fU 26 / 33

slide-27
SLIDE 27

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.C

  • Update the introduced scales.

for (int i=0; i<fK_.size(); i++) { fK_[i]=(1/sqrt(Cmi_.value()))*(pow(pow(cellVolume_[i],1.0/3.0)/ (pow(k_[i],3.0/2.0)/epsilon_[i]),2.0/3.0)); if (fK_[i]>1) fK_[i]=1; if (fK_[i]<0.1) fK_[i]=0.1; } CEps1_=1.4*(1.0+0.045/(sqrt(zitaU_))); C2U =CEps1_+ (fK_/fEpsilon_)*(CEps2_-CEps1_); // updating scales tU_=max(kU_/epsilon_, CTau_*sqrt(nu())/sqrt(epsilon_)); //time scale lU_=CL_*max(pow(kU_,3.0/2.0)/epsilon_, Ceta_*pow(nu()*nu()*nu()/epsilon_, 1.0/4.0)); Info << "Defining the new PANS k eps z f model" << endl; 27 / 33

slide-28
SLIDE 28

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.H

  • Name the new fields and model coefficients in the .H file

// Model coefficients dimensionedScalar Cmi_; dimensionedScalar Cmu_; dimensionedScalar C1_; \\new dimensionedScalar C2_; \\new dimensionedScalar CEps2_; dimensionedScalar sigmaEps_; dimensionedScalar fEpsilon_; dimensionedScalar CL_; \\new dimensionedScalar CTau_; \\new dimensionedScalar Ceta_; \\new // Fields volScalarField cellVolume_; volScalarField fK_; volScalarField C2U; volScalarField CEps1_; \\new volScalarField zitaU_; \\new volScalarField fU_; \\new volScalarField lU_; \\new volScalarField tU_; \\new volScalarField BfU_; \\new volScalarField k_; volScalarField epsilon_; volScalarField kU_; volScalarField epsilonU_; volScalarField nut_; 28 / 33

slide-29
SLIDE 29

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.H

  • Add the effective diffusivity for ζ.

//- Return the effective diffusivity for unresolved zita tmp<volScalarField> DzitaUEff() const { return tmp<volScalarField> ( new volScalarField("DzitaUEff", nut_/(fK_*fK_*sigmaEps_/fEpsilon_) + nu()) ); } 29 / 33

slide-30
SLIDE 30

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Low Reynolds wall functions

  • New boundary condition should be compiled and added to the

dynamic library.

  • Copy the folder

OpenFOAM-2.3.x/src/turbulenceModels/incompressible/ RAS/derivedFvPatchFields/ wallFunctions/epsilonWallFunctions in your own RAS folder

  • Change name of the folders and the files to

zitaUWallFunctions

  • Repeat the first two points for the fU equation boundary

condition.

  • Using the sed command substitute epsilon with zitaU or fU.
  • Now we have four new boundary condition for high and low

Reynolds number approach.

30 / 33

slide-31
SLIDE 31

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

kEpsZitaF.H

  • A small change in the wall function is done for f . fwall = −2νζ

y 2

if (yPlus > yPlusLam_) { fU[cellI] = w*Cmu75*pow(k[cellI], 1.5)/(kappa_*y[faceI]); } else { fU[cellI] = -2.0*zitaU[cellI]*nuw[faceI]/sqr(y[faceI]); } 31 / 33

slide-32
SLIDE 32

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Conclusion and Future Work

  • Two Hybrid PANS turbulence model has been implemented

and simulated on a test 2D backstep geometry.

  • Interesting results in the model behaviour are observed.

Improvement

  • Test the High Reynolds model on a 3D geometry.
  • Use an averaged fk from a previous RANS calculation.
  • Implement proper near wall functions for f and ζ parameters.
  • Test the Low Reynolds model on a 3D geometry.

32 / 33

slide-33
SLIDE 33

Outline Intro Turbulence models The test case PANS k − ε implementation k − ε − ζ − f Conclusion

Thank you

33 / 33