Commit b7d3d572 authored by sergio's avatar sergio
Browse files

ENH: Atmospheric Boundary Layer (ABL) profiles for epsilon and velocity

parent e3a0ed25
......@@ -44,6 +44,8 @@ derivedFvPatchFields/turbulentHeatFluxTemperature/turbulentHeatFluxTemperatureFv
derivedFvPatchFields/turbulentMixingLengthDissipationRateInlet/turbulentMixingLengthDissipationRateInletFvPatchScalarField.C
derivedFvPatchFields/turbulentMixingLengthFrequencyInlet/turbulentMixingLengthFrequencyInletFvPatchScalarField.C
derivedFvPatchFields/fixedShearStress/fixedShearStressFvPatchVectorField.C
derivedFvPatchFields/ABLInletEpsilon/ABLInletEpsilonFvPatchScalarField.C
derivedFvPatchFields/ABLInletVelocity/ABLInletVelocityFvPatchVectorField.C
backwardsCompatibility/wallFunctions/backwardsCompatibilityWallFunctions.C
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 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 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "ABLInletEpsilonFvPatchScalarField.H"
#include "addToRunTimeSelectionTable.H"
#include "fvPatchFieldMapper.H"
#include "volFields.H"
#include "surfaceFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
ABLInletEpsilonFvPatchScalarField::ABLInletEpsilonFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(p, iF),
Ustar_(0),
z_(pTraits<vector>::zero),
z0_(0),
kappa_(0.41)
{}
ABLInletEpsilonFvPatchScalarField::ABLInletEpsilonFvPatchScalarField
(
const ABLInletEpsilonFvPatchScalarField& ptf,
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
fixedValueFvPatchScalarField(ptf, p, iF, mapper),
Ustar_(ptf.Ustar_),
z_(ptf.z_),
z0_(ptf.z0_),
kappa_(ptf.kappa_)
{}
ABLInletEpsilonFvPatchScalarField::ABLInletEpsilonFvPatchScalarField
(
const fvPatch& p,
const DimensionedField<scalar, volMesh>& iF,
const dictionary& dict
)
:
fixedValueFvPatchScalarField(p, iF),
Ustar_(readScalar(dict.lookup("Ustar"))),
z_(dict.lookup("z")),
z0_(readScalar(dict.lookup("z0"))),
kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41))
{
if (mag(z_) < SMALL)
{
FatalErrorIn("ABLInletEpsilonFvPatchScalarField(dict)")
<< "z is not correct"
<< abort(FatalError);
}
z_ /= mag(z_);
evaluate();
}
ABLInletEpsilonFvPatchScalarField::ABLInletEpsilonFvPatchScalarField
(
const ABLInletEpsilonFvPatchScalarField& fcvpvf,
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(fcvpvf, iF),
Ustar_(fcvpvf.Ustar_),
z_(fcvpvf.z_),
z0_(fcvpvf.z0_),
kappa_(fcvpvf.kappa_)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void ABLInletEpsilonFvPatchScalarField::updateCoeffs()
{
const vectorField& c = patch().Cf();
scalarField coord = (c & z_);
scalarField::operator=(pow(Ustar_, 3.0)/(kappa_*(coord + z0_)));
}
// Write
void ABLInletEpsilonFvPatchScalarField::write(Ostream& os) const
{
fvPatchScalarField::write(os);
os.writeKeyword("Ustar")
<< Ustar_ << token::END_STATEMENT << nl;
os.writeKeyword("z")
<< z_ << token::END_STATEMENT << nl;
os.writeKeyword("z0")
<< z0_ << token::END_STATEMENT << nl;
os.writeKeyword("kappa")
<< kappa_ << token::END_STATEMENT << nl;
writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
makePatchTypeField(fvPatchScalarField, ABLInletEpsilonFvPatchScalarField);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace incompressible
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 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 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
ABLInletEpsilonFvPatchScalarField
Description
Boundary condition specifies a epsilon inlet for the atmospheric boundary
layer (ABL). This boundaty is to be used in conjunction with
ABLInletVelocity.
@verbatim
epsilon = Ustar^3 / (K(z + z0))
where:
Ustar is the frictional velocity
K is karman's constant
z is the verical coordinate
z0 is the surface roughness lenght
@endverbatim
Reference:
D.M. Hargreaves and N.G. Wright
"On the use of the k-epsilon model in commercial CFD software to model the
neutral atmospheric boundary layer"
Journal of Wind Engineering and Industrial Aerodynamics 95(2007) 355-369.
SourceFiles
ABLInletEpsilonFvPatchScalarField.C
\*---------------------------------------------------------------------------*/
#ifndef ABLInletEpsilonFvPatchScalarField_H
#define ABLInletEpsilonFvPatchScalarField_H
#include "fvPatchFields.H"
#include "fixedValueFvPatchFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
/*---------------------------------------------------------------------------*\
Class ABLInletEpsilonFvPatchScalarField Declaration
\*---------------------------------------------------------------------------*/
class ABLInletEpsilonFvPatchScalarField
:
public fixedValueFvPatchScalarField
{
// Private data
//- Peak velocity magnitude
scalar Ustar_;
//- Direction of the z-coordinate
vector z_;
//- Surface roughness lenght
scalar z0_;
//- Von Karman constant
scalar kappa_;
public:
//- Runtime type information
TypeName("ABLInletEpsilon");
// Constructors
//- Construct from patch and internal field
ABLInletEpsilonFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&
);
//- Construct from patch, internal field and dictionary
ABLInletEpsilonFvPatchScalarField
(
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const dictionary&
);
//- Construct by mapping given ABLInletEpsilonFvPatchScalarField
// onto a new patch
ABLInletEpsilonFvPatchScalarField
(
const ABLInletEpsilonFvPatchScalarField&,
const fvPatch&,
const DimensionedField<scalar, volMesh>&,
const fvPatchFieldMapper&
);
//- Construct and return a clone
virtual tmp<fvPatchScalarField> clone() const
{
return tmp<fvPatchScalarField>
(
new ABLInletEpsilonFvPatchScalarField(*this)
);
}
//- Construct as copy setting internal field reference
ABLInletEpsilonFvPatchScalarField
(
const ABLInletEpsilonFvPatchScalarField&,
const DimensionedField<scalar, volMesh>&
);
//- Construct and return a clone setting internal field reference
virtual tmp<fvPatchScalarField> clone
(
const DimensionedField<scalar, volMesh>& iF
) const
{
return tmp<fvPatchScalarField>
(
new ABLInletEpsilonFvPatchScalarField(*this, iF)
);
}
// Member functions
//- Return max value
scalar& Ustar()
{
return Ustar_;
}
//- Return z direction
vector& z()
{
return z_;
}
//- Update coefficients
virtual void updateCoeffs();
//- Write
virtual void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace incompressible
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 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 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "ABLInletVelocityFvPatchVectorField.H"
#include "addToRunTimeSelectionTable.H"
#include "fvPatchFieldMapper.H"
#include "volFields.H"
#include "surfaceFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace incompressible
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
ABLInletVelocityFvPatchVectorField::ABLInletVelocityFvPatchVectorField
(
const fvPatch& p,
const DimensionedField<vector, volMesh>& iF
)
:
fixedValueFvPatchVectorField(p, iF),
Ustar_(0),
n_(pTraits<vector>::zero),
z_(pTraits<vector>::zero),
z0_(0),
kappa_(0.41)
{}
ABLInletVelocityFvPatchVectorField::ABLInletVelocityFvPatchVectorField
(
const ABLInletVelocityFvPatchVectorField& ptf,
const fvPatch& p,
const DimensionedField<vector, volMesh>& iF,
const fvPatchFieldMapper& mapper
)
:
fixedValueFvPatchVectorField(ptf, p, iF, mapper),
Ustar_(ptf.Ustar_),
n_(ptf.n_),
z_(ptf.z_),
z0_(ptf.z0_),
kappa_(ptf.kappa_)
{}
ABLInletVelocityFvPatchVectorField::ABLInletVelocityFvPatchVectorField
(
const fvPatch& p,
const DimensionedField<vector, volMesh>& iF,
const dictionary& dict
)
:
fixedValueFvPatchVectorField(p, iF),
Ustar_(readScalar(dict.lookup("Ustar"))),
n_(dict.lookup("n")),
z_(dict.lookup("z")),
z0_(readScalar(dict.lookup("z0"))),
kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41))
{
if (mag(n_) < SMALL || mag(z_) < SMALL)
{
FatalErrorIn("ABLInletVelocityFvPatchVectorField(dict)")
<< "n or z given with zero size not correct"
<< abort(FatalError);
}
n_ /= mag(n_);
z_ /= mag(z_);
evaluate();
}
ABLInletVelocityFvPatchVectorField::ABLInletVelocityFvPatchVectorField
(
const ABLInletVelocityFvPatchVectorField& fcvpvf,
const DimensionedField<vector, volMesh>& iF
)
:
fixedValueFvPatchVectorField(fcvpvf, iF),
Ustar_(fcvpvf.Ustar_),
n_(fcvpvf.n_),
z_(fcvpvf.z_),
z0_(fcvpvf.z0_),
kappa_(fcvpvf.kappa_)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void ABLInletVelocityFvPatchVectorField::updateCoeffs()
{
const vectorField& c = patch().Cf();
scalarField coord = (c & z_);
vectorField::operator=(n_*(Ustar_/kappa_)*log((coord + z0_)/z0_));
}
// Write
void ABLInletVelocityFvPatchVectorField::write(Ostream& os) const
{
fvPatchVectorField::write(os);
os.writeKeyword("Ustar")
<< Ustar_ << token::END_STATEMENT << nl;
os.writeKeyword("z0")
<< z0_ << token::END_STATEMENT << nl;
os.writeKeyword("n")
<< n_ << token::END_STATEMENT << nl;
os.writeKeyword("z")
<< z_ << token::END_STATEMENT << nl;
os.writeKeyword("kappa")
<< kappa_ << token::END_STATEMENT << nl;
writeEntry("value", os);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
makePatchTypeField(fvPatchVectorField, ABLInletVelocityFvPatchVectorField);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace incompressible
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 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 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, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
Class
ABLInletVelocityFvPatchVectorField
Description
Boundary condition specifies a atmospheric boundary layer (ABL)
velocity inlet profile given the friction velocity value,
flow direction n and direction of the parabolic coordinate z.
@verbatim
U = (Ustar/K) ln((z + z0)/z0)
where:
Ustar is the frictional velocity
K is karman's constant
z0 is the surface roughness lenght
z is the verical coordinate
and:
Ustar = K Uref/ln((Zref + z0)/z0)
where:
Uref is the reference velocity at Zref
Zref is the reference height.
@endverbatim
Reference:
D.M. Hargreaves and N.G. Wright
"On the use of the k-epsilon model in commercial CFD software to model the
neutral atmospheric boundary layer"
Journal of Wind Engineering and Industrial Aerodynamics 95(2007) 355-369.
NOTE: D.M. Hargreaves and N.G. Wright recommend Gamma epsilon in the k-epsilon
model should be changed from 1.3 to 1.11 for consistency.
The roughness height (Er) is given by Er = 20 z0 following the same
reference
SourceFiles
ABLInletVelocityFvPatchVectorField.C
\*---------------------------------------------------------------------------*/
#ifndef ABLInletVelocityFvPatchVectorField_H
#define ABLInletVelocityFvPatchVectorField_H
#include "fvPatchFields.H"
#include "fixedValueFvPatchFields.H"