Previous: Example of creating new Up: Example of creating new Next: Testing new turbulence model
This is an automatically generated documentation by LaTeX2HTML utility. In case of any issue, please, contact us at info@cfdsupport.com.
Turbulence model modification
- We modify the turbulence model in following way
- Production of turbulent kinetic energy will be controlled using new constant
- Modified model production term:
(24.8) - New model constant will be read from constant/turbulenceProperties
- Copy standard model e.g. to directory models in your work directory:
# cd $WM_PROJECT_USER_DIR
# mkdir -p src/turbulenceModels/incompressible/RAS
# cd $WM_PROJECT_USER_DIR/src/turbulenceModels/incompressible/RAS
# cp -r $FOAM_SRC/TurbulenceModels/turbulenceModels/RAS/kEpsilon/ ./ - New model can be called e.g. kEpsilonCprod model
- Rename model:
# mv kEpsilon kEpsilonCprod
- Delete redundant file:
# cd kEpsilonCprod
# rm kEpsilon.dep - Rename model source files:
# rename 's/kEpsilon/kEpsilonCprod/g' *
# rename kEpsilon kEpsilonCprod * - Replace name of the model also in files:
# sed -i 's/kEpsilon/kEpsilonCprod/g' *
- Create directory Make with configuration files:
# cd ..
# mkdir Make - In directory Make create file files:
kEpsilonCprod/kEpsilonCprod.C LIB = $(FOAM_USER_LIBBIN)/libuserincompressibleRASModels
- In directory Make create file options:
EXE_INC = \ -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \ -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ -I$(LIB_SRC)/transportModels \ -I$(LIB_SRC)/transportModels/incompressible/lnInclude \ -I$(LIB_SRC)/finiteVolume/lnInclude \ -I$(LIB_SRC)/meshTools/lnInclude LIB_LIBS = \ -lturbulenceModels \ -lincompressibleTurbulenceModels \ -lfiniteVolume \ -lmeshTools
- Now we have just copied and renamed the model
- In file kEpsilonCprod.H delete lines 211 -213
#ifdef NoRepository # include "kEpsilonCprod.C" #endif
- Compile model with no changes:
# wmake
- New library is created in dictionary $FOAM_USER_LIBBIN
# ll $FOAM_USER_LIBBIN
- Now implement the new model
- Constant is read the same way as other model constants
- Modify file kEpsilonCprod.H in the following way
- Create variable for production control:
dimensionedScalar Cmu_; dimensionedScalar C1_; dimensionedScalar C2_; dimensionedScalar C3_; dimensionedScalar sigmak_; dimensionedScalar sigmaEps_; dimensionedScalar Cprod_;
- Modify the file kEpsilonCprod.C in the following way:
- Create a new variable for
Cprod_ ( dimensioned<scalar>::lookupOrAddToDict ( "Cprod", this->coeffDict_, 1.0 ) ),
- Modify function for reading model constants:
bool kEpsilonCprod<BasicTurbulenceModel>::read() { if (eddyViscosity<RASModel<BasicTurbulenceModel>>::read()) { Cmu_.readIfPresent(this->coeffDict()); C1_.readIfPresent(this->coeffDict()); C2_.readIfPresent(this->coeffDict()); C3_.readIfPresent(this->coeffDict()); sigmak_.readIfPresent(this->coeffDict()); sigmaEps_.readIfPresent(this->coeffDict()); Cprod_.readIfPresent(this->coeffDict()); return true; } else { return false; } }
- Modify function for computing production of turbulent kinetic energy:
Cprod_*nut.v()*(dev(twoSymm(tgradU().v())) && tgradU().v())
- Add modified turbulence model into turbulence models library:
} // End namespace RASModels } // End namespace Foam #include "addToRunTimeSelectionTable.H" #include "makeTurbulenceModel.H" #include "RASModel.H" #include "transportModel.H" #include "incompressibleTurbulenceModel.H" #include "IncompressibleTurbulenceModel.H" namespace Foam { typedef IncompressibleTurbulenceModel<transportModel> transportModelIncompressibleTurbulenceModel; typedef RASModel<transportModelIncompressibleTurbulenceModel> RAStransportModelIncompressibleTurbulenceModel; } makeTemplatedTurbulenceModel(transportModelIncompressibleTurbulenceModel, RAS, kEpsilonCprod)
- Compile the new turbulence model:
# wmake