surfaceInterpolationScheme.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 "geometricOneField.H"
30 #include "coupledFvPatchField.H"
31 
32 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
33 
34 template<class Type>
37 (
38  const fvMesh& mesh,
39  Istream& schemeData
40 )
41 {
42  if (schemeData.eof())
43  {
45  (
46  schemeData
47  ) << "Discretisation scheme not specified"
48  << endl << endl
49  << "Valid schemes are :" << endl
50  << MeshConstructorTablePtr_->sortedToc()
51  << exit(FatalIOError);
52  }
53 
54  const word schemeName(schemeData);
55 
56  if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
57  {
58  InfoInFunction << "Discretisation scheme = " << schemeName << endl;
59  }
60 
61  typename MeshConstructorTable::iterator constructorIter =
62  MeshConstructorTablePtr_->find(schemeName);
63 
64  if (constructorIter == MeshConstructorTablePtr_->end())
65  {
67  (
68  schemeData
69  ) << "Unknown discretisation scheme "
70  << schemeName << nl << nl
71  << "Valid schemes are :" << endl
72  << MeshConstructorTablePtr_->sortedToc()
73  << exit(FatalIOError);
74  }
75 
76  return constructorIter()(mesh, schemeData);
77 }
78 
79 
80 template<class Type>
83 (
84  const fvMesh& mesh,
85  const surfaceScalarField& faceFlux,
86  Istream& schemeData
87 )
88 {
89  if (schemeData.eof())
90  {
92  (
93  schemeData
94  ) << "Discretisation scheme not specified"
95  << endl << endl
96  << "Valid schemes are :" << endl
97  << MeshConstructorTablePtr_->sortedToc()
98  << exit(FatalIOError);
99  }
100 
101  const word schemeName(schemeData);
102 
103  if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
104  {
106  << "Discretisation scheme = " << schemeName << endl;
107  }
108 
109  typename MeshFluxConstructorTable::iterator constructorIter =
110  MeshFluxConstructorTablePtr_->find(schemeName);
111 
112  if (constructorIter == MeshFluxConstructorTablePtr_->end())
113  {
115  (
116  schemeData
117  ) << "Unknown discretisation scheme "
118  << schemeName << nl << nl
119  << "Valid schemes are :" << endl
120  << MeshFluxConstructorTablePtr_->sortedToc()
121  << exit(FatalIOError);
122  }
123 
124  return constructorIter()(mesh, faceFlux, schemeData);
125 }
126 
127 
128 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
129 
130 template<class Type>
132 {}
133 
134 
135 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
136 
137 template<class Type>
140 (
142  const tmp<surfaceScalarField>& tlambdas,
143  const tmp<surfaceScalarField>& tys
144 )
145 {
146  if (surfaceInterpolation::debug)
147  {
149  << "Interpolating "
150  << vf.type() << " "
151  << vf.name()
152  << " from cells to faces "
153  "without explicit correction"
154  << endl;
155  }
156 
157  const surfaceScalarField& lambdas = tlambdas();
158  const surfaceScalarField& ys = tys();
159 
160  const Field<Type>& vfi = vf;
161  const scalarField& lambda = lambdas;
162  const scalarField& y = ys;
163 
164  const fvMesh& mesh = vf.mesh();
165  const labelUList& P = mesh.owner();
166  const labelUList& N = mesh.neighbour();
167 
169  (
171  (
172  IOobject
173  (
174  "interpolate("+vf.name()+')',
175  vf.instance(),
176  vf.db()
177  ),
178  mesh,
179  vf.dimensions()
180  )
181  );
183 
184  Field<Type>& sfi = sf.primitiveFieldRef();
185 
186  for (label fi=0; fi<P.size(); fi++)
187  {
188  sfi[fi] = lambda[fi]*vfi[P[fi]] + y[fi]*vfi[N[fi]];
189  }
190 
191 
192  // Interpolate across coupled patches using given lambdas and ys
194  Boundary& sfbf = sf.boundaryFieldRef();
195 
196  forAll(lambdas.boundaryField(), pi)
197  {
198  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
199  const fvsPatchScalarField& pY = ys.boundaryField()[pi];
200 
201  if (vf.boundaryField()[pi].coupled())
202  {
203  sfbf[pi] =
204  pLambda*vf.boundaryField()[pi].patchInternalField()
205  + pY*vf.boundaryField()[pi].patchNeighbourField();
206  }
207  else
208  {
209  sfbf[pi] = vf.boundaryField()[pi];
210  }
211  }
212 
213  tlambdas.clear();
214  tys.clear();
215 
216  return tsf;
217 }
218 
219 
220 template<class Type>
221 template<class SFType>
222 Foam::tmp
223 <
225  <
226  typename Foam::innerProduct<typename SFType::value_type, Type>::type,
227  Foam::fvsPatchField,
229  >
230 >
232 (
233  const SFType& Sf,
235  const tmp<surfaceScalarField>& tlambdas
236 )
237 {
238  if (surfaceInterpolation::debug)
239  {
241  << "Interpolating "
242  << vf.type() << " "
243  << vf.name()
244  << " from cells to faces "
245  "without explicit correction"
246  << endl;
247  }
248 
249  typedef typename Foam::innerProduct<typename SFType::value_type, Type>::type
250  RetType;
251 
252  const surfaceScalarField& lambdas = tlambdas();
253 
254  const Field<Type>& vfi = vf;
255  const scalarField& lambda = lambdas;
256 
257  const fvMesh& mesh = vf.mesh();
258  const labelUList& P = mesh.owner();
259  const labelUList& N = mesh.neighbour();
260 
262  (
264  (
265  IOobject
266  (
267  "interpolate("+vf.name()+')',
268  vf.instance(),
269  vf.db()
270  ),
271  mesh,
272  Sf.dimensions()*vf.dimensions()
273  )
274  );
276 
277  Field<RetType>& sfi = sf.primitiveFieldRef();
278 
279  const typename SFType::Internal& Sfi = Sf();
280 
281  for (label fi=0; fi<P.size(); fi++)
282  {
283  sfi[fi] = Sfi[fi] & (lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]]);
284  }
285 
286  // Interpolate across coupled patches using given lambdas
287 
289  Boundary& sfbf = sf.boundaryFieldRef();
290 
291  forAll(lambdas.boundaryField(), pi)
292  {
293  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
294  const typename SFType::Patch& pSf = Sf.boundaryField()[pi];
295  fvsPatchField<RetType>& psf = sfbf[pi];
296 
297  if (vf.boundaryField()[pi].coupled())
298  {
299  psf =
300  pSf
301  & (
302  pLambda*vf.boundaryField()[pi].patchInternalField()
303  + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField()
304  );
305  }
306  else
307  {
308  psf = pSf & vf.boundaryField()[pi];
309  }
310  }
311 
312  tlambdas.clear();
313 
314  return tsf;
315 }
316 
317 
318 template<class Type>
321 (
323  const tmp<surfaceScalarField>& tlambdas
324 )
325 {
326  return dotInterpolate(geometricOneField(), vf, tlambdas);
327 }
328 
329 
330 template<class Type>
331 Foam::tmp
332 <
334  <
335  typename Foam::innerProduct<Foam::vector, Type>::type,
336  Foam::fvsPatchField,
338  >
339 >
341 (
342  const surfaceVectorField& Sf,
344 ) const
345 {
346  if (surfaceInterpolation::debug)
347  {
349  << "Interpolating "
350  << vf.type() << " "
351  << vf.name()
352  << " from cells to faces"
353  << endl;
354  }
355 
356  tmp
357  <
359  <
360  typename Foam::innerProduct<Foam::vector, Type>::type,
363  >
364  > tsf = dotInterpolate(Sf, vf, weights(vf));
365 
366  if (corrected())
367  {
368  tsf.ref() += Sf & correction(vf);
369  }
370 
371  return tsf;
372 }
373 
374 
375 template<class Type>
376 Foam::tmp
377 <
379  <
380  typename Foam::innerProduct<Foam::vector, Type>::type,
381  Foam::fvsPatchField,
383  >
384 >
386 (
387  const surfaceVectorField& Sf,
389 ) const
390 {
391  tmp
392  <
394  <
395  typename Foam::innerProduct<Foam::vector, Type>::type,
398  >
399  > tSfDotinterpVf = dotInterpolate(Sf, tvf());
400  tvf.clear();
401  return tSfDotinterpVf;
402 }
403 
404 
405 template<class Type>
408 (
410 ) const
411 {
412  if (surfaceInterpolation::debug)
413  {
415  << "Interpolating "
416  << vf.type() << " "
417  << vf.name()
418  << " from cells to faces"
419  << endl;
420  }
421 
423  = interpolate(vf, weights(vf));
424 
425  if (corrected())
426  {
427  tsf.ref() += correction(vf);
428  }
429 
430  return tsf;
431 }
432 
433 
434 template<class Type>
437 (
439 ) const
440 {
442  = interpolate(tvf());
443  tvf.clear();
444  return tinterpVf;
445 }
446 
447 
448 // ************************************************************************* //
Foam::surfaceFields.
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
const word & name() const
Return name.
Definition: IOobject.H:297
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:230
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Mesh data needed to do the Finite Volume discretisation.
Definition: surfaceMesh.H:47
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:288
scalar y
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Pre-declare SubField and related Field type.
Definition: Field.H:57
A class for handling words, derived from string.
Definition: word.H:59
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
bool eof() const
Return true if end of input seen.
Definition: IOstream.H:339
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:61
static tmp< GeometricField< typename innerProduct< typename SFType::value_type, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const SFType &Sf, const GeometricField< Type, fvPatchField, volMesh > &vf, const tmp< surfaceScalarField > &tlambdas)
Return the face-interpolate of the given cell field.
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static const char nl
Definition: Ostream.H:265
const Mesh & mesh() const
Return mesh.
volScalarField sf(fieldObject, mesh)
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.
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.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:282
const fileName & instance() const
Definition: IOobject.H:392
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
label N
Definition: createFields.H:22
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:361
label size() const
Return the number of elements in the UList.
Definition: UListI.H:299
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Abstract base class for surface interpolation schemes.
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
IOerror FatalIOError
#define InfoInFunction
Report an information message using Foam::Info.