75#ifndef cellCoBlended_H
76#define cellCoBlended_H
80#include "surfaceInterpolate.H"
144 Co1_(readScalar(is)),
149 Co2_(readScalar(is)),
159 if (Co1_ < 0 || Co2_ < 0 || Co1_ >= Co2_)
162 <<
"coefficients = " << Co1_ <<
" and " << Co2_
163 <<
" should be > 0 and Co2 > Co1"
178 Co1_(readScalar(is)),
183 Co2_(readScalar(is)),
190 if (Co1_ < 0 || Co2_ < 0 || Co1_ >= Co2_)
193 <<
"coefficients = " << Co1_ <<
" and " << Co2_
194 <<
" should be > 0 and Co2 > Co1"
216 mesh.objectRegistry::template lookupObject<volScalarField>
224 <<
"dimensions of faceFlux are not correct"
238 extrapolatedCalculatedFvPatchScalarField::typeName
254 vf.
name() +
"BlendingFactor",
280 bf*tScheme1_().weights(vf)
281 + (scalar(1) - bf)*tScheme2_().weights(vf);
296 bf*tScheme1_().interpolate(vf)
297 + (scalar(1) - bf)*tScheme2_().interpolate(vf);
304 return tScheme1_().corrected() || tScheme2_().corrected();
const dimensionSet & dimensions() const
Return dimensions.
const Field< Type > & field() const
Return field.
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const word & name() const noexcept
Return the object name.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
scalar deltaTValue() const noexcept
Return time step value.
static word timeName(const scalar t, const int precision=precision_)
Base class for blended schemes to provide access to the blending factor surface field.
Two-scheme cell-based Courant number based blending differencing scheme.
virtual bool corrected() const
Return true if this scheme uses an explicit correction.
tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the interpolation weighting factors.
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > correction(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the explicit correction to the face-interpolate.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-interpolate of the given cell field.
virtual tmp< surfaceScalarField > blendingFactor(const GeometricField< Type, fvPatchField, volMesh > &vf) const
Return the face-based blending factor.
TypeName("cellCoBlended")
Runtime type information.
cellCoBlended(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &is)
Construct from mesh, faceFlux and Istream.
Mesh data needed to do the Finite Volume discretisation.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
const Time & time() const
Return the top-level database.
Abstract base class for surface interpolation schemes.
const fvMesh & mesh() const
Return mesh reference.
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
tmp< GeometricField< Type, fvPatchField, volMesh > > surfaceSum(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
const dimensionSet dimless
Dimensionless.
const dimensionSet dimArea(sqr(dimLength))
const dimensionSet dimVelocity
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimDensity
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.