limitedSurfaceInterpolationScheme.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) 2011 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 
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "coupledFvPatchField.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
37 
38 template<class Type>
39 tmp<limitedSurfaceInterpolationScheme<Type> >
41 (
42  const fvMesh& mesh,
43  Istream& schemeData
44 )
45 {
46  if (surfaceInterpolation::debug)
47  {
48  Info<< "limitedSurfaceInterpolationScheme<Type>::"
49  "New(const fvMesh&, Istream&)"
50  " : constructing limitedSurfaceInterpolationScheme<Type>"
51  << endl;
52  }
53 
54  if (schemeData.eof())
55  {
57  (
58  "limitedSurfaceInterpolationScheme<Type>::"
59  "New(const fvMesh&, Istream&)",
60  schemeData
61  ) << "Discretisation scheme not specified"
62  << endl << endl
63  << "Valid schemes are :" << endl
64  << MeshConstructorTablePtr_->sortedToc()
65  << exit(FatalIOError);
66  }
67 
68  const word schemeName(schemeData);
69 
70  typename MeshConstructorTable::iterator constructorIter =
71  MeshConstructorTablePtr_->find(schemeName);
72 
73  if (constructorIter == MeshConstructorTablePtr_->end())
74  {
76  (
77  "limitedSurfaceInterpolationScheme<Type>::"
78  "New(const fvMesh&, Istream&)",
79  schemeData
80  ) << "Unknown discretisation scheme "
81  << schemeName << nl << nl
82  << "Valid schemes are :" << endl
83  << MeshConstructorTablePtr_->sortedToc()
84  << exit(FatalIOError);
85  }
86 
87  return constructorIter()(mesh, schemeData);
88 }
89 
90 
91 // Return weighting factors for scheme given by name in dictionary
92 template<class Type>
95 (
96  const fvMesh& mesh,
97  const surfaceScalarField& faceFlux,
98  Istream& schemeData
99 )
100 {
101  if (surfaceInterpolation::debug)
102  {
103  Info<< "limitedSurfaceInterpolationScheme<Type>::New"
104  "(const fvMesh&, const surfaceScalarField&, Istream&) : "
105  "constructing limitedSurfaceInterpolationScheme<Type>"
106  << endl;
107  }
108 
109  if (schemeData.eof())
110  {
112  (
113  "limitedSurfaceInterpolationScheme<Type>::New"
114  "(const fvMesh&, const surfaceScalarField&, Istream&)",
115  schemeData
116  ) << "Discretisation scheme not specified"
117  << endl << endl
118  << "Valid schemes are :" << endl
119  << MeshConstructorTablePtr_->sortedToc()
120  << exit(FatalIOError);
121  }
122 
123  const word schemeName(schemeData);
124 
125  typename MeshFluxConstructorTable::iterator constructorIter =
126  MeshFluxConstructorTablePtr_->find(schemeName);
127 
128  if (constructorIter == MeshFluxConstructorTablePtr_->end())
129  {
131  (
132  "limitedSurfaceInterpolationScheme<Type>::New"
133  "(const fvMesh&, const surfaceScalarField&, Istream&)",
134  schemeData
135  ) << "Unknown discretisation scheme "
136  << schemeName << nl << nl
137  << "Valid schemes are :" << endl
138  << MeshFluxConstructorTablePtr_->sortedToc()
139  << exit(FatalIOError);
140  }
141 
142  return constructorIter()(mesh, faceFlux, schemeData);
143 }
144 
145 
146 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
147 
148 template<class Type>
150 {}
151 
152 
153 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
154 
155 template<class Type>
157 (
159  const surfaceScalarField& CDweights,
160  tmp<surfaceScalarField> tLimiter
161 ) const
162 {
163  // Note that here the weights field is initialised as the limiter
164  // from which the weight is calculated using the limiter value
165  surfaceScalarField& Weights = tLimiter();
166 
167  scalarField& pWeights = Weights.internalField();
168 
169  forAll(pWeights, face)
170  {
171  pWeights[face] =
172  pWeights[face]*CDweights[face]
173  + (1.0 - pWeights[face])*pos(faceFlux_[face]);
174  }
175 
176  surfaceScalarField::GeometricBoundaryField& bWeights =
177  Weights.boundaryField();
178 
179  forAll(bWeights, patchI)
180  {
181  scalarField& pWeights = bWeights[patchI];
182 
183  const scalarField& pCDweights = CDweights.boundaryField()[patchI];
184  const scalarField& pFaceFlux = faceFlux_.boundaryField()[patchI];
185 
186  forAll(pWeights, face)
187  {
188  pWeights[face] =
189  pWeights[face]*pCDweights[face]
190  + (1.0 - pWeights[face])*pos(pFaceFlux[face]);
191  }
192  }
193 
194  return tLimiter;
195 }
196 
197 template<class Type>
199 (
201 ) const
202 {
203  return this->weights
204  (
205  phi,
207  this->limiter(phi)
208  );
209 }
210 
211 template<class Type>
214 (
216 ) const
217 {
218  return faceFlux_*this->interpolate(phi);
219 }
220 
221 
222 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
223 
224 } // End namespace Foam
225 
226 // ************************************************************************* //
const surfaceScalarField & weights() const
Return reference to linear difference weighting factors.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
A class for handling words, derived from string.
Definition: word.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
InternalField & internalField()
Return internal field.
messageStream Info
dynamicFvMesh & mesh
tmp< surfaceScalarField > weights(const GeometricField< Type, fvPatchField, volMesh > &, const surfaceScalarField &CDweights, tmp< surfaceScalarField > tLimiter) const
Return the interpolation weighting factors for the given field,.
Namespace for OpenFOAM.
bool eof() const
Return true if end of input seen.
Definition: IOstream.H:339
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOerror FatalIOError
#define forAll(list, i)
Definition: UList.H:421
virtual tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > flux(const GeometricField< Type, fvPatchField, volMesh > &) const
Return the interpolation weighting factors.
dimensionedScalar pos(const dimensionedScalar &ds)
static tmp< limitedSurfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
void limiter(scalarField &allLambda, const RdeltaTType &rDeltaT, const RhoType &rho, const volScalarField &psi, const surfaceScalarField &phiBD, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su, const scalar psiMax, const scalar psiMin)
#define FatalIOErrorIn(functionName, ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:325
A class for managing temporary objects.
Definition: PtrList.H:118