meshToMeshTemplates.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) 2012-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 "fvMesh.H"
27 #include "volFields.H"
29 #include "calculatedFvPatchField.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  //- Helper class for list
37  template<class Type>
38  class ListPlusEqOp
39  {
40  public:
41  void operator()(List<Type>& x, const List<Type> y) const
42  {
43  if (y.size())
44  {
45  if (x.size())
46  {
47  label sz = x.size();
48  x.setSize(sz + y.size());
49  forAll(y, i)
50  {
51  x[sz++] = y[i];
52  }
53  }
54  else
55  {
56  x = y;
57  }
58  }
59  }
60  };
61 }
62 
63 
64 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
65 
66 template<class Type>
67 void Foam::meshToMesh::add
68 (
69  UList<Type>& fld,
70  const label offset
71 ) const
72 {
73  forAll(fld, i)
74  {
75  fld[i] += offset;
76  }
77 }
78 
79 
80 template<class Type, class CombineOp>
82 (
83  const UList<Type>& srcField,
84  const CombineOp& cop,
85  List<Type>& result
86 ) const
87 {
88  if (result.size() != tgtToSrcCellAddr_.size())
89  {
91  (
92  "void Foam::meshToMesh::mapSrcToTgt"
93  "("
94  "const UList<Type>&, "
95  "const CombineOp&, "
96  "List<Type>&"
97  ") const"
98  ) << "Supplied field size is not equal to target mesh size" << nl
99  << " source mesh = " << srcToTgtCellAddr_.size() << nl
100  << " target mesh = " << tgtToSrcCellAddr_.size() << nl
101  << " supplied field = " << result.size()
102  << abort(FatalError);
103  }
104 
106 
107  if (singleMeshProc_ == -1)
108  {
109  const mapDistribute& map = srcMapPtr_();
110 
111  List<Type> work(srcField);
112  map.distribute(work);
113 
114  forAll(result, cellI)
115  {
116  const labelList& srcAddress = tgtToSrcCellAddr_[cellI];
117  const scalarList& srcWeight = tgtToSrcCellWght_[cellI];
118 
119  if (srcAddress.size())
120  {
121 // result[cellI] = pTraits<Type>::zero;
122  result[cellI] *= (1.0 - sum(srcWeight));
123  forAll(srcAddress, i)
124  {
125  label srcI = srcAddress[i];
126  scalar w = srcWeight[i];
127  cbop(result[cellI], cellI, work[srcI], w);
128  }
129  }
130  }
131  }
132  else
133  {
134  forAll(result, cellI)
135  {
136  const labelList& srcAddress = tgtToSrcCellAddr_[cellI];
137  const scalarList& srcWeight = tgtToSrcCellWght_[cellI];
138 
139  if (srcAddress.size())
140  {
141 // result[cellI] = pTraits<Type>::zero;
142  result[cellI] *= (1.0 - sum(srcWeight));
143  forAll(srcAddress, i)
144  {
145  label srcI = srcAddress[i];
146  scalar w = srcWeight[i];
147  cbop(result[cellI], cellI, srcField[srcI], w);
148  }
149  }
150  }
151  }
152 }
153 
154 
155 template<class Type, class CombineOp>
157 (
158  const Field<Type>& srcField,
159  const CombineOp& cop
160 ) const
161 {
162  tmp<Field<Type> > tresult
163  (
164  new Field<Type>
165  (
166  tgtToSrcCellAddr_.size(),
168  )
169  );
170 
171  mapSrcToTgt(srcField, cop, tresult());
172 
173  return tresult;
174 }
175 
176 
177 template<class Type, class CombineOp>
179 (
180  const tmp<Field<Type> >& tsrcField,
181  const CombineOp& cop
182 ) const
183 {
184  return mapSrcToTgt(tsrcField(), cop);
185 }
186 
187 
188 template<class Type>
190 (
191  const Field<Type>& srcField
192 ) const
193 {
194  return mapSrcToTgt(srcField, plusEqOp<Type>());
195 }
196 
197 
198 template<class Type>
200 (
201  const tmp<Field<Type> >& tsrcField
202 ) const
203 {
204  return mapSrcToTgt(tsrcField());
205 }
206 
207 
208 template<class Type, class CombineOp>
210 (
211  const UList<Type>& tgtField,
212  const CombineOp& cop,
213  List<Type>& result
214 ) const
215 {
216  if (result.size() != srcToTgtCellAddr_.size())
217  {
219  (
220  "void Foam::meshToMesh::mapTgtToSrc"
221  "("
222  "const UList<Type>&, "
223  "const CombineOp&, "
224  "List<Type>&"
225  ") const"
226  ) << "Supplied field size is not equal to source mesh size" << nl
227  << " source mesh = " << srcToTgtCellAddr_.size() << nl
228  << " target mesh = " << tgtToSrcCellAddr_.size() << nl
229  << " supplied field = " << result.size()
230  << abort(FatalError);
231  }
232 
234 
235  if (singleMeshProc_ == -1)
236  {
237  const mapDistribute& map = tgtMapPtr_();
238 
239  List<Type> work(tgtField);
240  map.distribute(work);
241 
242  forAll(result, cellI)
243  {
244  const labelList& tgtAddress = srcToTgtCellAddr_[cellI];
245  const scalarList& tgtWeight = srcToTgtCellWght_[cellI];
246 
247  if (tgtAddress.size())
248  {
249  result[cellI] *= (1.0 - sum(tgtWeight));
250  forAll(tgtAddress, i)
251  {
252  label tgtI = tgtAddress[i];
253  scalar w = tgtWeight[i];
254  cbop(result[cellI], cellI, work[tgtI], w);
255  }
256  }
257  }
258  }
259  else
260  {
261  forAll(result, cellI)
262  {
263  const labelList& tgtAddress = srcToTgtCellAddr_[cellI];
264  const scalarList& tgtWeight = srcToTgtCellWght_[cellI];
265 
266  if (tgtAddress.size())
267  {
268  result[cellI] *= (1.0 - sum(tgtWeight));
269  forAll(tgtAddress, i)
270  {
271  label tgtI = tgtAddress[i];
272  scalar w = tgtWeight[i];
273  cbop(result[cellI], cellI, tgtField[tgtI], w);
274  }
275  }
276  }
277  }
278 }
279 
280 
281 template<class Type, class CombineOp>
283 (
284  const Field<Type>& tgtField,
285  const CombineOp& cop
286 ) const
287 {
288  tmp<Field<Type> > tresult
289  (
290  new Field<Type>
291  (
292  srcToTgtCellAddr_.size(),
294  )
295  );
296 
297  mapTgtToSrc(tgtField, cop, tresult());
298 
299  return tresult;
300 }
301 
302 
303 template<class Type, class CombineOp>
305 (
306  const tmp<Field<Type> >& ttgtField,
307  const CombineOp& cop
308 ) const
309 {
310  return mapTgtToSrc(ttgtField(), cop);
311 }
312 
313 
314 template<class Type>
316 (
317  const Field<Type>& tgtField
318 ) const
319 {
320  return mapTgtToSrc(tgtField, plusEqOp<Type>());
321 }
322 
323 
324 template<class Type>
326 (
327  const tmp<Field<Type> >& ttgtField
328 ) const
329 {
330  return mapTgtToSrc(ttgtField(), plusEqOp<Type>());
331 }
332 
333 
334 template<class Type, class CombineOp>
336 (
338  const CombineOp& cop,
340 ) const
341 {
342  mapSrcToTgt(field, cop, result.internalField());
343 
344  const PtrList<AMIPatchToPatchInterpolation>& AMIList = patchAMIs();
345 
346  forAll(AMIList, i)
347  {
348  label srcPatchI = srcPatchID_[i];
349  label tgtPatchI = tgtPatchID_[i];
350 
351  const fvPatchField<Type>& srcField = field.boundaryField()[srcPatchI];
352  fvPatchField<Type>& tgtField = result.boundaryField()[tgtPatchI];
353 
354  // 2.3 does not do distributed mapping yet so only do if
355  // running on single processor
356  if (AMIList[i].singlePatchProc() != -1)
357  {
358  // Clone and map (since rmap does not do general mapping)
359  tmp<fvPatchField<Type> > tnewTgt
360  (
362  (
363  srcField,
364  tgtField.patch(),
365  result.dimensionedInternalField(),
367  (
368  AMIList[i].tgtAddress(),
369  AMIList[i].tgtWeights()
370  )
371  )
372  );
373 
374  // Transfer all mapped quantities (value and e.g. gradient) onto
375  // tgtField. Value will get overwritten below.
376  tgtField.rmap(tnewTgt(), identity(tgtField.size()));
377  }
378 
379  tgtField == pTraits<Type>::zero;
380 
381  AMIList[i].interpolateToTarget
382  (
383  srcField,
385  tgtField,
387  );
388  }
389 
390  forAll(cuttingPatches_, i)
391  {
392  label patchI = cuttingPatches_[i];
393  fvPatchField<Type>& pf = result.boundaryField()[patchI];
394  pf == pf.patchInternalField();
395  }
396 }
397 
398 
399 template<class Type, class CombineOp>
402 (
404  const CombineOp& cop
405 ) const
406 {
408 
409  const fvMesh& tgtMesh = static_cast<const fvMesh&>(tgtRegion_);
410 
411  const fvBoundaryMesh& tgtBm = tgtMesh.boundary();
412  const typename fieldType::GeometricBoundaryField& srcBfld =
413  field.boundaryField();
414 
415  PtrList<fvPatchField<Type> > tgtPatchFields(tgtBm.size());
416 
417  // constuct tgt boundary patch types as copy of 'field' boundary types
418  // note: this will provide place holders for fields with additional
419  // entries, but these values will need to be reset
420  forAll(tgtPatchID_, i)
421  {
422  label srcPatchI = srcPatchID_[i];
423  label tgtPatchI = tgtPatchID_[i];
424 
425  if (!tgtPatchFields.set(tgtPatchI))
426  {
427  tgtPatchFields.set
428  (
429  tgtPatchI,
431  (
432  srcBfld[srcPatchI],
433  tgtMesh.boundary()[tgtPatchI],
436  (
437  labelList(tgtMesh.boundary()[tgtPatchI].size(), -1)
438  )
439  )
440  );
441  }
442  }
443 
444  // Any unset tgtPatchFields become calculated
445  forAll(tgtPatchFields, tgtPatchI)
446  {
447  if (!tgtPatchFields.set(tgtPatchI))
448  {
449  // Note: use factory New method instead of direct generation of
450  // calculated so we keep constraints
451  tgtPatchFields.set
452  (
453  tgtPatchI,
455  (
457  tgtMesh.boundary()[tgtPatchI],
459  )
460  );
461  }
462  }
463 
464  tmp<fieldType> tresult
465  (
466  new fieldType
467  (
468  IOobject
469  (
470  type() + ":interpolate(" + field.name() + ")",
471  tgtMesh.time().timeName(),
472  tgtMesh,
475  ),
476  tgtMesh,
477  field.dimensions(),
479  tgtPatchFields
480  )
481  );
482 
483  mapSrcToTgt(field, cop, tresult());
484 
485  return tresult;
486 }
487 
488 
489 template<class Type, class CombineOp>
492 (
494  const CombineOp& cop
495 ) const
496 {
497  return mapSrcToTgt(tfield(), cop);
498 }
499 
500 
501 template<class Type>
504 (
506 ) const
507 {
508  return mapSrcToTgt(field, plusEqOp<Type>());
509 }
510 
511 
512 template<class Type>
515 (
517 ) const
518 {
519  return mapSrcToTgt(tfield(), plusEqOp<Type>());
520 }
521 
522 
523 template<class Type, class CombineOp>
525 (
527  const CombineOp& cop,
529 ) const
530 {
531  mapTgtToSrc(field, cop, result.internalField());
532 
533  const PtrList<AMIPatchToPatchInterpolation>& AMIList = patchAMIs();
534 
535  forAll(AMIList, i)
536  {
537  label srcPatchI = srcPatchID_[i];
538  label tgtPatchI = tgtPatchID_[i];
539 
540  fvPatchField<Type>& srcField = result.boundaryField()[srcPatchI];
541  const fvPatchField<Type>& tgtField = field.boundaryField()[tgtPatchI];
542 
543  // 2.3 does not do distributed mapping yet so only do if
544  // running on single processor
545  if (AMIList[i].singlePatchProc() != -1)
546  {
547  // Clone and map (since rmap does not do general mapping)
548  tmp<fvPatchField<Type> > tnewSrc
549  (
551  (
552  tgtField,
553  srcField.patch(),
554  result.dimensionedInternalField(),
556  (
557  AMIList[i].srcAddress(),
558  AMIList[i].srcWeights()
559  )
560  )
561  );
562 
563  // Transfer all mapped quantities (value and e.g. gradient) onto
564  // srcField. Value will get overwritten below
565  srcField.rmap(tnewSrc(), identity(srcField.size()));
566  }
567 
568  srcField == pTraits<Type>::zero;
569 
570  AMIList[i].interpolateToSource
571  (
572  tgtField,
574  srcField,
576  );
577  }
578 
579  forAll(cuttingPatches_, i)
580  {
581  label patchI = cuttingPatches_[i];
582  fvPatchField<Type>& pf = result.boundaryField()[patchI];
583  pf == pf.patchInternalField();
584  }
585 }
586 
587 
588 template<class Type, class CombineOp>
591 (
593  const CombineOp& cop
594 ) const
595 {
597 
598  const fvMesh& srcMesh = static_cast<const fvMesh&>(srcRegion_);
599 
600  const fvBoundaryMesh& srcBm = srcMesh.boundary();
601  const typename fieldType::GeometricBoundaryField& tgtBfld =
602  field.boundaryField();
603 
604  PtrList<fvPatchField<Type> > srcPatchFields(srcBm.size());
605 
606  // constuct src boundary patch types as copy of 'field' boundary types
607  // note: this will provide place holders for fields with additional
608  // entries, but these values will need to be reset
609  forAll(srcPatchID_, i)
610  {
611  label srcPatchI = srcPatchID_[i];
612  label tgtPatchI = tgtPatchID_[i];
613 
614  if (!srcPatchFields.set(tgtPatchI))
615  {
616  srcPatchFields.set
617  (
618  srcPatchI,
620  (
621  tgtBfld[srcPatchI],
622  srcMesh.boundary()[tgtPatchI],
625  (
626  labelList(srcMesh.boundary()[srcPatchI].size(), -1)
627  )
628  )
629  );
630  }
631  }
632 
633  // Any unset srcPatchFields become calculated
634  forAll(srcPatchFields, srcPatchI)
635  {
636  if (!srcPatchFields.set(srcPatchI))
637  {
638  // Note: use factory New method instead of direct generation of
639  // calculated so we keep constraints
640  srcPatchFields.set
641  (
642  srcPatchI,
644  (
646  srcMesh.boundary()[srcPatchI],
648  )
649  );
650  }
651  }
652 
653  tmp<fieldType> tresult
654  (
655  new fieldType
656  (
657  IOobject
658  (
659  type() + ":interpolate(" + field.name() + ")",
660  srcMesh.time().timeName(),
661  srcMesh,
664  ),
665  srcMesh,
666  field.dimensions(),
668  srcPatchFields
669  )
670  );
671 
672  mapTgtToSrc(field, cop, tresult());
673 
674  return tresult;
675 }
676 
677 
678 template<class Type, class CombineOp>
681 (
683  const CombineOp& cop
684 ) const
685 {
686  return mapTgtToSrc(tfield(), cop);
687 }
688 
689 
690 template<class Type>
693 (
695 ) const
696 {
697  return mapTgtToSrc(field, plusEqOp<Type>());
698 }
699 
700 
701 template<class Type>
704 (
706 ) const
707 {
708  return mapTgtToSrc(tfield(), plusEqOp<Type>());
709 }
710 
711 
712 // ************************************************************************* //
This boundary condition is not designed to be evaluated; it is assmued that the value is assigned via...
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
static void distribute(const Pstream::commsTypes commsType, const List< labelPair > &schedule, const label constructSize, const labelListList &subMap, const labelListList &constructMap, List< T > &, const int tag=UPstream::msgType())
Distribute data. Note:schedule only used for Pstream::scheduled.
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.
void mapTgtToSrc(const UList< Type > &tgtFld, const CombineOp &cop, List< Type > &result) const
Map field from tgt to src mesh with defined operation.
Class containing processor-to-processor mapping information.
Foam::fvBoundaryMesh.
void mapSrcToTgt(const UList< Type > &srcFld, const CombineOp &cop, List< Type > &result) const
Map field from src to tgt mesh with defined operation.
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
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
InternalField & internalField()
Return internal field.
label nCells() const
DimensionedInternalField & dimensionedInternalField()
Return dimensioned internal field.
Namespace for OpenFOAM.
FieldMapper with weighted mapping.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:291
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:215
#define forAll(list, i)
Definition: UList.H:421
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const dimensionSet & dimensions() const
Return dimensions.
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
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const word & name() const
Return name.
Definition: IOobject.H:260
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Generic GeometricField class.
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
scalar y
List< label > labelList
A List of labels.
Definition: labelList.H:56
void operator()(List< Type > &x, const List< Type > y) const
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:300
A class for managing temporary objects.
Definition: PtrList.H:118
static const DimensionedField< Type, GeoMesh > & null()
Return a null DimensionedField.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:552