surfaceInterpolate.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-2016 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 "surfaceInterpolate.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 template<class Type>
33 (
34  const surfaceScalarField& faceFlux,
35  Istream& streamData
36 )
37 {
39  (
40  faceFlux.mesh(),
41  faceFlux,
42  streamData
43  );
44 }
45 
46 
47 template<class Type>
49 (
50  const surfaceScalarField& faceFlux,
51  const word& name
52 )
53 {
55  (
56  faceFlux.mesh(),
57  faceFlux,
58  faceFlux.mesh().interpolationScheme(name)
59  );
60 }
61 
62 
63 template<class Type>
65 (
66  const fvMesh& mesh,
67  Istream& streamData
68 )
69 {
71  (
72  mesh,
73  streamData
74  );
75 }
76 
77 
78 template<class Type>
80 (
81  const fvMesh& mesh,
82  const word& name
83 )
84 {
86  (
87  mesh,
88  mesh.interpolationScheme(name)
89  );
90 }
91 
92 
93 template<class Type>
96 (
97  const GeometricField<Type, fvPatchField, volMesh>& vf,
98  const surfaceScalarField& faceFlux,
99  Istream& schemeData
100 )
101 {
102  if (surfaceInterpolation::debug)
103  {
105  << "interpolating GeometricField<Type, fvPatchField, volMesh> "
106  << vf.name() << endl;
107  }
108 
109  return scheme<Type>(faceFlux, schemeData)().interpolate(vf);
110 }
111 
112 
113 template<class Type>
116 (
117  const GeometricField<Type, fvPatchField, volMesh>& vf,
118  const surfaceScalarField& faceFlux,
119  const word& name
120 )
121 {
122  if (surfaceInterpolation::debug)
123  {
125  << "interpolating GeometricField<Type, fvPatchField, volMesh> "
126  << vf.name() << " using " << name << endl;
127  }
128 
129  return scheme<Type>(faceFlux, name)().interpolate(vf);
130 }
131 
132 template<class Type>
135 (
136  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf,
137  const surfaceScalarField& faceFlux,
138  const word& name
139 )
140 {
141  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf =
142  interpolate(tvf(), faceFlux, name);
143 
144  tvf.clear();
145 
146  return tsf;
147 }
148 
149 template<class Type>
152 (
153  const GeometricField<Type, fvPatchField, volMesh>& vf,
154  const tmp<surfaceScalarField>& tFaceFlux,
155  const word& name
156 )
157 {
158  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf =
159  interpolate(vf, tFaceFlux(), name);
160 
161  tFaceFlux.clear();
162 
163  return tsf;
164 }
165 
166 template<class Type>
169 (
170  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf,
171  const tmp<surfaceScalarField>& tFaceFlux,
172  const word& name
173 )
174 {
175  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf =
176  interpolate(tvf(), tFaceFlux(), name);
177 
178  tvf.clear();
179  tFaceFlux.clear();
180 
181  return tsf;
182 }
183 
184 
185 template<class Type>
188 (
189  const GeometricField<Type, fvPatchField, volMesh>& vf,
190  Istream& schemeData
191 )
192 {
193  if (surfaceInterpolation::debug)
194  {
196  << "interpolating GeometricField<Type, fvPatchField, volMesh> "
197  << vf.name() << endl;
198  }
199 
200  return scheme<Type>(vf.mesh(), schemeData)().interpolate(vf);
201 }
202 
203 template<class Type>
206 (
207  const GeometricField<Type, fvPatchField, volMesh>& vf,
208  const word& name
209 )
210 {
211  if (surfaceInterpolation::debug)
212  {
214  << "interpolating GeometricField<Type, fvPatchField, volMesh> "
215  << vf.name() << " using " << name
216  << endl;
217  }
218 
219  return scheme<Type>(vf.mesh(), name)().interpolate(vf);
220 }
221 
222 template<class Type>
225 (
226  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf,
227  const word& name
228 )
229 {
230  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf =
231  interpolate(tvf(), name);
232 
233  tvf.clear();
234 
235  return tsf;
236 }
237 
238 
239 template<class Type>
242 (
243  const GeometricField<Type, fvPatchField, volMesh>& vf
244 )
245 {
246  if (surfaceInterpolation::debug)
247  {
249  << "interpolating GeometricField<Type, fvPatchField, volMesh> "
250  << vf.name() << " using run-time selected scheme"
251  << endl;
252  }
253 
254  return interpolate(vf, "interpolate(" + vf.name() + ')');
255 }
256 
257 
258 template<class Type>
261 (
262  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf
263 )
264 {
265  tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tsf =
266  interpolate(tvf());
267  tvf.clear();
268  return tsf;
269 }
270 
271 
272 template<class Type>
275 (
276  const FieldField<fvPatchField, Type>& fvpff
277 )
278 {
279  FieldField<fvsPatchField, Type>* fvspffPtr
280  (
281  new FieldField<fvsPatchField, Type>(fvpff.size())
282  );
283 
284  forAll(*fvspffPtr, patchi)
285  {
286  fvspffPtr->set
287  (
288  patchi,
289  fvsPatchField<Type>::NewCalculatedType(fvpff[patchi].patch()).ptr()
290  );
291  (*fvspffPtr)[patchi] = fvpff[patchi];
292  }
293 
294  return tmp<FieldField<fvsPatchField, Type>>(fvspffPtr);
295 }
296 
297 
298 template<class Type>
301 (
302  const tmp<FieldField<fvPatchField, Type>>& tfvpff
303 )
304 {
305  tmp<FieldField<fvsPatchField, Type>> tfvspff = interpolate(tfvpff());
306  tfvpff.clear();
307  return tfvspff;
308 }
309 
310 
311 template<class Type>
312 Foam::tmp
313 <
315  <
316  typename Foam::innerProduct<Foam::vector, Type>::type,
317  Foam::fvsPatchField,
319  >
320 >
322 (
323  const surfaceVectorField& Sf,
324  const GeometricField<Type, fvPatchField, volMesh>& vf
325 )
326 {
327  if (surfaceInterpolation::debug)
328  {
330  << "interpolating GeometricField<Type, fvPatchField, volMesh> "
331  << vf.name() << " using run-time selected scheme"
332  << endl;
333  }
334 
335  return scheme<Type>
336  (
337  vf.mesh(),
338  "dotInterpolate(" + Sf.name() + ',' + vf.name() + ')'
339  )().dotInterpolate(Sf, vf);
340 }
341 
342 
343 template<class Type>
344 Foam::tmp
345 <
347  <
348  typename Foam::innerProduct<Foam::vector, Type>::type,
349  Foam::fvsPatchField,
351  >
352 >
354 (
355  const surfaceVectorField& Sf,
356  const tmp<GeometricField<Type, fvPatchField, volMesh>>& tvf
357 )
358 {
359  tmp
360  <
361  GeometricField
362  <
363  typename Foam::innerProduct<Foam::vector, Type>::type,
364  fvsPatchField,
365  surfaceMesh
366  >
367  > tsf = dotInterpolate(Sf, tvf());
368  tvf.clear();
369  return tsf;
370 }
371 
372 
373 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
Mesh data needed to do the Finite Volume discretisation.
Definition: surfaceMesh.H:47
static tmp< GeometricField< typename innerProduct< vector, Type >::type, fvsPatchField, surfaceMesh > > dotInterpolate(const surfaceVectorField &Sf, const GeometricField< Type, fvPatchField, volMesh > &tvf)
Interpolate field onto faces.
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
static tmp< surfaceInterpolationScheme< Type > > scheme(const surfaceScalarField &faceFlux, Istream &schemeData)
Return weighting factors for scheme given from Istream.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Generic GeometricField class.
tmp< surfaceScalarField > interpolate(const RhoType &rho)
Surface Interpolation.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
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.
label patchi
A class for managing temporary objects.
Definition: PtrList.H:54
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
#define InfoInFunction
Report an information message using Foam::Info.