Commit 83422839 authored by henry's avatar henry
Browse files

Added DDES and IDDES SGS models for incompressible LES.

parent 4e4f9e0a
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "SpalartAllmarasDDES.H"
#include "addToRunTimeSelectionTable.H"
#include "wallDist.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
namespace LESModels
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(SpalartAllmarasDDES, 0);
addToRunTimeSelectionTable(LESModel, SpalartAllmarasDDES, dictionary);
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
tmp<volScalarField> SpalartAllmarasDDES::rd
(
const volScalarField& visc,
const volScalarField& S
) const
{
volScalarField d = wallDist(mesh_).y();
tmp<volScalarField> trd
(
new volScalarField
(
min
(
visc
/(
max
(
S,
dimensionedScalar("SMALL", S.dimensions(), SMALL)
)*sqr(kappa_*d)
+ dimensionedScalar
(
"ROOTVSMALL",
dimensionSet(0, 2 , -1, 0, 0),
ROOTVSMALL
)
), scalar(10.0)
)
)
);
return trd;
}
tmp<volScalarField> SpalartAllmarasDDES::fd(const volScalarField& S)
{
return 1.0 - tanh(pow3(8.0*rd(nuSgs_ + transport().nu(), S)));
}
void SpalartAllmarasDDES::dTildaUpdate(const volScalarField& S)
{
dTilda_ =
wallDist(mesh_).y()
- fd(S)*max
(
dimensionedScalar("zero", dimLength, 0.0),
wallDist(mesh_).y() - CDES_*delta()
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
SpalartAllmarasDDES::SpalartAllmarasDDES
(
const volVectorField& U,
const surfaceScalarField& phi,
transportModel& transport
)
:
SpalartAllmaras(U, phi, transport, typeName)
{}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace LESModels
} // End namespace incompressible
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::LESmodels::SpalartAllmarasDDES
Description
SpalartAllmaras DDES LES turbulence model for incompressible flows
Reference:
P.R. Spalart, S. Deck, S., M.L.Shur, K.D. Squires, M.Kh Strelets, and
A. Travin. `A new version of detached-eddy simulation, resistant to
ambiguous grid densities'. Theor. Comp. Fluid Dyn., 20:181-195, 2006.
SourceFiles
SpalartAllmarasDDES.C
\*---------------------------------------------------------------------------*/
#ifndef SpalartAllmarasDDES_H
#define SpalartAllmarasDDES_H
#include "SpalartAllmaras.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
namespace LESModels
{
/*---------------------------------------------------------------------------*\
Class SpalartAllmarasDDES Declaration
\*---------------------------------------------------------------------------*/
class SpalartAllmarasDDES
:
public SpalartAllmaras
{
// Private member functions
// Disallow default bitwise copy construct and assignment
SpalartAllmarasDDES(const SpalartAllmarasDDES&);
SpalartAllmarasDDES& operator=(const SpalartAllmarasDDES&);
protected:
// Protected member functions
tmp<volScalarField> fd(const volScalarField& S);
tmp<volScalarField> rd
(
const volScalarField& visc,
const volScalarField& S
) const;
//- Length scale calculation
virtual void dTildaUpdate(const volScalarField& S);
public:
//- Runtime type information
TypeName("SpalartAllmarasDDES");
// Constructors
//- Constructor from components
SpalartAllmarasDDES
(
const volVectorField& U,
const surfaceScalarField& phi,
transportModel& transport
);
//- Destructor
virtual ~SpalartAllmarasDDES()
{}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace LESModels
} // End namespace incompressible
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "IDDESDelta.H"
#include "addToRunTimeSelectionTable.H"
#include "wallDist.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(IDDESDelta, 0);
addToRunTimeSelectionTable(LESdelta, IDDESDelta, dictionary);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::IDDESDelta::calcDelta()
{
const Vector<label>& directions = mesh().directions();
label nD = (directions.nComponents + cmptSum(directions))/2;
// - Init hwn as wall distant.
volScalarField hwn = wallDist(mesh()).y();
scalar deltamaxTmp = 0.;
const cellList& cells = mesh().cells();
forAll(cells,cellI)
{
scalar deltaminTmp = 1.e10;
const labelList& cFaces = mesh().cells()[cellI];
const point& centrevector = mesh().cellCentres()[cellI];
forAll(cFaces, cFaceI)
{
label faceI = cFaces[cFaceI];
const point& facevector = mesh().faceCentres()[faceI];
scalar tmp = mag(facevector-centrevector);
if (tmp > deltamaxTmp)
{
deltamaxTmp = tmp;
}
if (tmp < deltaminTmp)
{
deltaminTmp = tmp;
}
}
hwn[cellI] = 2.0*deltaminTmp;
}
dimensionedScalar deltamax("deltamax",dimLength,2.0*deltamaxTmp);
if (nD == 3)
{
delta_.internalField() =
deltaCoeff_
*min(max(max(cw_*wallDist(mesh()).y(),cw_*deltamax),hwn),deltamax);
}
else if (nD == 2)
{
WarningIn("IDDESDelta::calcDelta()")
<< "Case is 2D, LES is not strictly applicable\n"
<< endl;
delta_.internalField() =
deltaCoeff_
*min(max(max(cw_*wallDist(mesh()).y(),cw_*deltamax),hwn),deltamax);
}
else
{
FatalErrorIn("IDDESDelta::calcDelta()")
<< "Case is not 3D or 2D, LES is not applicable"
<< exit(FatalError);
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::IDDESDelta::IDDESDelta
(
const word& name,
const fvMesh& mesh,
const dictionary& dd
)
:
LESdelta(name, mesh),
deltaCoeff_(readScalar(dd.subDict(type() + "Coeffs").lookup("deltaCoeff"))),
cw_(0)
{
dd.subDict(type() + "Coeffs").readIfPresent("cw", cw_);
calcDelta();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::IDDESDelta::read(const dictionary& dd)
{
dd.subDict(type() + "Coeffs").lookup("deltaCoeff") >> deltaCoeff_;
calcDelta();
}
void Foam::IDDESDelta::correct()
{
if (mesh_.changing())
{
calcDelta();
}
}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
Foam::IDDESDelta
Description
IDDESDelta used by the IDDES (improved low Re Spalart-Allmaras DES model)
The min and max delta are calculated using the double distance of the min or
max from the face centre to the cell centre.
SourceFiles
IDDESDelta.C
\*---------------------------------------------------------------------------*/
#ifndef IDDESDeltaDelta_H
#define IDDESDeltaDelta_H
#include "LESdelta.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class IDDESDelta Declaration
\*---------------------------------------------------------------------------*/
class IDDESDelta
:
public LESdelta
{
// Private data
scalar deltaCoeff_;
scalar cw_;
// Private Member Functions
//- Disallow default bitwise copy construct and assignment
IDDESDelta(const IDDESDelta&);
void operator=(const IDDESDelta&);
//- Calculate the delta values
void calcDelta();
public:
//- Runtime type information
TypeName("IDDESDelta");
// Constructors
//- Construct from name, mesh and IOdictionary
IDDESDelta
(
const word& name,
const fvMesh& mesh,
const dictionary&
);
// Destructor
~IDDESDelta()
{}
// Member Functions
//- Read the LESdelta dictionary
void read(const dictionary&);
// Correct values
void correct();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "SpalartAllmarasIDDES.H"
#include "addToRunTimeSelectionTable.H"
#include "wallDist.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
namespace LESModels
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(SpalartAllmarasIDDES, 0);
addToRunTimeSelectionTable(LESModel, SpalartAllmarasIDDES, dictionary);
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
tmp<volScalarField> SpalartAllmarasIDDES::alpha() const
{
return
0.25
- wallDist(mesh_).y()
/dimensionedScalar("hMax", dimLength, max(cmptMax(delta())));
}
tmp<volScalarField> SpalartAllmarasIDDES::ft(const volScalarField& S) const
{
return tanh(pow3(sqr(ct_)*r(nuSgs_, S)));
}
tmp<volScalarField> SpalartAllmarasIDDES::fl(const volScalarField& S) const
{
return tanh(pow(sqr(cl_)*r(transport().nu(), S), 10));
}
tmp<volScalarField> SpalartAllmarasIDDES::rd
(
const volScalarField& visc,
const volScalarField& S
) const
{
volScalarField d = wallDist(mesh_).y();
tmp<volScalarField> trd
(
new volScalarField
(
min
(
visc
/(
max
(
S,
dimensionedScalar("SMALL", S.dimensions(), SMALL)
)*sqr(kappa_*d)
+ dimensionedScalar
(
"ROOTVSMALL",
dimensionSet(0, 2 , -1, 0, 0),
ROOTVSMALL
)
), scalar(10.0)
)
)
);
return trd;