blendingFactorTemplates.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "gaussConvectionScheme.H"
27 #include "blendedSchemeBase.H"
28 #include "fvcCellReduce.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 template<class Type>
33 Foam::volScalarField& Foam::blendingFactor::factor
34 (
36 )
37 {
38  const word fieldName = "blendingFactor:" + field.name();
39 
40  if (!obr_.foundObject<volScalarField>(fieldName))
41  {
42  const fvMesh& mesh = refCast<const fvMesh>(obr_);
43 
44  volScalarField* factorPtr =
45  new volScalarField
46  (
47  IOobject
48  (
49  fieldName,
50  mesh.time().timeName(),
51  mesh,
52  IOobject::NO_READ,
53  IOobject::NO_WRITE
54  ),
55  mesh,
56  dimensionedScalar("0", dimless, 0.0),
57  zeroGradientFvPatchScalarField::typeName
58  );
59 
60  obr_.store(factorPtr);
61  }
62 
63  return
64  const_cast<volScalarField&>
65  (
66  obr_.lookupObject<volScalarField>(fieldName)
67  );
68 }
69 
70 
71 template<class Type>
72 void Foam::blendingFactor::calc()
73 {
75 
76  if (!obr_.foundObject<fieldType>(fieldName_))
77  {
78  return;
79  }
80 
81  const fvMesh& mesh = refCast<const fvMesh>(obr_);
82 
83  const fieldType& field = mesh.lookupObject<fieldType>(fieldName_);
84 
85  const word divScheme("div(" + phiName_ + ',' + fieldName_ + ')');
86  ITstream& its = mesh.divScheme(divScheme);
87 
88  const surfaceScalarField& phi =
89  mesh.lookupObject<surfaceScalarField>(phiName_);
90 
92  fv::convectionScheme<Type>::New(mesh, phi, its);
93 
95  refCast<const fv::gaussConvectionScheme<Type> >(cs());
96 
97  const surfaceInterpolationScheme<Type>& interpScheme =
98  gcs.interpScheme();
99 
100  if (!isA<blendedSchemeBase<Type> >(interpScheme))
101  {
102  FatalErrorIn("void Foam::blendingFactor::execute()")
103  << interpScheme.typeName << " is not a blended scheme"
104  << exit(FatalError);
105  }
106 
107  // retrieve the face-based blending factor
108  const blendedSchemeBase<Type>& blendedScheme =
109  refCast<const blendedSchemeBase<Type> >(interpScheme);
110  const surfaceScalarField factorf(blendedScheme.blendingFactor(field));
111 
112  // convert into vol field whose values represent the local face maxima
113  volScalarField& factor = this->factor(field);
114  factor = fvc::cellReduce(factorf, maxEqOp<scalar>());
115  factor.correctBoundaryConditions();
116 }
117 
118 
119 // ************************************************************************* //
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const surfaceInterpolationScheme< Type > & interpScheme() const
Input token stream.
Definition: ITstream.H:49
Construct a volume field from a surface field using a combine operator.
bool foundObject(const word &name) const
Is the named Type found?
ITstream & divScheme(const word &name) const
Definition: fvSchemes.C:421
tmp< GeometricField< Type, fvPatchField, volMesh > > cellReduce(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf, const CombineOp &cop)
Definition: fvcCellReduce.C:46
A class for handling words, derived from string.
Definition: word.H:59
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dynamicFvMesh & mesh
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
static tmp< convectionScheme< Type > > New(const fvMesh &mesh, const surfaceScalarField &faceFlux, Istream &schemeData)
Return a pointer to a new convectionScheme created on freestore.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
Abstract base class for surface interpolation schemes.
Basic second-order convection using face-gradients and Gauss&#39; theorem.
void store()
Transfer ownership of this object to its registry.
Definition: regIOobjectI.H:34
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const word & name() const
Return name.
Definition: IOobject.H:260
surfaceScalarField & phi
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Base class for blended schemes to provide access to the blending factor surface field.
error FatalError
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
virtual tmp< surfaceScalarField > blendingFactor(const GeometricField< Type, fvPatchField, volMesh > &vf) const =0
Return the face-based blending factor.
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.
Definition: PtrList.H:118