meshToMesh0Templates.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 "meshToMesh0.H"
27 #include "volFields.H"
28 #include "interpolationCellPoint.H"
29 #include "SubField.H"
30 #include "mixedFvPatchField.H"
31 
32 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
33 
34 template<class Type, class CombineOp>
36 (
37  Field<Type>& toF,
38  const Field<Type>& fromVf,
39  const labelList& adr,
40  const CombineOp& cop
41 ) const
42 {
43  // Direct mapping of nearest-cell values
44 
45  forAll(toF, celli)
46  {
47  if (adr[celli] != -1)
48  {
49  cop(toF[celli], fromVf[adr[celli]]);
50  }
51  }
52 
53  //toF.map(fromVf, adr);
54 }
55 
56 
57 template<class Type, class CombineOp>
59 (
60  Field<Type>& toF,
62  const labelListList& adr,
63  const scalarListList& weights,
64  const CombineOp& cop
65 ) const
66 {
67  // Inverse volume weighted interpolation
68  forAll(toF, celli)
69  {
70  const labelList& overlapCells = adr[celli];
71  const scalarList& w = weights[celli];
72 
73  Type f = Zero;
74  forAll(overlapCells, i)
75  {
76  label fromCelli = overlapCells[i];
77  f += fromVf[fromCelli]*w[i];
78  cop(toF[celli], f);
79  }
80  }
81 }
82 
83 
84 template<class Type, class CombineOp>
86 (
87  Field<Type>& toF,
89  const labelList& adr,
90  const scalarListList& weights,
91  const CombineOp& cop
92 ) const
93 {
94  // Inverse distance weighted interpolation
95 
96  // get reference to cellCells
97  const labelListList& cc = fromMesh_.cellCells();
98 
99  forAll(toF, celli)
100  {
101  if (adr[celli] != -1)
102  {
103  const labelList& neighbours = cc[adr[celli]];
104  const scalarList& w = weights[celli];
105 
106  Type f = fromVf[adr[celli]]*w[0];
107 
108  for (label ni = 1; ni < w.size(); ni++)
109  {
110  f += fromVf[neighbours[ni - 1]]*w[ni];
111  }
112 
113  cop(toF[celli], f);
114  }
115  }
116 }
117 
118 
119 template<class Type, class CombineOp>
121 (
122  Field<Type>& toF,
124  const labelList& adr,
125  const vectorField& centres,
126  const CombineOp& cop
127 ) const
128 {
129  // Cell-Point interpolation
130  interpolationCellPoint<Type> interpolator(fromVf);
131 
132  forAll(toF, celli)
133  {
134  if (adr[celli] != -1)
135  {
136  cop
137  (
138  toF[celli],
139  interpolator.interpolate
140  (
141  centres[celli],
142  adr[celli]
143  )
144  );
145  }
146  }
147 }
148 
149 
150 template<class Type, class CombineOp>
152 (
153  Field<Type>& toF,
155  meshToMesh0::order ord,
156  const CombineOp& cop
157 ) const
158 {
159  if (fromVf.mesh() != fromMesh_)
160  {
162  << "the argument field does not correspond to the right mesh. "
163  << "Field size: " << fromVf.size()
164  << " mesh size: " << fromMesh_.nCells()
165  << exit(FatalError);
166  }
167 
168  if (toF.size() != toMesh_.nCells())
169  {
171  << "the argument field does not correspond to the right mesh. "
172  << "Field size: " << toF.size()
173  << " mesh size: " << toMesh_.nCells()
174  << exit(FatalError);
175  }
176 
177  switch(ord)
178  {
179  case MAP:
180  mapField(toF, fromVf, cellAddressing_, cop);
181  break;
182 
183  case INTERPOLATE:
184  {
185  interpolateField
186  (
187  toF,
188  fromVf,
189  cellAddressing_,
190  inverseDistanceWeights(),
191  cop
192  );
193  break;
194  }
195  case CELL_POINT_INTERPOLATE:
196  {
197  interpolateField
198  (
199  toF,
200  fromVf,
201  cellAddressing_,
202  toMesh_.cellCentres(),
203  cop
204  );
205 
206  break;
207  }
208  case CELL_VOLUME_WEIGHT:
209  {
210  const labelListList& cellToCell = cellToCellAddressing();
211  const scalarListList& invVolWeights = inverseVolumeWeights();
212 
213  interpolateField
214  (
215  toF,
216  fromVf,
217  cellToCell,
218  invVolWeights,
219  cop
220  );
221  break;
222  }
223  default:
225  << "unknown interpolation scheme " << ord
226  << exit(FatalError);
227  }
228 }
229 
230 
231 template<class Type, class CombineOp>
233 (
234  Field<Type>& toF,
236  meshToMesh0::order ord,
237  const CombineOp& cop
238 ) const
239 {
240  interpolateInternalField(toF, tfromVf(), ord, cop);
241  tfromVf.clear();
242 }
243 
244 
245 template<class Type, class CombineOp>
247 (
250  meshToMesh0::order ord,
251  const CombineOp& cop
252 ) const
253 {
254  interpolateInternalField(toVf, fromVf, ord, cop);
255 
257  Boundary& toVfBf = toVf.boundaryFieldRef();
258 
259  forAll(toMesh_.boundaryMesh(), patchi)
260  {
261  const fvPatch& toPatch = toMesh_.boundary()[patchi];
262 
263  if (cuttingPatches_.found(toPatch.name()))
264  {
265  switch(ord)
266  {
267  case MAP:
268  {
269  mapField
270  (
271  toVfBf[patchi],
272  fromVf,
273  boundaryAddressing_[patchi],
274  cop
275  );
276  break;
277  }
278 
279  case INTERPOLATE:
280  {
281  interpolateField
282  (
283  toVfBf[patchi],
284  fromVf,
285  boundaryAddressing_[patchi],
286  toPatch.Cf(),
287  cop
288  );
289  break;
290  }
291 
292  case CELL_POINT_INTERPOLATE:
293  {
294  interpolateField
295  (
296  toVfBf[patchi],
297  fromVf,
298  boundaryAddressing_[patchi],
299  toPatch.Cf(),
300  cop
301  );
302  break;
303  }
304  case CELL_VOLUME_WEIGHT:
305  {
306  break;
307  }
308 
309  default:
311  << "unknown interpolation scheme " << ord
312  << exit(FatalError);
313  }
314 
315  if (isA<mixedFvPatchField<Type>>(toVfBf[patchi]))
316  {
317  refCast<mixedFvPatchField<Type>>
318  (
319  toVfBf[patchi]
320  ).refValue() = toVfBf[patchi];
321  }
322  }
323  else if
324  (
325  patchMap_.found(toPatch.name())
326  && fromMeshPatches_.found(patchMap_.find(toPatch.name())())
327  )
328  {
329  mapField
330  (
331  toVfBf[patchi],
332  fromVf.boundaryField()
333  [
334  fromMeshPatches_.find(patchMap_.find(toPatch.name())())()
335  ],
336  boundaryAddressing_[patchi],
337  cop
338  );
339  }
340  }
341 }
342 
343 
344 template<class Type, class CombineOp>
346 (
349  meshToMesh0::order ord,
350  const CombineOp& cop
351 ) const
352 {
353  interpolate(toVf, tfromVf(), ord, cop);
354  tfromVf.clear();
355 }
356 
357 
358 template<class Type, class CombineOp>
361 (
363  meshToMesh0::order ord,
364  const CombineOp& cop
365 ) const
366 {
367  // Create and map the internal-field values
368  Field<Type> internalField(toMesh_.nCells());
369  interpolateInternalField(internalField, fromVf, ord, cop);
370 
371  // check whether both meshes have got the same number
372  // of boundary patches
373  if (fromMesh_.boundary().size() != toMesh_.boundary().size())
374  {
376  << "Incompatible meshes: different number of boundaries, "
377  "only internal field may be interpolated"
378  << exit(FatalError);
379  }
380 
381  // Create and map the patch field values
382  PtrList<fvPatchField<Type>> patchFields
383  (
384  boundaryAddressing_.size()
385  );
386 
387  forAll(boundaryAddressing_, patchi)
388  {
389  patchFields.set
390  (
391  patchi,
393  (
394  fromVf.boundaryField()[patchi],
395  toMesh_.boundary()[patchi],
398  (
399  boundaryAddressing_[patchi]
400  )
401  )
402  );
403  }
404 
405 
406  // Create the complete field from the pieces
408  (
410  (
411  IOobject
412  (
413  "interpolated(" + fromVf.name() + ')',
414  toMesh_.time().timeName(),
415  toMesh_,
416  IOobject::NO_READ,
417  IOobject::NO_WRITE
418  ),
419  toMesh_,
420  fromVf.dimensions(),
421  internalField,
422  patchFields
423  )
424  );
425 
426  return ttoF;
427 }
428 
429 
430 template<class Type, class CombineOp>
433 (
435  meshToMesh0::order ord,
436  const CombineOp& cop
437 ) const
438 {
440  interpolate(tfromVf(), ord, cop);
441  tfromVf.clear();
442 
443  return tint;
444 }
445 
446 
447 // ************************************************************************* //
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
#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:291
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
void mapField(Field< Type > &, const Field< Type > &, const labelList &adr, const CombineOp &cop) const
Map field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
This boundary condition provides a base class for &#39;mixed&#39; type boundary conditions, i.e. conditions that mix fixed value and patch-normal gradient conditions.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
order
Enumeration specifying required accuracy.
Definition: meshToMesh0.H:142
Generic GeometricField class.
void interpolateField(Field< Type > &, const GeometricField< Type, fvPatchField, volMesh > &, const labelList &adr, const scalarListList &weights, const CombineOp &cop) const
Interpolate field using inverse-distance weights.
void interpolateInternalField(Field< Type > &, const GeometricField< Type, fvPatchField, volMesh > &, order=INTERPOLATE, const CombineOp &cop=eqOp< Type >()) const
Interpolate internal volume field.
const word & name() const
Return name.
Definition: fvPatch.H:149
A topoSetSource to select the cells from another cellSet.
Definition: cellToCell.H:48
const dimensionSet & dimensions() const
Return dimensions.
Pre-declare SubField and related Field type.
Definition: Field.H:57
Type interpolate(const cellPointWeight &cpw) const
Interpolate field for the given cellPointWeight.
const vectorField & Cf() const
Return face centres.
Definition: fvPatch.C:99
static const zero Zero
Definition: zero.H:91
const Mesh & mesh() const
Return mesh.
labelList f(nPoints)
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
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
Patch-field interpolation class.
Definition: meshToMesh0.H:178
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void interpolate(GeometricField< Type, fvPatchField, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &, order=INTERPOLATE, const CombineOp &cop=eqOp< Type >()) const
Interpolate volume field.