surfaceInterpolationScheme.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-2015 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 Description
25  Abstract base class for surface interpolation schemes.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 #include "volFields.H"
31 #include "surfaceFields.H"
32 #include "coupledFvPatchField.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
40 
41 // Return weighting factors for scheme given by name in dictionary
42 template<class Type>
43 tmp<surfaceInterpolationScheme<Type> > surfaceInterpolationScheme<Type>::New
44 (
45  const fvMesh& mesh,
46  Istream& schemeData
47 )
48 {
49  if (schemeData.eof())
50  {
52  (
53  "surfaceInterpolationScheme<Type>::New(const fvMesh&, Istream&)",
54  schemeData
55  ) << "Discretisation scheme not specified"
56  << endl << endl
57  << "Valid schemes are :" << endl
58  << MeshConstructorTablePtr_->sortedToc()
59  << exit(FatalIOError);
60  }
61 
62  const word schemeName(schemeData);
63 
64  if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
65  {
66  Info<< "surfaceInterpolationScheme<Type>::New"
67  "(const fvMesh&, Istream&)"
68  " : discretisation scheme = "
69  << schemeName
70  << endl;
71  }
72 
73  typename MeshConstructorTable::iterator constructorIter =
74  MeshConstructorTablePtr_->find(schemeName);
75 
76  if (constructorIter == MeshConstructorTablePtr_->end())
77  {
79  (
80  "surfaceInterpolationScheme<Type>::New(const fvMesh&, Istream&)",
81  schemeData
82  ) << "Unknown discretisation scheme "
83  << schemeName << nl << nl
84  << "Valid schemes are :" << endl
85  << MeshConstructorTablePtr_->sortedToc()
86  << exit(FatalIOError);
87  }
88 
89  return constructorIter()(mesh, schemeData);
90 }
91 
92 
93 // Return weighting factors for scheme given by name in dictionary
94 template<class Type>
96 (
97  const fvMesh& mesh,
98  const surfaceScalarField& faceFlux,
99  Istream& schemeData
100 )
101 {
102  if (schemeData.eof())
103  {
105  (
106  "surfaceInterpolationScheme<Type>::New"
107  "(const fvMesh&, const surfaceScalarField&, Istream&)",
108  schemeData
109  ) << "Discretisation scheme not specified"
110  << endl << endl
111  << "Valid schemes are :" << endl
112  << MeshConstructorTablePtr_->sortedToc()
113  << exit(FatalIOError);
114  }
115 
116  const word schemeName(schemeData);
117 
118  if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
119  {
120  Info<< "surfaceInterpolationScheme<Type>::New"
121  "(const fvMesh&, const surfaceScalarField&, Istream&)"
122  " : discretisation scheme = "
123  << schemeName
124  << endl;
125  }
126 
127  typename MeshFluxConstructorTable::iterator constructorIter =
128  MeshFluxConstructorTablePtr_->find(schemeName);
129 
130  if (constructorIter == MeshFluxConstructorTablePtr_->end())
131  {
133  (
134  "surfaceInterpolationScheme<Type>::New"
135  "(const fvMesh&, const surfaceScalarField&, Istream&)",
136  schemeData
137  ) << "Unknown discretisation scheme "
138  << schemeName << nl << nl
139  << "Valid schemes are :" << endl
140  << MeshFluxConstructorTablePtr_->sortedToc()
141  << exit(FatalIOError);
142  }
143 
144  return constructorIter()(mesh, faceFlux, schemeData);
145 }
146 
147 
148 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
149 
150 template<class Type>
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
157 //- Return the face-interpolate of the given cell field
158 // with the given owner and neighbour weighting factors
159 template<class Type>
162 (
164  const tmp<surfaceScalarField>& tlambdas,
165  const tmp<surfaceScalarField>& tys
166 )
167 {
168  if (surfaceInterpolation::debug)
169  {
170  Info<< "surfaceInterpolationScheme<Type>::uncorrectedInterpolate"
171  "(const GeometricField<Type, fvPatchField, volMesh>&, "
172  "const tmp<surfaceScalarField>&, "
173  "const tmp<surfaceScalarField>&) : "
174  "interpolating "
175  << vf.type() << " "
176  << vf.name()
177  << " from cells to faces "
178  "without explicit correction"
179  << endl;
180  }
181 
182  const surfaceScalarField& lambdas = tlambdas();
183  const surfaceScalarField& ys = tys();
184 
185  const Field<Type>& vfi = vf.internalField();
186  const scalarField& lambda = lambdas.internalField();
187  const scalarField& y = ys.internalField();
188 
189  const fvMesh& mesh = vf.mesh();
190  const labelUList& P = mesh.owner();
191  const labelUList& N = mesh.neighbour();
192 
194  (
196  (
197  IOobject
198  (
199  "interpolate("+vf.name()+')',
200  vf.instance(),
201  vf.db()
202  ),
203  mesh,
204  vf.dimensions()
205  )
206  );
208 
209  Field<Type>& sfi = sf.internalField();
210 
211  for (label fi=0; fi<P.size(); fi++)
212  {
213  sfi[fi] = lambda[fi]*vfi[P[fi]] + y[fi]*vfi[N[fi]];
214  }
215 
216 
217  // Interpolate across coupled patches using given lambdas and ys
218 
219  forAll(lambdas.boundaryField(), pi)
220  {
221  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
222  const fvsPatchScalarField& pY = ys.boundaryField()[pi];
223 
224  if (vf.boundaryField()[pi].coupled())
225  {
226  sf.boundaryField()[pi] =
227  pLambda*vf.boundaryField()[pi].patchInternalField()
228  + pY*vf.boundaryField()[pi].patchNeighbourField();
229  }
230  else
231  {
232  sf.boundaryField()[pi] = vf.boundaryField()[pi];
233  }
234  }
235 
236  tlambdas.clear();
237  tys.clear();
238 
239  return tsf;
240 }
241 
242 
243 //- Return the face-interpolate of the given cell field
244 // with the given weighting factors
245 template<class Type>
248 (
250  const tmp<surfaceScalarField>& tlambdas
251 )
252 {
253  if (surfaceInterpolation::debug)
254  {
255  Info<< "surfaceInterpolationScheme<Type>::interpolate"
256  "(const GeometricField<Type, fvPatchField, volMesh>&, "
257  "const tmp<surfaceScalarField>&) : "
258  "interpolating "
259  << vf.type() << " "
260  << vf.name()
261  << " from cells to faces "
262  "without explicit correction"
263  << endl;
264  }
265 
266  const surfaceScalarField& lambdas = tlambdas();
267 
268  const Field<Type>& vfi = vf.internalField();
269  const scalarField& lambda = lambdas.internalField();
270 
271  const fvMesh& mesh = vf.mesh();
272  const labelUList& P = mesh.owner();
273  const labelUList& N = mesh.neighbour();
274 
276  (
278  (
279  IOobject
280  (
281  "interpolate("+vf.name()+')',
282  vf.instance(),
283  vf.db()
284  ),
285  mesh,
286  vf.dimensions()
287  )
288  );
290 
291  Field<Type>& sfi = sf.internalField();
292 
293  for (label fi=0; fi<P.size(); fi++)
294  {
295  sfi[fi] = lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]];
296  }
297 
298  // Interpolate across coupled patches using given lambdas
299 
300  forAll(lambdas.boundaryField(), pi)
301  {
302  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
303 
304  if (vf.boundaryField()[pi].coupled())
305  {
306  tsf().boundaryField()[pi] =
307  pLambda*vf.boundaryField()[pi].patchInternalField()
308  + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField();
309  }
310  else
311  {
312  sf.boundaryField()[pi] = vf.boundaryField()[pi];
313  }
314  }
315 
316  tlambdas.clear();
317 
318  return tsf;
319 }
320 
321 
322 //- Return the face-interpolate of the given cell field
323 // with explicit correction
324 template<class Type>
327 (
329 ) const
330 {
331  if (surfaceInterpolation::debug)
332  {
333  Info<< "surfaceInterpolationScheme<Type>::interpolate"
334  "(const GeometricField<Type, fvPatchField, volMesh>&) : "
335  "interpolating "
336  << vf.type() << " "
337  << vf.name()
338  << " from cells to faces"
339  << endl;
340  }
341 
343  = interpolate(vf, weights(vf));
344 
345  if (corrected())
346  {
347  tsf() += correction(vf);
348  }
349 
350  return tsf;
351 }
352 
353 
354 //- Return the face-interpolate of the given cell field
355 // with explicit correction
356 template<class Type>
359 (
361 ) const
362 {
364  = interpolate(tvf());
365  tvf.clear();
366  return tinterpVf;
367 }
368 
369 
370 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
371 
372 } // End namespace Foam
373 
374 // ************************************************************************* //
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Foam::surfaceFields.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:174
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
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const fileName & instance() const
Definition: IOobject.H:337
InternalField & internalField()
Return internal field.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:282
messageStream Info
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &, const tmp< surfaceScalarField > &, const tmp< surfaceScalarField > &)
Return the face-interpolate of the given cell field.
const Mesh & mesh() const
Return mesh.
dynamicFvMesh & mesh
volScalarField sf(fieldObject, mesh)
Namespace for OpenFOAM.
bool eof() const
Return true if end of input seen.
Definition: IOstream.H:339
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOerror FatalIOError
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
Abstract base class for surface interpolation schemes.
label N
Definition: createFields.H:22
#define forAll(list, i)
Definition: UList.H:421
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
const dimensionSet & dimensions() const
Return dimensions.
Pre-declare SubField and related Field type.
Definition: Field.H:57
const word & name() const
Return name.
Definition: IOobject.H:260
scalar y
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:288
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
#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
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65