fvSurfaceMapper.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-2019 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 Description
25  FV surface mapper.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "fvSurfaceMapper.H"
30 #include "fvMesh.H"
31 #include "mapPolyMesh.H"
32 #include "faceMapper.H"
33 
34 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
35 
36 void Foam::fvSurfaceMapper::calcAddressing() const
37 {
38  if
39  (
40  directAddrPtr_
41  || interpolationAddrPtr_
42  || weightsPtr_
43  || insertedObjectLabelsPtr_
44  )
45  {
47  << "Addressing already calculated"
48  << abort(FatalError);
49  }
50 
51  // Mapping
52 
53  const label oldNInternal = faceMap_.nOldInternalFaces();
54 
55  // Assemble the maps
56  if (direct())
57  {
58  // Direct mapping - slice to size
59  directAddrPtr_ =
60  new labelList
61  (
63  (
64  faceMap_.directAddressing(),
65  mesh_.nInternalFaces()
66  )
67  );
68  labelList& addr = *directAddrPtr_;
69 
70  // Adjust for creation of an internal face from a boundary face
71  forAll(addr, facei)
72  {
73  if (addr[facei] > oldNInternal)
74  {
75  addr[facei] = 0;
76  }
77  }
78  }
79  else
80  {
81  // Interpolative mapping - slice to size
82  interpolationAddrPtr_ =
83  new labelListList
84  (
86  (
87  faceMap_.addressing(),
88  mesh_.nInternalFaces()
89  )
90  );
91  labelListList& addr = *interpolationAddrPtr_;
92 
93  weightsPtr_ =
94  new scalarListList
95  (
97  (
98  faceMap_.weights(),
99  mesh_.nInternalFaces()
100  )
101  );
102  scalarListList& w = *weightsPtr_;
103 
104  // Adjust for creation of an internal face from a boundary face
105  forAll(addr, facei)
106  {
107  if (max(addr[facei]) >= oldNInternal)
108  {
109  addr[facei] = labelList(1, label(0));
110  w[facei] = scalarList(1, 1.0);
111  }
112  }
113  }
114 
115  // Inserted objects
116 
117  // If there are, assemble the labels
118  if (insertedObjects())
119  {
120  const labelList& insFaces = faceMap_.insertedObjectLabels();
121 
122  insertedObjectLabelsPtr_ = new labelList(insFaces.size());
123  labelList& ins = *insertedObjectLabelsPtr_;
124 
125  label nIns = 0;
126 
127  forAll(insFaces, facei)
128  {
129  // If the face is internal, keep it here
130  if (insFaces[facei] < mesh_.nInternalFaces())
131  {
132  ins[nIns] = insFaces[facei];
133  nIns++;
134  }
135  }
136 
137  ins.setSize(nIns);
138  }
139  else
140  {
141  // No inserted objects
142  insertedObjectLabelsPtr_ = new labelList(0);
143  }
144 }
145 
146 
147 void Foam::fvSurfaceMapper::clearOut()
148 {
149  deleteDemandDrivenData(directAddrPtr_);
150  deleteDemandDrivenData(interpolationAddrPtr_);
151  deleteDemandDrivenData(weightsPtr_);
152 
153  deleteDemandDrivenData(insertedObjectLabelsPtr_);
154 }
155 
156 
157 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
158 
160 (
161  const fvMesh& mesh,
162  const faceMapper& fMapper
163 )
164 :
165  mesh_(mesh),
166  faceMap_(fMapper),
167  directAddrPtr_(nullptr),
168  interpolationAddrPtr_(nullptr),
169  weightsPtr_(nullptr),
170  insertedObjectLabelsPtr_(nullptr)
171 {}
172 
173 
174 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
175 
177 {
178  clearOut();
179 }
180 
181 
182 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
183 
185 {
186  if (!direct())
187  {
189  << "Requested direct addressing for an interpolative mapper."
190  << abort(FatalError);
191  }
192 
193  if (!directAddrPtr_)
194  {
195  calcAddressing();
196  }
197 
198  return *directAddrPtr_;
199 }
200 
201 
203 {
204  if (direct())
205  {
207  << "Requested interpolative addressing for a direct mapper."
208  << abort(FatalError);
209  }
210 
211  if (!interpolationAddrPtr_)
212  {
213  calcAddressing();
214  }
215 
216  return *interpolationAddrPtr_;
217 }
218 
219 
221 {
222  if (direct())
223  {
225  << "Requested interpolative weights for a direct mapper."
226  << abort(FatalError);
227  }
228 
229  if (!weightsPtr_)
230  {
231  calcAddressing();
232  }
233 
234  return *weightsPtr_;
235 }
236 
237 
239 {
240  if (!insertedObjectLabelsPtr_)
241  {
242  calcAddressing();
243  }
244 
245  return *insertedObjectLabelsPtr_;
246 }
247 
248 
249 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
250 
251 
252 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
253 
254 
255 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
256 
257 
258 // ************************************************************************* //
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#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
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual bool insertedObjects() const
Are there any inserted faces.
This object provides mapping and fill-in information for face data between the two meshes after the t...
Definition: faceMapper.H:55
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
label nInternalFaces() const
fvSurfaceMapper(const fvMesh &mesh, const faceMapper &fMapper)
Construct from components.
virtual bool direct() const
Is the mapping direct.
virtual const labelList & insertedObjectLabels() const
Return list of inserted faces.
Definition: faceMapper.C:357
SubList< label > subList
Declare type of subList.
Definition: List.H:199
virtual ~fvSurfaceMapper()
Destructor.
virtual const labelListList & addressing() const
Return interpolated addressing.
Definition: faceMapper.C:321
virtual const labelUList & directAddressing() const
Return direct addressing.
List< scalarList > scalarListList
Definition: scalarList.H:51
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
virtual const scalarListList & weights() const
Return interpolaion weights.
Definition: faceMapper.C:339
List< label > labelList
A List of labels.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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 const labelList & insertedObjectLabels() const
Return list of inserted faces.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual const scalarListList & weights() const
Return interpolaion weights.
virtual label nOldInternalFaces() const
Return number of old internalFaces.
Definition: faceMapper.C:382
void deleteDemandDrivenData(DataPtr &dataPtr)
virtual const labelUList & directAddressing() const
Return direct addressing.
Definition: faceMapper.C:295
virtual const labelListList & addressing() const
Return interpolated addressing.