Skip to content
Snippets Groups Projects
Commit 619ddc23 authored by Andrew Heather's avatar Andrew Heather
Browse files

INT: Integration of Upstream CFD's DeltaOmegaTilde delta function

- Initial code supplied by Marian Fuchs, Upstream CFD GmbH
- Code cleaned/refactored/integrated by OpenCFD
parent 4fc34c8a
Branches
Tags
1 merge request!560Integration of grey area turbulence models from Upstream CFD
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 Upstream CFD GmbH
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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 3 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, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "DeltaOmegaTildeDelta.H"
#include "fvcCurl.H"
#include "maxDeltaxyz.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace LESModels
{
defineTypeNameAndDebug(DeltaOmegaTildeDelta, 0);
addToRunTimeSelectionTable(LESdelta, DeltaOmegaTildeDelta, dictionary);
}
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::LESModels::DeltaOmegaTildeDelta::calcDelta()
{
const fvMesh& mesh = turbulenceModel_.mesh();
const volVectorField& U0 = turbulenceModel_.U();
const volVectorField vorticity(fvc::curl(U0));
const volVectorField nvecvort
(
vorticity
/ (
max
(
mag(vorticity),
dimensionedScalar("SMALL", dimless/dimTime, SMALL)
)
)
);
const cellList& cells = mesh.cells();
const vectorField& cellCentres = mesh.cellCentres();
const vectorField& faceCentres = mesh.faceCentres();
forAll(cells, celli)
{
const labelList& cFaces = cells[celli];
const point& cc = cellCentres[celli];
const vector& nv = nvecvort[celli];
scalar deltaMaxTmp = 0.0;
for (const label facei : cFaces)
{
const point& fc = faceCentres[facei];
scalar tmp = 2.0*mag(nv ^ (fc - cc));
if (tmp > deltaMaxTmp)
{
deltaMaxTmp = tmp;
}
}
delta_[celli] = deltaCoeff_*deltaMaxTmp;
}
const label nD = mesh.nGeometricD();
if (nD == 2)
{
WarningInFunction
<< "Case is 2D, LES is not strictly applicable" << nl
<< endl;
}
else if (nD != 3)
{
FatalErrorInFunction
<< "Case must be either 2D or 3D" << exit(FatalError);
}
// Handle coupled boundaries
delta_.correctBoundaryConditions();
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::LESModels::DeltaOmegaTildeDelta::DeltaOmegaTildeDelta
(
const word& name,
const turbulenceModel& turbulence,
const dictionary& dict
)
:
LESdelta(name, turbulence),
hmaxPtr_(nullptr),
deltaCoeff_
(
dict.optionalSubDict(type() + "Coeffs").getOrDefault<scalar>
(
"deltaCoeff",
1.035
)
),
requireUpdate_
(
dict.optionalSubDict(type() + "Coeffs").getOrDefault<bool>
(
"requireUpdate",
true
)
)
{
if (dict.optionalSubDict(type() + "Coeffs").found("hmax"))
{
// User-defined hmax
hmaxPtr_ =
LESdelta::New
(
IOobject::groupName("hmax", turbulence.U().group()),
turbulence,
dict.optionalSubDict("hmaxCoeffs"),
"hmax"
);
}
else
{
Info<< "Employing " << maxDeltaxyz::typeName << " for hmax" << endl;
hmaxPtr_.reset
(
new maxDeltaxyz
(
IOobject::groupName("hmax", turbulence.U().group()),
turbulence,
dict.optionalSubDict("hmaxCoeffs")
)
);
}
calcDelta();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::LESModels::DeltaOmegaTildeDelta::read(const dictionary& dict)
{
const dictionary& coeffsDict(dict.optionalSubDict(type() + "Coeffs"));
coeffsDict.readIfPresent<scalar>("deltaCoeff", deltaCoeff_);
coeffsDict.readIfPresent<bool>("requireUpdate", requireUpdate_);
calcDelta();
}
void Foam::LESModels::DeltaOmegaTildeDelta::correct()
{
if (turbulenceModel_.mesh().changing() && requireUpdate_)
{
hmaxPtr_->correct();
}
calcDelta();
}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2022 Upstream CFD GmbH
Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
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 3 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, see <http://www.gnu.org/licenses/>.
Class
Foam::DeltaOmegaTildeDelta
Description
Delta formulation that accounts for the orientation of the vorticity
vector. In "2D-regions" (i.e. xy-plane), delta is of the order
max(delta_x,delta_y), so that the influence of delta_z is discarded.
This can help to accelerate the transition from RANS to LES in hybrid
RANS/LES simulations.
Reference:
\verbatim
Shur, M. L., Spalart, P. R., Strelets, M. K., & Travin, A. K. (2015).
An enhanced version of DES with rapid transition from RANS to LES in
separated flows.
Flow, Turbulence and Combustion, 95, 709-737, 2015.
\endverbatim
SourceFiles
DeltaOmegaTildeDelta.C
\*---------------------------------------------------------------------------*/
#ifndef LESModels_DeltaOmegaTildeDelta_H
#define LESModels_DeltaOmegaTildeDelta_H
#include "LESdelta.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace LESModels
{
/*---------------------------------------------------------------------------*\
Class DeltaOmegaTildeDelta Declaration
\*---------------------------------------------------------------------------*/
class DeltaOmegaTildeDelta
:
public LESdelta
{
// Private data
//- Run-time selectable delta for hmax
// Defaults to the maxDeltaxyz model if not supplied
autoPtr<LESdelta> hmaxPtr_;
//- Model coefficient
scalar deltaCoeff_;
//- Flag to indicate whether hmax requires updating
bool requireUpdate_;
// Private Member Functions
//- No copy construct
DeltaOmegaTildeDelta(const DeltaOmegaTildeDelta&) = delete;
//- No copy assignment
void operator=(const DeltaOmegaTildeDelta&) = delete;
//- Calculate the delta values
void calcDelta();
public:
//- Runtime type information
TypeName("DeltaOmegaTilde");
// Constructors
//- Construct from name, turbulenceModel and dictionary
DeltaOmegaTildeDelta
(
const word& name,
const turbulenceModel& turbulence,
const dictionary&
);
//- Destructor
virtual ~DeltaOmegaTildeDelta() = default;
// Member Functions
//- Read the LESdelta dictionary
virtual void read(const dictionary&);
//- Return the hmax delta field
const volScalarField& hmax() const
{
return hmaxPtr_();
}
// Correct values
void correct();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace LESModels
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
......@@ -10,7 +10,7 @@ $(LESdelta)/smoothDelta/smoothDelta.C
$(LESdelta)/maxDeltaxyz/maxDeltaxyz.C
$(LESdelta)/IDDESDelta/IDDESDelta.C
$(LESdelta)/maxDeltaxyzCubeRootLESDelta/maxDeltaxyzCubeRootLESDelta.C
$(LESdelta)/DeltaOmegaTildeDelta/DeltaOmegaTildeDelta.C
LESfilters = LES/LESfilters
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment