meshToMeshTemplates.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2012-2018 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  << "Supplied field size is not equal to target mesh size" << nl
92  << " source mesh = " << srcToTgtCellAddr_.size() << nl
93  << " target mesh = " << tgtToSrcCellAddr_.size() << nl
94  << " supplied field = " << result.size()
95  << abort(FatalError);
96  }
97 
99 
100  if (singleMeshProc_ == -1)
101  {
102  const mapDistribute& map = srcMapPtr_();
103 
104  List<Type> work(srcField);
105  map.distribute(work);
106 
107  forAll(result, celli)
108  {
109  const labelList& srcAddress = tgtToSrcCellAddr_[celli];
110  const scalarList& srcWeight = tgtToSrcCellWght_[celli];
111 
112  if (srcAddress.size())
113  {
114 // result[celli] = Zero;
115  result[celli] *= (1.0 - sum(srcWeight));
116  forAll(srcAddress, i)
117  {
118  label srcI = srcAddress[i];
119  scalar w = srcWeight[i];
120  cbop(result[celli], celli, work[srcI], w);
121  }
122  }
123  }
124  }
125  else
126  {
127  forAll(result, celli)
128  {
129  const labelList& srcAddress = tgtToSrcCellAddr_[celli];
130  const scalarList& srcWeight = tgtToSrcCellWght_[celli];
131 
132  if (srcAddress.size())
133  {
134 // result[celli] = Zero;
135  result[celli] *= (1.0 - sum(srcWeight));
136  forAll(srcAddress, i)
137  {
138  label srcI = srcAddress[i];
139  scalar w = srcWeight[i];
140  cbop(result[celli], celli, srcField[srcI], w);
141  }
142  }
143  }
144  }
145 }
146 
147 
148 template<class Type, class CombineOp>
150 (
151  const Field<Type>& srcField,
152  const CombineOp& cop
153 ) const
154 {
155  tmp<Field<Type>> tresult
156  (
157  new Field<Type>
158  (
159  tgtToSrcCellAddr_.size(),
160  Zero
161  )
162  );
163 
164  mapSrcToTgt(srcField, cop, tresult.ref());
165 
166  return tresult;
167 }
168 
169 
170 template<class Type, class CombineOp>
172 (
173  const tmp<Field<Type>>& tsrcField,
174  const CombineOp& cop
175 ) const
176 {
177  return mapSrcToTgt(tsrcField(), cop);
178 }
179 
180 
181 template<class Type>
183 (
184  const Field<Type>& srcField
185 ) const
186 {
187  return mapSrcToTgt(srcField, plusEqOp<Type>());
188 }
189 
190 
191 template<class Type>
193 (
194  const tmp<Field<Type>>& tsrcField
195 ) const
196 {
197  return mapSrcToTgt(tsrcField());
198 }
199 
200 
201 template<class Type, class CombineOp>
203 (
204  const UList<Type>& tgtField,
205  const CombineOp& cop,
206  List<Type>& result
207 ) const
208 {
209  if (result.size() != srcToTgtCellAddr_.size())
210  {
212  << "Supplied field size is not equal to source mesh size" << nl
213  << " source mesh = " << srcToTgtCellAddr_.size() << nl
214  << " target mesh = " << tgtToSrcCellAddr_.size() << nl
215  << " supplied field = " << result.size()
216  << abort(FatalError);
217  }
218 
220 
221  if (singleMeshProc_ == -1)
222  {
223  const mapDistribute& map = tgtMapPtr_();
224 
225  List<Type> work(tgtField);
226  map.distribute(work);
227 
228  forAll(result, celli)
229  {
230  const labelList& tgtAddress = srcToTgtCellAddr_[celli];
231  const scalarList& tgtWeight = srcToTgtCellWght_[celli];
232 
233  if (tgtAddress.size())
234  {
235  result[celli] *= (1.0 - sum(tgtWeight));
236  forAll(tgtAddress, i)
237  {
238  label tgtI = tgtAddress[i];
239  scalar w = tgtWeight[i];
240  cbop(result[celli], celli, work[tgtI], w);
241  }
242  }
243  }
244  }
245  else
246  {
247  forAll(result, celli)
248  {
249  const labelList& tgtAddress = srcToTgtCellAddr_[celli];
250  const scalarList& tgtWeight = srcToTgtCellWght_[celli];
251 
252  if (tgtAddress.size())
253  {
254  result[celli] *= (1.0 - sum(tgtWeight));
255  forAll(tgtAddress, i)
256  {
257  label tgtI = tgtAddress[i];
258  scalar w = tgtWeight[i];
259  cbop(result[celli], celli, tgtField[tgtI], w);
260  }
261  }
262  }
263  }
264 }
265 
266 
267 template<class Type, class CombineOp>
269 (
270  const Field<Type>& tgtField,
271  const CombineOp& cop
272 ) const
273 {
274  tmp<Field<Type>> tresult
275  (
276  new Field<Type>
277  (
278  srcToTgtCellAddr_.size(),
279  Zero
280  )
281  );
282 
283  mapTgtToSrc(tgtField, cop, tresult.ref());
284 
285  return tresult;
286 }
287 
288 
289 template<class Type, class CombineOp>
291 (
292  const tmp<Field<Type>>& ttgtField,
293  const CombineOp& cop
294 ) const
295 {
296  return mapTgtToSrc(ttgtField(), cop);
297 }
298 
299 
300 template<class Type>
302 (
303  const Field<Type>& tgtField
304 ) const
305 {
306  return mapTgtToSrc(tgtField, plusEqOp<Type>());
307 }
308 
309 
310 template<class Type>
312 (
313  const tmp<Field<Type>>& ttgtField
314 ) const
315 {
316  return mapTgtToSrc(ttgtField(), plusEqOp<Type>());
317 }
318 
319 
320 template<class Type, class CombineOp>
321 void Foam::meshToMesh::mapAndOpSrcToTgt
322 (
323  const AMIInterpolation& AMI,
324  const Field<Type>& srcField,
325  Field<Type>& tgtField,
326  const CombineOp& cop
327 ) const
328 {
329  tgtField = pTraits<Type>::zero;
330 
332  (
333  srcField,
335  tgtField,
337  );
338 }
339 
340 
341 template<class Type, class CombineOp>
343 (
345  const CombineOp& cop,
347 ) const
348 {
349  mapSrcToTgt(field, cop, result.primitiveFieldRef());
350 
351  const PtrList<AMIInterpolation>& AMIList = patchAMIs();
352 
354  Boundary& resultBf = result.boundaryFieldRef();
355 
356  forAll(AMIList, i)
357  {
358  label srcPatchi = srcPatchID_[i];
359  label tgtPatchi = tgtPatchID_[i];
360 
361  const fvPatchField<Type>& srcField = field.boundaryField()[srcPatchi];
362  fvPatchField<Type>& tgtField = resultBf[tgtPatchi];
363 
364  // Clone and map (since rmap does not do general mapping)
365  tmp<fvPatchField<Type>> tnewTgt
366  (
368  (
369  srcField,
370  tgtField.patch(),
371  result(),
373  (
374  AMIList[i].singlePatchProc(),
375  (
376  AMIList[i].singlePatchProc() == -1
377  ? &AMIList[i].srcMap()
378  : nullptr
379  ),
380  AMIList[i].tgtAddress(),
381  AMIList[i].tgtWeights()
382  )
383  )
384  );
385 
386  // Transfer all mapped quantities (value and e.g. gradient) onto
387  // tgtField. Value will get overwritten below.
388  tgtField.rmap(tnewTgt(), identity(tgtField.size()));
389 
390  // Override value to account for CombineOp (note: is dummy template
391  // specialisation for plusEqOp)
392  mapAndOpSrcToTgt(AMIList[i], srcField, tgtField, cop);
393  }
394 
395  forAll(cuttingPatches_, i)
396  {
397  label patchi = cuttingPatches_[i];
398  fvPatchField<Type>& pf = resultBf[patchi];
399  pf == pf.patchInternalField();
400  }
401 }
402 
403 
404 template<class Type, class CombineOp>
407 (
409  const CombineOp& cop
410 ) const
411 {
413 
414  const fvMesh& tgtMesh = static_cast<const fvMesh&>(tgtRegion_);
415 
416  const fvBoundaryMesh& tgtBm = tgtMesh.boundary();
417  const typename fieldType::Boundary& srcBfld =
418  field.boundaryField();
419 
420  PtrList<fvPatchField<Type>> tgtPatchFields(tgtBm.size());
421 
422  // construct tgt boundary patch types as copy of 'field' boundary types
423  // note: this will provide place holders for fields with additional
424  // entries, but these values will need to be reset
425  forAll(tgtPatchID_, i)
426  {
427  label srcPatchi = srcPatchID_[i];
428  label tgtPatchi = tgtPatchID_[i];
429 
430  if (!tgtPatchFields.set(tgtPatchi))
431  {
432  tgtPatchFields.set
433  (
434  tgtPatchi,
436  (
437  srcBfld[srcPatchi],
438  tgtMesh.boundary()[tgtPatchi],
441  (
442  labelList(tgtMesh.boundary()[tgtPatchi].size(), -1)
443  )
444  )
445  );
446  }
447  }
448 
449  // Any unset tgtPatchFields become calculated
450  forAll(tgtPatchFields, tgtPatchi)
451  {
452  if (!tgtPatchFields.set(tgtPatchi))
453  {
454  // Note: use factory New method instead of direct generation of
455  // calculated so we keep constraints
456  tgtPatchFields.set
457  (
458  tgtPatchi,
460  (
462  tgtMesh.boundary()[tgtPatchi],
464  )
465  );
466  }
467  }
468 
469  tmp<fieldType> tresult
470  (
471  new fieldType
472  (
473  IOobject
474  (
475  type() + ":interpolate(" + field.name() + ")",
476  tgtMesh.time().timeName(),
477  tgtMesh,
480  ),
481  tgtMesh,
482  field.dimensions(),
483  Field<Type>(tgtMesh.nCells(), Zero),
484  tgtPatchFields
485  )
486  );
487 
488  mapSrcToTgt(field, cop, tresult.ref());
489 
490  return tresult;
491 }
492 
493 
494 template<class Type, class CombineOp>
497 (
499  const CombineOp& cop
500 ) const
501 {
502  return mapSrcToTgt(tfield(), cop);
503 }
504 
505 
506 template<class Type>
509 (
511 ) const
512 {
513  return mapSrcToTgt(field, plusEqOp<Type>());
514 }
515 
516 
517 template<class Type>
520 (
522 ) const
523 {
524  return mapSrcToTgt(tfield(), plusEqOp<Type>());
525 }
526 
527 
528 template<class Type, class CombineOp>
529 void Foam::meshToMesh::mapAndOpTgtToSrc
530 (
531  const AMIInterpolation& AMI,
532  Field<Type>& srcField,
533  const Field<Type>& tgtField,
534  const CombineOp& cop
535 ) const
536 {
537  srcField = pTraits<Type>::zero;
538 
540  (
541  tgtField,
543  srcField,
545  );
546 }
547 
548 
549 template<class Type, class CombineOp>
551 (
553  const CombineOp& cop,
555 ) const
556 {
557  mapTgtToSrc(field, cop, result.primitiveFieldRef());
558 
559  const PtrList<AMIInterpolation>& AMIList = patchAMIs();
560 
561  forAll(AMIList, i)
562  {
563  label srcPatchi = srcPatchID_[i];
564  label tgtPatchi = tgtPatchID_[i];
565 
566  fvPatchField<Type>& srcField = result.boundaryFieldRef()[srcPatchi];
567  const fvPatchField<Type>& tgtField = field.boundaryField()[tgtPatchi];
568 
569 
570  // Clone and map (since rmap does not do general mapping)
571  tmp<fvPatchField<Type>> tnewSrc
572  (
574  (
575  tgtField,
576  srcField.patch(),
577  result(),
579  (
580  AMIList[i].singlePatchProc(),
581  (
582  AMIList[i].singlePatchProc() == -1
583  ? &AMIList[i].tgtMap()
584  : nullptr
585  ),
586  AMIList[i].srcAddress(),
587  AMIList[i].srcWeights()
588  )
589  )
590  );
591 
592  // Transfer all mapped quantities (value and e.g. gradient) onto
593  // srcField. Value will get overwritten below
594  srcField.rmap(tnewSrc(), identity(srcField.size()));
595 
596 
597  // Override value to account for CombineOp (could be dummy for
598  // plusEqOp)
599  mapAndOpTgtToSrc(AMIList[i], srcField, tgtField, cop);
600  }
601 
602  forAll(cuttingPatches_, i)
603  {
604  label patchi = cuttingPatches_[i];
606  pf == pf.patchInternalField();
607  }
608 }
609 
610 
611 template<class Type, class CombineOp>
614 (
616  const CombineOp& cop
617 ) const
618 {
620 
621  const fvMesh& srcMesh = static_cast<const fvMesh&>(srcRegion_);
622 
623  const fvBoundaryMesh& srcBm = srcMesh.boundary();
624  const typename fieldType::Boundary& tgtBfld =
625  field.boundaryField();
626 
627  PtrList<fvPatchField<Type>> srcPatchFields(srcBm.size());
628 
629  // construct src boundary patch types as copy of 'field' boundary types
630  // note: this will provide place holders for fields with additional
631  // entries, but these values will need to be reset
632  forAll(srcPatchID_, i)
633  {
634  label srcPatchi = srcPatchID_[i];
635  label tgtPatchi = tgtPatchID_[i];
636 
637  if (!srcPatchFields.set(tgtPatchi))
638  {
639  srcPatchFields.set
640  (
641  srcPatchi,
643  (
644  tgtBfld[srcPatchi],
645  srcMesh.boundary()[tgtPatchi],
648  (
649  labelList(srcMesh.boundary()[srcPatchi].size(), -1)
650  )
651  )
652  );
653  }
654  }
655 
656  // Any unset srcPatchFields become calculated
657  forAll(srcPatchFields, srcPatchi)
658  {
659  if (!srcPatchFields.set(srcPatchi))
660  {
661  // Note: use factory New method instead of direct generation of
662  // calculated so we keep constraints
663  srcPatchFields.set
664  (
665  srcPatchi,
667  (
669  srcMesh.boundary()[srcPatchi],
671  )
672  );
673  }
674  }
675 
676  tmp<fieldType> tresult
677  (
678  new fieldType
679  (
680  IOobject
681  (
682  type() + ":interpolate(" + field.name() + ")",
683  srcMesh.time().timeName(),
684  srcMesh,
687  ),
688  srcMesh,
689  field.dimensions(),
690  Field<Type>(srcMesh.nCells(), Zero),
691  srcPatchFields
692  )
693  );
694 
695  mapTgtToSrc(field, cop, tresult.ref());
696 
697  return tresult;
698 }
699 
700 
701 template<class Type, class CombineOp>
704 (
706  const CombineOp& cop
707 ) const
708 {
709  return mapTgtToSrc(tfield(), cop);
710 }
711 
712 
713 template<class Type>
716 (
718 ) const
719 {
720  return mapTgtToSrc(field, plusEqOp<Type>());
721 }
722 
723 
724 template<class Type>
727 (
729 ) const
730 {
731  return mapTgtToSrc(tfield(), plusEqOp<Type>());
732 }
733 
734 
735 // ************************************************************************* //
void operator()(List< Type > &x, const List< Type > y) const
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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:303
FieldMapper with weighted mapping from (optionally remote) quantities.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
label nCells() const
Traits class for primitives.
Definition: pTraits.H:50
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
Generic GeometricField class.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
scalar y
const dimensionSet & dimensions() const
Return dimensions.
static const DimensionedField< Type, GeoMesh > & null()
Return a null DimensionedField.
Pre-declare SubField and related Field type.
Definition: Field.H:56
This boundary condition is not designed to be evaluated; it is assumed that the value is assigned via...
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:199
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:97
void interpolateToTarget(const UList< Type > &fld, const CombineOp &cop, List< Type > &result, const UList< Type > &defaultValues=UList< Type >::null()) const
Interpolate from source to target with supplied op.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:337
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
static const char nl
Definition: Ostream.H:260
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:174
void interpolateToSource(const UList< Type > &fld, const CombineOp &cop, List< Type > &result, const UList< Type > &defaultValues=UList< Type >::null()) const
Interpolate from target to source with supplied op.
Internal & ref()
Return a reference to the dimensioned internal field.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
void setSize(const label)
Reset size of List.
Definition: List.C:281
label patchi
Class containing processor-to-processor mapping information.
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:70
Foam::fvBoundaryMesh.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
A class for managing temporary objects.
Definition: PtrList.H:53
void mapSrcToTgt(const UList< Type > &srcFld, const CombineOp &cop, List< Type > &result) const
Map field from src to tgt mesh with defined operation.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.
void mapTgtToSrc(const UList< Type > &tgtFld, const CombineOp &cop, List< Type > &result) const
Map field from tgt to src mesh with defined operation.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540