fvFieldDecomposer.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 "fvFieldDecomposer.H"
27 
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 Foam::labelList Foam::fvFieldDecomposer::patchFieldDecomposer::alignAddressing
32 (
33  const labelUList& addressingSlice,
34  const label addressingOffset
35 ) const
36 {
37  labelList addressing(addressingSlice.size());
38 
39  forAll(addressing, i)
40  {
41  // Subtract one to align addressing.
42  addressing[i] = addressingSlice[i] - (addressingOffset + 1);
43  }
44 
45  return addressing;
46 }
47 
48 
50 (
51  const labelUList& addressingSlice,
52  const label addressingOffset
53 )
54 :
55  labelList(alignAddressing(addressingSlice, addressingOffset)),
56  directFvPatchFieldMapper(static_cast<const labelList&>(*this))
57 {}
58 
59 
60 Foam::labelList Foam::fvFieldDecomposer::processorVolPatchFieldDecomposer::
61 alignAddressing
62 (
63  const fvMesh& mesh,
64  const labelUList& addressingSlice
65 ) const
66 {
67  labelList addressing(addressingSlice.size());
68 
69  const labelList& own = mesh.faceOwner();
70  const labelList& neighb = mesh.faceNeighbour();
71 
72  forAll(addressing, i)
73  {
74  // Subtract one to align addressing.
75  label ai = mag(addressingSlice[i]) - 1;
76 
77  if (ai < neighb.size())
78  {
79  // This is a regular face. it has been an internal face
80  // of the original mesh and now it has become a face
81  // on the parallel boundary.
82  // Give face the value of the neighbour.
83 
84  if (addressingSlice[i] >= 0)
85  {
86  // I have the owner so use the neighbour value
87  addressing[i] = neighb[ai];
88  }
89  else
90  {
91  addressing[i] = own[ai];
92  }
93  }
94  else
95  {
96  // This is a face that used to be on a cyclic boundary
97  // but has now become a parallel patch face. I cannot
98  // do the interpolation properly (I would need to look
99  // up the different (face) list of data), so I will
100  // just grab the value from the owner cell
101 
102  addressing[i] = own[ai];
103  }
104  }
105 
106  return addressing;
107 }
108 
109 
112 (
113  const fvMesh& mesh,
114  const labelUList& addressingSlice
115 )
116 :
117  labelList(alignAddressing(mesh, addressingSlice)),
118  directFvPatchFieldMapper(static_cast<const labelList&>(*this))
119 {}
120 
121 
123 (
124  const fvMesh& completeMesh,
125  const fvMesh& procMesh,
126  const labelList& faceAddressing,
127  const labelList& cellAddressing,
128  const labelList& boundaryAddressing
129 )
130 :
131  completeMesh_(completeMesh),
132  procMesh_(procMesh),
133  faceAddressing_(faceAddressing),
134  cellAddressing_(cellAddressing),
135  boundaryAddressing_(boundaryAddressing),
136  patchFieldDecomposerPtrs_
137  (
138  procMesh_.boundary().size(),
139  static_cast<patchFieldDecomposer*>(nullptr)
140  ),
141  processorVolPatchFieldDecomposerPtrs_
142  (
143  procMesh_.boundary().size(),
144  static_cast<processorVolPatchFieldDecomposer*>(nullptr)
145  )
146 {
147  forAll(boundaryAddressing_, patchi)
148  {
149  const fvPatch& procPatch = procMesh.boundary()[patchi];
150 
151  label fromPatchi = boundaryAddressing_[patchi];
152  if (fromPatchi < 0 && isA<processorCyclicFvPatch>(procPatch))
153  {
154  const label referPatchi =
155  refCast<const processorCyclicPolyPatch>
156  (procPatch.patch()).referPatchID();
157  fromPatchi = boundaryAddressing_[referPatchi];
158  }
159 
160  if (fromPatchi >= 0)
161  {
162  patchFieldDecomposerPtrs_[patchi] = new patchFieldDecomposer
163  (
164  procMesh_.boundary()[patchi].patchSlice(faceAddressing_),
165  completeMesh_.boundaryMesh()[fromPatchi].start()
166  );
167  }
168 
169  if (boundaryAddressing_[patchi] < 0)
170  {
171  processorVolPatchFieldDecomposerPtrs_[patchi] =
173  (
174  completeMesh_,
175  procMesh_.boundary()[patchi].patchSlice(faceAddressing_)
176  );
177  }
178  }
179 }
180 
181 
182 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
183 
185 {
186  forAll(patchFieldDecomposerPtrs_, patchi)
187  {
188  if (patchFieldDecomposerPtrs_[patchi])
189  {
190  delete patchFieldDecomposerPtrs_[patchi];
191  }
192  }
193 
194  forAll(processorVolPatchFieldDecomposerPtrs_, patchi)
195  {
196  if (processorVolPatchFieldDecomposerPtrs_[patchi])
197  {
198  delete processorVolPatchFieldDecomposerPtrs_[patchi];
199  }
200  }
201 }
202 
203 // ************************************************************************* //
#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
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1175
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
directFvPatchFieldMapper(const labelUList &addressing)
Construct given addressing.
UList< label > labelUList
Definition: UList.H:65
processorVolPatchFieldDecomposer(const fvMesh &mesh, const labelUList &addressingSlice)
Construct given addressing.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:138
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1169
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
fvFieldDecomposer(const fvMesh &completeMesh, const fvMesh &procMesh, const labelList &faceAddressing, const labelList &cellAddressing, const labelList &boundaryAddressing)
Construct from components.
virtual const labelUList & addressing() const
Access to the direct map addressing.
Processor patch field decomposer class. Maps either owner or.
patchFieldDecomposer(const labelUList &addressingSlice, const label addressingOffset)
Construct given addressing.
label patchi
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
~fvFieldDecomposer()
Destructor.
dimensioned< scalar > mag(const dimensioned< Type > &)
label size() const
Return the number of elements in the UList.
Definition: UListI.H:311
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540