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-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 \*---------------------------------------------------------------------------*/
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 = pTraits<Type>::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  (
163  "meshToMesh0::interpolateInternalField(Field<Type>&, "
164  "const GeometricField<Type, fvPatchField, volMesh>&, "
165  "meshToMesh0::order, const CombineOp&) const"
166  ) << "the argument field does not correspond to the right mesh. "
167  << "Field size: " << fromVf.size()
168  << " mesh size: " << fromMesh_.nCells()
169  << exit(FatalError);
170  }
171 
172  if (toF.size() != toMesh_.nCells())
173  {
175  (
176  "meshToMesh0::interpolateInternalField(Field<Type>&, "
177  "const GeometricField<Type, fvPatchField, volMesh>&, "
178  "meshToMesh0::order, const CombineOp&) const"
179  ) << "the argument field does not correspond to the right mesh. "
180  << "Field size: " << toF.size()
181  << " mesh size: " << toMesh_.nCells()
182  << exit(FatalError);
183  }
184 
185  switch(ord)
186  {
187  case MAP:
188  mapField(toF, fromVf, cellAddressing_, cop);
189  break;
190 
191  case INTERPOLATE:
192  {
193  interpolateField
194  (
195  toF,
196  fromVf,
197  cellAddressing_,
198  inverseDistanceWeights(),
199  cop
200  );
201  break;
202  }
203  case CELL_POINT_INTERPOLATE:
204  {
205  interpolateField
206  (
207  toF,
208  fromVf,
209  cellAddressing_,
210  toMesh_.cellCentres(),
211  cop
212  );
213 
214  break;
215  }
216  case CELL_VOLUME_WEIGHT:
217  {
218  const labelListList& cellToCell = cellToCellAddressing();
219  const scalarListList& invVolWeights = inverseVolumeWeights();
220 
221  interpolateField
222  (
223  toF,
224  fromVf,
225  cellToCell,
226  invVolWeights,
227  cop
228  );
229  break;
230  }
231  default:
233  (
234  "meshToMesh0::interpolateInternalField(Field<Type>&, "
235  "const GeometricField<Type, fvPatchField, volMesh>&, "
236  "meshToMesh0::order, const CombineOp&) const"
237  ) << "unknown interpolation scheme " << ord
238  << exit(FatalError);
239  }
240 }
241 
242 
243 template<class Type, class CombineOp>
245 (
246  Field<Type>& toF,
248  meshToMesh0::order ord,
249  const CombineOp& cop
250 ) const
251 {
252  interpolateInternalField(toF, tfromVf(), ord, cop);
253  tfromVf.clear();
254 }
255 
256 
257 template<class Type, class CombineOp>
259 (
262  meshToMesh0::order ord,
263  const CombineOp& cop
264 ) const
265 {
266  interpolateInternalField(toVf, fromVf, ord, cop);
267 
268  forAll(toMesh_.boundaryMesh(), patchi)
269  {
270  const fvPatch& toPatch = toMesh_.boundary()[patchi];
271 
272  if (cuttingPatches_.found(toPatch.name()))
273  {
274  switch(ord)
275  {
276  case MAP:
277  {
278  mapField
279  (
280  toVf.boundaryField()[patchi],
281  fromVf,
282  boundaryAddressing_[patchi],
283  cop
284  );
285  break;
286  }
287 
288  case INTERPOLATE:
289  {
290  interpolateField
291  (
292  toVf.boundaryField()[patchi],
293  fromVf,
294  boundaryAddressing_[patchi],
295  toPatch.Cf(),
296  cop
297  );
298  break;
299  }
300 
301  case CELL_POINT_INTERPOLATE:
302  {
303  interpolateField
304  (
305  toVf.boundaryField()[patchi],
306  fromVf,
307  boundaryAddressing_[patchi],
308  toPatch.Cf(),
309  cop
310  );
311  break;
312  }
313  case CELL_VOLUME_WEIGHT:
314  {
315  // Do nothing
316  break;
317  }
318 
319  default:
321  (
322  "meshToMesh0::interpolate("
323  "GeometricField<Type, fvPatchField, volMesh>&, "
324  "const GeometricField<Type, fvPatchField, volMesh>&, "
325  "meshToMesh0::order, const CombineOp&) const"
326  ) << "unknown interpolation scheme " << ord
327  << exit(FatalError);
328  }
329 
331  {
332  refCast<mixedFvPatchField<Type> >
333  (
334  toVf.boundaryField()[patchi]
335  ).refValue() = toVf.boundaryField()[patchi];
336  }
337  }
338  else if
339  (
340  patchMap_.found(toPatch.name())
341  && fromMeshPatches_.found(patchMap_.find(toPatch.name())())
342  )
343  {
344  /*
345  toVf.boundaryField()[patchi].map
346  (
347  fromVf.boundaryField()
348  [
349  fromMeshPatches_.find(patchMap_.find(toPatch.name())())()
350  ],
351  boundaryAddressing_[patchi]
352  );
353  */
354 
355  mapField
356  (
357  toVf.boundaryField()[patchi],
358  fromVf.boundaryField()
359  [
360  fromMeshPatches_.find(patchMap_.find(toPatch.name())())()
361  ],
362  boundaryAddressing_[patchi],
363  cop
364  );
365  }
366  }
367 }
368 
369 
370 template<class Type, class CombineOp>
372 (
375  meshToMesh0::order ord,
376  const CombineOp& cop
377 ) const
378 {
379  interpolate(toVf, tfromVf(), ord, cop);
380  tfromVf.clear();
381 }
382 
383 
384 template<class Type, class CombineOp>
387 (
389  meshToMesh0::order ord,
390  const CombineOp& cop
391 ) const
392 {
393  // Create and map the internal-field values
394  Field<Type> internalField(toMesh_.nCells());
395  interpolateInternalField(internalField, fromVf, ord, cop);
396 
397  // check whether both meshes have got the same number
398  // of boundary patches
399  if (fromMesh_.boundary().size() != toMesh_.boundary().size())
400  {
402  (
403  "meshToMesh0::interpolate"
404  "(const GeometricField<Type, fvPatchField, volMesh>&,"
405  "meshToMesh0::order, const CombineOp&) const"
406  ) << "Incompatible meshes: different number of boundaries, "
407  "only internal field may be interpolated"
408  << exit(FatalError);
409  }
410 
411  // Create and map the patch field values
412  PtrList<fvPatchField<Type> > patchFields
413  (
414  boundaryAddressing_.size()
415  );
416 
417  forAll(boundaryAddressing_, patchI)
418  {
419  patchFields.set
420  (
421  patchI,
423  (
424  fromVf.boundaryField()[patchI],
425  toMesh_.boundary()[patchI],
428  (
429  boundaryAddressing_[patchI]
430  )
431  )
432  );
433  }
434 
435 
436  // Create the complete field from the pieces
438  (
440  (
441  IOobject
442  (
443  "interpolated(" + fromVf.name() + ')',
444  toMesh_.time().timeName(),
445  toMesh_,
446  IOobject::NO_READ,
447  IOobject::NO_WRITE
448  ),
449  toMesh_,
450  fromVf.dimensions(),
452  patchFields
453  )
454  );
455 
456  return ttoF;
457 }
458 
459 
460 template<class Type, class CombineOp>
463 (
465  meshToMesh0::order ord,
466  const CombineOp& cop
467 ) const
468 {
470  interpolate(tfromVf(), ord, cop);
471  tfromVf.clear();
472 
473  return tint;
474 }
475 
476 
477 // ************************************************************************* //
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.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
labelList f(nPoints)
const vectorField & Cf() const
Return face centres.
Definition: fvPatch.C:99
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
void interpolateInternalField(Field< Type > &, const GeometricField< Type, fvPatchField, volMesh > &, order=INTERPOLATE, const CombineOp &cop=eqOp< Type >()) const
Interpolate internal volume field.
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
tmp< surfaceScalarField > interpolate(const RhoType &rho)
const Mesh & mesh() const
Return mesh.
const word & name() const
Return name.
Definition: fvPatch.H:149
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.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define forAll(list, i)
Definition: UList.H:421
label patchi
void mapField(Field< Type > &, const Field< Type > &, const labelList &adr, const CombineOp &cop) const
Map field.
const dimensionSet & dimensions() const
Return dimensions.
Type interpolate(const cellPointWeight &cpw) const
Interpolate field for the given cellPointWeight.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
const word & name() const
Return name.
Definition: IOobject.H:260
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Given cell centre values and point (vertex) values decompose into tetrahedra and linear interpolate w...
Generic GeometricField class.
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
void interpolate(GeometricField< Type, fvPatchField, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &, order=INTERPOLATE, const CombineOp &cop=eqOp< Type >()) const
Interpolate volume field.
order
Enumeration specifying required accuracy.
Definition: meshToMesh0.H:142
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Patch-field interpolation class.
Definition: meshToMesh0.H:178
bool isA(const Type &t)
Check if a dynamic_cast to typeid is possible.
Definition: typeInfo.H:134
A class for managing temporary objects.
Definition: PtrList.H:118
conserve internalField()+
A topoSetSource to select the cells from another cellSet.
Definition: cellToCell.H:48