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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2019 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "surfaceInterpolationScheme.H"
30 #include "volFields.H"
31 #include "surfaceFields.H"
32 #include "geometricOneField.H"
33 #include "coupledFvPatchField.H"
34 
35 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
36 
37 template<class Type>
38 Foam::tmp<Foam::surfaceInterpolationScheme<Type>>
39 Foam::surfaceInterpolationScheme<Type>::New
40 (
41  const fvMesh& mesh,
42  Istream& schemeData
43 )
44 {
45  if (schemeData.eof())
46  {
47  FatalIOErrorInFunction(schemeData)
48  << "Discretisation scheme not specified\n\n"
49  << "Valid schemes:\n"
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  auto cstrIter = MeshConstructorTablePtr_->cfind(schemeName);
62 
63  if (!cstrIter.found())
64  {
65  FatalIOErrorInLookup
66  (
67  schemeData,
68  "discretisation",
69  schemeName,
70  *MeshConstructorTablePtr_
71  ) << exit(FatalError);
72  }
73 
74  return cstrIter()(mesh, schemeData);
75 }
76 
77 
78 template<class Type>
79 Foam::tmp<Foam::surfaceInterpolationScheme<Type>>
80 Foam::surfaceInterpolationScheme<Type>::New
81 (
82  const fvMesh& mesh,
83  const surfaceScalarField& faceFlux,
84  Istream& schemeData
85 )
86 {
87  if (schemeData.eof())
88  {
89  FatalIOErrorInFunction(schemeData)
90  << "Discretisation scheme not specified"
91  << endl << endl
92  << "Valid schemes are :" << endl
93  << MeshConstructorTablePtr_->sortedToc()
94  << exit(FatalIOError);
95  }
96 
97  const word schemeName(schemeData);
98 
99  if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
100  {
101  InfoInFunction
102  << "Discretisation scheme = " << schemeName << endl;
103  }
104 
105  auto cstrIter = MeshFluxConstructorTablePtr_->cfind(schemeName);
106 
107  if (!cstrIter.found())
108  {
109  FatalIOErrorInLookup
110  (
111  schemeData,
112  "discretisation",
113  schemeName,
114  *MeshFluxConstructorTablePtr_
115  ) << exit(FatalIOError);
116  }
117 
118  return cstrIter()(mesh, faceFlux, schemeData);
119 }
120 
121 
122 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
123 
124 template<class Type>
125 Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
126 Foam::surfaceInterpolationScheme<Type>::interpolate
127 (
128  const GeometricField<Type, fvPatchField, volMesh>& vf,
129  const tmp<surfaceScalarField>& tlambdas,
130  const tmp<surfaceScalarField>& tys
131 )
132 {
133  if (surfaceInterpolation::debug)
134  {
135  InfoInFunction
136  << "Interpolating "
137  << vf.type() << " "
138  << vf.name()
139  << " from cells to faces without explicit correction"
140  << endl;
141  }
142 
143  const surfaceScalarField& lambdas = tlambdas();
144  const surfaceScalarField& ys = tys();
145 
146  const Field<Type>& vfi = vf;
147  const scalarField& lambda = lambdas;
148  const scalarField& y = ys;
149 
150  const fvMesh& mesh = vf.mesh();
151  const labelUList& P = mesh.owner();
152  const labelUList& N = mesh.neighbour();
153 
154  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf
155  (
156  new GeometricField<Type, fvsPatchField, surfaceMesh>
157  (
158  IOobject
159  (
160  "interpolate("+vf.name()+')',
161  vf.instance(),
162  vf.db()
163  ),
164  mesh,
165  vf.dimensions()
166  )
167  );
168  GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsf.ref();
169 
170  Field<Type>& sfi = sf.primitiveFieldRef();
171 
172  for (label fi=0; fi<P.size(); fi++)
173  {
174  sfi[fi] = lambda[fi]*vfi[P[fi]] + y[fi]*vfi[N[fi]];
175  }
176 
177 
178  // Interpolate across coupled patches using given lambdas and ys
179  typename GeometricField<Type, fvsPatchField, surfaceMesh>::
180  Boundary& sfbf = sf.boundaryFieldRef();
181 
182  forAll(lambdas.boundaryField(), pi)
183  {
184  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
185  const fvsPatchScalarField& pY = ys.boundaryField()[pi];
186 
187  if (vf.boundaryField()[pi].coupled())
188  {
189  sfbf[pi] =
190  pLambda*vf.boundaryField()[pi].patchInternalField()
191  + pY*vf.boundaryField()[pi].patchNeighbourField();
192  }
193  else
194  {
195  sfbf[pi] = vf.boundaryField()[pi];
196  }
197  }
198 
199  tlambdas.clear();
200  tys.clear();
201 
202  return tsf;
203 }
204 
205 
206 template<class Type>
207 template<class SFType>
208 Foam::tmp
209 <
210  Foam::GeometricField
211  <
212  typename Foam::innerProduct<typename SFType::value_type, Type>::type,
213  Foam::fvsPatchField,
214  Foam::surfaceMesh
215  >
216 >
217 Foam::surfaceInterpolationScheme<Type>::dotInterpolate
218 (
219  const SFType& Sf,
220  const GeometricField<Type, fvPatchField, volMesh>& vf,
221  const tmp<surfaceScalarField>& tlambdas
222 )
223 {
224  if (surfaceInterpolation::debug)
225  {
226  InfoInFunction
227  << "Interpolating "
228  << vf.type() << " "
229  << vf.name()
230  << " from cells to faces without explicit correction"
231  << endl;
232  }
233 
234  typedef typename Foam::innerProduct<typename SFType::value_type, Type>::type
235  RetType;
236 
237  const surfaceScalarField& lambdas = tlambdas();
238 
239  const Field<Type>& vfi = vf;
240  const scalarField& lambda = lambdas;
241 
242  const fvMesh& mesh = vf.mesh();
243  const labelUList& P = mesh.owner();
244  const labelUList& N = mesh.neighbour();
245 
246  tmp<GeometricField<RetType, fvsPatchField, surfaceMesh>> tsf
247  (
248  new GeometricField<RetType, fvsPatchField, surfaceMesh>
249  (
250  IOobject
251  (
252  "interpolate("+vf.name()+')',
253  vf.instance(),
254  vf.db()
255  ),
256  mesh,
257  Sf.dimensions()*vf.dimensions()
258  )
259  );
260  GeometricField<RetType, fvsPatchField, surfaceMesh>& sf = tsf.ref();
261 
262  Field<RetType>& sfi = sf.primitiveFieldRef();
263 
264  const typename SFType::Internal& Sfi = Sf();
265 
266  for (label fi=0; fi<P.size(); fi++)
267  {
268  sfi[fi] = Sfi[fi] & (lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]]);
269  }
270 
271  // Interpolate across coupled patches using given lambdas
272 
273  typename GeometricField<RetType, fvsPatchField, surfaceMesh>::
274  Boundary& sfbf = sf.boundaryFieldRef();
275 
276  forAll(lambdas.boundaryField(), pi)
277  {
278  const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
279  const typename SFType::Patch& pSf = Sf.boundaryField()[pi];
280  fvsPatchField<RetType>& psf = sfbf[pi];
281 
282  if (vf.boundaryField()[pi].coupled())
283  {
284  psf =
285  pSf
286  & (
287  pLambda*vf.boundaryField()[pi].patchInternalField()
288  + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField()
289  );
290  }
291  else
292  {
293  psf = pSf & vf.boundaryField()[pi];
294  }
295  }
296 
297  tlambdas.clear();
298 
299 // tsf.ref().oriented() = Sf.oriented();
300 
301  return tsf;
302 }
303 
304 
305 template<class Type>
306 Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
307 Foam::surfaceInterpolationScheme<Type>::interpolate
308 (
309  const GeometricField<Type, fvPatchField, volMesh>& vf,
310  const tmp<surfaceScalarField>& tlambdas
311 )
312 {
313  return dotInterpolate(geometricOneField(), vf, tlambdas);
314 }
315 
316 
317 template<class Type>
318 Foam::tmp
319 <
320  Foam::GeometricField
321  <
322  typename Foam::innerProduct<Foam::vector, Type>::type,
323  Foam::fvsPatchField,
324  Foam::surfaceMesh
325  >
326 >
327 Foam::surfaceInterpolationScheme<Type>::dotInterpolate
328 (
329  const surfaceVectorField& Sf,
330  const GeometricField<Type, fvPatchField, volMesh>& vf
331 ) const
332 {
333  if (surfaceInterpolation::debug)
334  {
335  InfoInFunction
336  << "Interpolating "
337  << vf.type() << " "
338  << vf.name()
339  << " from cells to faces"
340  << endl;
341  }
342 
343  tmp
344  <
345  GeometricField
346  <
347  typename Foam::innerProduct<Foam::vector, Type>::type,
348  fvsPatchField,
349  surfaceMesh
350  >
351  > tsf = dotInterpolate(Sf, vf, weights(vf));
352 
353  tsf.ref().oriented() = Sf.oriented();
354 
355  if (corrected())
356  {
357  tsf.ref() += Sf & correction(vf);
358  }
359 
360  return tsf;
361 }
362 
363 
364 template<class Type>
365 Foam::tmp
366 <
367  Foam::GeometricField
368  <
369  typename Foam::innerProduct<Foam::vector, Type>::type,
370  Foam::fvsPatchField,
371  Foam::surfaceMesh
372  >
373 >
374 Foam::surfaceInterpolationScheme<Type>::dotInterpolate
375 (
376  const surfaceVectorField& Sf,
377  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf
378 ) const
379 {
380  tmp
381  <
382  GeometricField
383  <
384  typename Foam::innerProduct<Foam::vector, Type>::type,
385  fvsPatchField,
386  surfaceMesh
387  >
388  > tSfDotinterpVf = dotInterpolate(Sf, tvf());
389 
390  tvf.clear();
391  return tSfDotinterpVf;
392 }
393 
394 
395 template<class Type>
396 Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
397 Foam::surfaceInterpolationScheme<Type>::interpolate
398 (
399  const GeometricField<Type, fvPatchField, volMesh>& vf
400 ) const
401 {
402  if (surfaceInterpolation::debug)
403  {
404  InfoInFunction
405  << "Interpolating "
406  << vf.type() << " "
407  << vf.name()
408  << " from cells to faces"
409  << endl;
410  }
411 
412  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf
413  = interpolate(vf, weights(vf));
414 
415  if (corrected())
416  {
417  tsf.ref() += correction(vf);
418  }
419 
420  return tsf;
421 }
422 
423 
424 template<class Type>
425 Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh>>
426 Foam::surfaceInterpolationScheme<Type>::interpolate
427 (
428  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf
429 ) const
430 {
431  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tinterpVf
432  = interpolate(tvf());
433  tvf.clear();
434  return tinterpVf;
435 }
436 
437 
438 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
volFields.H
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:104
Foam::surfaceInterpolationScheme::New
static tmp< surfaceInterpolationScheme< Type > > New(const fvMesh &mesh, Istream &schemeData)
Return new tmp interpolation scheme.
Definition: surfaceInterpolationScheme.C:40
Foam::surfaceInterpolationScheme::dotInterpolate
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.
InfoInFunction
#define InfoInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:316
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
Foam::tmp::clear
void clear() const noexcept
Definition: tmpI.H:325
Foam::fvc::dotInterpolate
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:59
Foam::surfaceInterpolationScheme::interpolate
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.
Definition: surfaceInterpolationScheme.C:127
Foam::IOstream::eof
bool eof() const
Return true if end of input seen.
Definition: IOstream.H:222
Foam::fvsPatchScalarField
fvsPatchField< scalar > fvsPatchScalarField
Definition: fvsPatchFieldsFwd.H:45
Foam::surfaceMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: surfaceMesh.H:49
Foam::fvsPatchField
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:68
Foam::geometricOneField
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Definition: geometricOneField.H:54
Foam::FatalIOError
IOerror FatalIOError
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:337
surfaceFields.H
Foam::surfaceFields.
Foam::correction
tmp< fvMatrix< Type > > correction(const fvMatrix< Type > &)
Return the correction form of the given matrix.
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:290
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:258
FatalIOErrorInLookup
#define FatalIOErrorInLookup(ios, lookupTag, lookupName, lookupTable)
Report an error message using Foam::FatalIOError.
Definition: error.H:380
Foam::label
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:62
Foam::Field
Generic templated field type.
Definition: Field.H:63
Foam::Istream
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:61
Foam::FatalError
error FatalError
lambda
dimensionedScalar lambda("lambda", dimTime/sqr(dimLength), laminarTransport)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:84
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::GeometricField::primitiveFieldRef
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
Definition: GeometricField.C:735
geometricOneField.H
Foam::constant::mathematical::pi
constexpr scalar pi(M_PI)
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition: surfaceFieldsFwd.H:54
Foam::GeometricField::ref
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
Definition: GeometricField.C:718
Foam::GeometricField::boundaryFieldRef
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Definition: GeometricField.C:752
Foam::UList< label >
Foam::surfaceInterpolationScheme
Abstract base class for surface interpolation schemes.
Definition: surfaceInterpolationScheme.H:57
Foam::UList::size
void size(const label n) noexcept
Override size to be inconsistent with allocated storage.
Definition: UListI.H:360
FatalIOErrorInFunction
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:375
N
const Vector< label > N(dict.get< Vector< label >>("N"))
coupledFvPatchField.H
Foam::labelUList
UList< label > labelUList
A UList of labels.
Definition: UList.H:79
Foam::surfaceVectorField
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
Definition: surfaceFieldsFwd.H:57
Foam::GeometricField< scalar, fvsPatchField, surfaceMesh >
Foam::fac::interpolate
static tmp< GeometricField< Type, faePatchField, edgeMesh > > interpolate(const GeometricField< Type, faPatchField, areaMesh > &tvf, const edgeScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
surfaceInterpolationScheme.H
Foam::innerProduct::type
typeOfRank< typename pTraits< arg1 >::cmptType, direction(pTraits< arg1 >::rank)+direction(pTraits< arg2 >::rank) - 2 >::type type
Definition: products.H:149
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
y
scalar y
Definition: LISASMDCalcMethod1.H:14