fvMeshSubsetInterpolate.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) 2011-2021 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 "fvMeshSubset.H"
27 #include "surfaceFields.H"
28 #include "internalFvsPatchField.H"
30 #include "internalFvPatchFields.H"
32 #include "directPointPatchFieldMapper.H"
33 #include "flipOp.H"
34 
35 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
36 
37 template<class Type>
40 (
42  const fvMesh& sMesh,
43  const labelList& patchMap,
44  const labelList& cellMap,
45  const labelList& faceMap
46 )
47 {
48  // 1. Create the complete field with dummy patch fields
49  PtrList<fvPatchField<Type>> patchFields(patchMap.size());
50 
51  forAll(patchFields, patchi)
52  {
53  // Set the first one by hand as it corresponds to the
54  // exposed internal faces. Additional interpolation can be put here
55  // as necessary.
56  if (patchMap[patchi] == -1)
57  {
58  patchFields.set
59  (
60  patchi,
62  (
63  sMesh.boundary()[patchi],
65  )
66  );
67  }
68  else
69  {
70  patchFields.set
71  (
72  patchi,
74  (
76  sMesh.boundary()[patchi],
78  )
79  );
80  }
81  }
82 
84  (
86  (
87  IOobject
88  (
89  "subset"+vf.name(),
90  sMesh.time().timeName(),
91  sMesh,
92  IOobject::NO_READ,
93  IOobject::NO_WRITE,
94  false
95  ),
96  sMesh,
97  vf.dimensions(),
98  Field<Type>(vf.primitiveField(), cellMap),
99  patchFields
100  )
101  );
103 
104 
105  // 2. Change the fvPatchFields to the correct type using a mapper
106  // constructor (with reference to the now correct internal field)
107 
109  Boundary& bf = resF.boundaryFieldRef();
110 
111  forAll(bf, patchi)
112  {
113  if (patchMap[patchi] != -1)
114  {
115  // Construct addressing
116  const fvPatch& subPatch = sMesh.boundary()[patchi];
117  const fvPatch& basePatch = vf.mesh().boundary()[patchMap[patchi]];
118  const label baseStart = basePatch.start();
119  const label baseSize = basePatch.size();
120 
121  labelList directAddressing(subPatch.size());
122 
123  forAll(directAddressing, i)
124  {
125  const label baseFacei = faceMap[subPatch.start() + i];
126 
127  if (baseFacei >= baseStart && baseFacei < baseStart+baseSize)
128  {
129  directAddressing[i] = baseFacei-baseStart;
130  }
131  else
132  {
133  // Mapped from internal face. Do what? Leave up to
134  // fvPatchField
135  directAddressing[i] = -1;
136  }
137  }
138 
139  bf.set
140  (
141  patchi,
143  (
144  vf.boundaryField()[patchMap[patchi]],
145  subPatch,
146  resF(),
147  directFvPatchFieldMapper(directAddressing)
148  )
149  );
150  }
151  }
152 
153  return tresF;
154 }
155 
156 
157 template<class Type>
160 (
162 ) const
163 {
164  return interpolate
165  (
166  vf,
167  subMesh(),
168  patchMap(),
169  cellMap(),
170  faceMap()
171  );
172 }
173 
174 
175 template<class Type>
178 (
180  const fvMesh& sMesh,
181  const labelList& patchMap,
182  const labelList& cellMap,
183  const labelList& faceMap
184 )
185 {
186  const bool negateIfFlipped = isFlux(sf);
187 
188  // 1. Create the complete field with dummy patch fields
189  PtrList<fvsPatchField<Type>> patchFields(patchMap.size());
190 
191  forAll(patchFields, patchi)
192  {
193  // Set the first one by hand as it corresponds to the
194  // exposed internal faces. Additional interpolation can be put here
195  // as necessary.
196  if (patchMap[patchi] == -1)
197  {
198  patchFields.set
199  (
200  patchi,
202  (
203  sMesh.boundary()[patchi],
205  )
206  );
207  }
208  else
209  {
210  patchFields.set
211  (
212  patchi,
214  (
216  sMesh.boundary()[patchi],
218  )
219  );
220  }
221  }
222 
223  // Create the complete field from the pieces
225  (
227  (
228  IOobject
229  (
230  "subset"+sf.name(),
231  sMesh.time().timeName(),
232  sMesh,
233  IOobject::NO_READ,
234  IOobject::NO_WRITE,
235  false
236  ),
237  sMesh,
238  sf.dimensions(),
240  (
241  sf.primitiveField(),
243  (
244  faceMap,
245  sMesh.nInternalFaces()
246  )
247  ),
248  patchFields
249  )
250  );
252 
253  // 2. Change the fvsPatchFields to the correct type using a mapper
254  // constructor (with reference to the now correct internal field)
255 
257  Boundary& bf = resF.boundaryFieldRef();
258 
259  forAll(bf, patchi)
260  {
261  const fvPatch& subPatch = sMesh.boundary()[patchi];
262 
263  if (patchMap[patchi] != -1)
264  {
265  // Construct addressing
266  const fvPatch& basePatch = sf.mesh().boundary()[patchMap[patchi]];
267  const label baseStart = basePatch.start();
268  const label baseSize = basePatch.size();
269 
270  labelList directAddressing(subPatch.size());
271 
272  forAll(directAddressing, i)
273  {
274  label baseFacei = faceMap[subPatch.start() + i];
275 
276  if (baseFacei >= baseStart && baseFacei < baseStart+baseSize)
277  {
278  directAddressing[i] = baseFacei-baseStart;
279  }
280  else
281  {
282  // Mapped from internal face. Do what? Leave up to
283  // patchField. This would require also to pass in
284  // original internal field so for now do as postprocessing
285  directAddressing[i] = -1;
286  }
287  }
288 
289  bf.set
290  (
291  patchi,
293  (
294  sf.boundaryField()[patchMap[patchi]],
295  subPatch,
296  resF(),
297  directFvPatchFieldMapper(directAddressing)
298  )
299  );
300  }
301 
302  // Map internal face values onto the patch elected to hold
303  // the exposed faces
304  fvsPatchField<Type>& pfld = bf[patchi];
305  const labelUList& fc = bf[patchi].patch().faceCells();
306  const labelList& own = sf.mesh().faceOwner();
307 
308  forAll(pfld, i)
309  {
310  const label baseFacei = faceMap[subPatch.start() + i];
311 
312  if (baseFacei < sf.primitiveField().size())
313  {
314  Type val = sf.internalField()[baseFacei];
315 
316  if (cellMap[fc[i]] == own[baseFacei] || !negateIfFlipped)
317  {
318  pfld[i] = val;
319  }
320  else
321  {
322  pfld[i] = flipOp()(val);
323  }
324  }
325  else
326  {
327  // Exposed face from other patch.
328  // Only possible in case of a coupled boundary
329  const label patchi = sf.mesh().boundaryMesh().whichPatch
330  (
331  baseFacei
332  );
333  const fvPatch& otherPatch = sf.mesh().boundary()[patchi];
334  const label patchFacei =
335  otherPatch.patch().whichFace(baseFacei);
336 
337  pfld[i] = sf.boundaryField()[patchi][patchFacei];
338  }
339  }
340  }
341 
342  return tresF;
343 }
344 
345 
346 template<class Type>
349 (
351 ) const
352 {
353  return interpolate
354  (
355  sf,
356  subMesh(),
357  patchMap(),
358  cellMap(),
359  faceMap()
360  );
361 }
362 
363 
364 template<class Type>
367 (
369  const pointMesh& sMesh,
370  const labelList& patchMap,
371  const labelList& pointMap
372 )
373 {
374  // 1. Create the complete field with dummy patch fields
375  PtrList<pointPatchField<Type>> patchFields(patchMap.size());
376 
377  forAll(patchFields, patchi)
378  {
379  // Set the first one by hand as it corresponds to the
380  // exposed internal faces. Additional interpolation can be put here
381  // as necessary.
382  if (patchMap[patchi] == -1)
383  {
384  patchFields.set
385  (
386  patchi,
388  (
389  sMesh.boundary()[patchi],
391  )
392  );
393  }
394  else
395  {
396  patchFields.set
397  (
398  patchi,
400  (
402  sMesh.boundary()[patchi],
404  )
405  );
406  }
407  }
408 
409  // Create the complete field from the pieces
411  (
413  (
414  IOobject
415  (
416  "subset"+pf.name(),
417  sMesh.time().timeName(),
418  sMesh.thisDb(),
419  IOobject::NO_READ,
420  IOobject::NO_WRITE,
421  false
422  ),
423  sMesh,
424  pf.dimensions(),
425  Field<Type>(pf.primitiveField(), pointMap),
426  patchFields
427  )
428  );
430 
431 
432  // 2. Change the pointPatchFields to the correct type using a mapper
433  // constructor (with reference to the now correct internal field)
434 
436  Boundary& bf = resF.boundaryFieldRef();
437 
438  forAll(bf, patchi)
439  {
440  // Set the first one by hand as it corresponds to the
441  // exposed internal faces. Additional interpolation can be put here
442  // as necessary.
443  if (patchMap[patchi] != -1)
444  {
445  // Construct addressing
446  const pointPatch& basePatch =
447  pf.mesh().boundary()[patchMap[patchi]];
448 
449  const labelList& meshPoints = basePatch.meshPoints();
450 
451  // Make addressing from mesh to patch point
452  Map<label> meshPointMap(2*meshPoints.size());
453  forAll(meshPoints, localI)
454  {
455  meshPointMap.insert(meshPoints[localI], localI);
456  }
457 
458  // Find which subpatch points originate from which patch point
459  const pointPatch& subPatch = sMesh.boundary()[patchi];
460  const labelList& subMeshPoints = subPatch.meshPoints();
461 
462  // If mapped from outside patch leave handling up to patchField
463  labelList directAddressing(subPatch.size(), -1);
464 
465  forAll(subMeshPoints, localI)
466  {
467  // Get mesh point on original mesh.
468  label meshPointi = pointMap[subMeshPoints[localI]];
469 
470  Map<label>::const_iterator iter = meshPointMap.find(meshPointi);
471 
472  if (iter != meshPointMap.end())
473  {
474  directAddressing[localI] = iter();
475  }
476  }
477 
478  bf.set
479  (
480  patchi,
482  (
483  pf.boundaryField()[patchMap[patchi]],
484  subPatch,
485  resF(),
486  directPointPatchFieldMapper(directAddressing)
487  )
488  );
489  }
490  }
491 
492  return tresF;
493 }
494 
495 
496 template<class Type>
499 (
501 ) const
502 {
503  return interpolate
504  (
505  sf,
506  pointMesh::New(subMesh()), // subsetted point mesh
507  patchMap(),
508  pointMap()
509  );
510 }
511 
512 
513 template<class Type>
516 (
518  const fvMesh& sMesh,
519  const labelList& cellMap
520 )
521 {
522  // Create the complete field from the pieces
524  (
526  (
527  IOobject
528  (
529  "subset"+df.name(),
530  sMesh.time().timeName(),
531  sMesh,
532  IOobject::NO_READ,
533  IOobject::NO_WRITE
534  ),
535  sMesh,
536  df.dimensions(),
537  Field<Type>(df, cellMap)
538  )
539  );
540 
541  return tresF;
542 }
543 
544 
545 template<class Type>
548 (
550 ) const
551 {
552  return interpolate(df, subMesh(), cellMap());
553 }
554 
555 
556 // ************************************************************************* //
Foam::surfaceFields.
Class containing functor to negate primitives. Dummy for all other types.
Definition: flipOp.H:50
#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
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
virtual label start() const
Return start label of this patch in the polyMesh face list.
Definition: fvPatch.H:150
label nInternalFaces() const
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
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
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
A calculated boundary condition for pointField.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
Constraint pointPatchField to hold values for internal face exposed by sub-setting.
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
Abstract base class for point-mesh patch fields.
const dimensionSet & dimensions() const
Return dimensions.
A List obtained as a section of another List.
Definition: SubList.H:53
Constraint fvsPatchField to hold values for internal face exposed by sub-setting. ...
Pre-declare SubField and related Field type.
Definition: Field.H:56
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:138
This boundary condition is not designed to be evaluated; it is assumed that the value is assigned via...
const objectRegistry & thisDb() const
Return database. For now is its polyMesh.
Definition: pointMesh.H:115
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
virtual label size() const
Return size.
Definition: fvPatch.H:156
const Mesh & mesh() const
Return mesh.
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.
Internal & ref()
Return a reference to the dimensioned internal field.
bool isFlux(const DimensionedField< Type, surfaceMesh > &df)
Check if surfaceField is a flux.
Definition: surfaceFields.H:55
static tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &, const fvMesh &sMesh, const labelList &patchMap, const labelList &cellMap, const labelList &faceMap)
Map volume field.
label patchi
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:56
const Time & time() const
Return time.
Definition: IOobject.C:318
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
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Foam::calculatedFvsPatchField.
A class for managing temporary objects.
Definition: PtrList.H:53
const pointBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: pointMesh.H:103
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
virtual const labelList & meshPoints() const =0
Return mesh points.
Constraint fvPatchField to hold values for internal face exposed by sub-setting.
virtual label size() const =0
Return size.
label whichFace(const label l) const
Return label of face in patch from global face label.
Definition: polyPatch.H:389
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540