All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
pointMapper.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-2022 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 "pointMapper.H"
27 #include "demandDrivenData.H"
28 #include "pointMesh.H"
29 #include "polyTopoChangeMap.H"
30 
31 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
32 
33 void Foam::pointMapper::calcAddressing() const
34 {
35  if
36  (
37  directAddrPtr_
38  || interpolationAddrPtr_
39  || weightsPtr_
40  || insertedPointLabelsPtr_
41  )
42  {
44  << "Addressing already calculated."
45  << abort(FatalError);
46  }
47 
48  if (direct())
49  {
50  // Direct addressing, no weights
51 
52  directAddrPtr_ = new labelList(map_.pointMap());
53  labelList& directAddr = *directAddrPtr_;
54 
55  // Not necessary to resize the list as there are no retired points
56  // directAddr.setSize(pMesh_.size());
57 
58  insertedPointLabelsPtr_ = new labelList(pMesh_.size());
59  labelList& insertedPoints = *insertedPointLabelsPtr_;
60 
61  label nInsertedPoints = 0;
62 
63  forAll(directAddr, pointi)
64  {
65  if (directAddr[pointi] < 0)
66  {
67  // Found inserted point
68  directAddr[pointi] = 0;
69  insertedPoints[nInsertedPoints] = pointi;
70  nInsertedPoints++;
71  }
72  }
73 
74  insertedPoints.setSize(nInsertedPoints);
75  }
76  else
77  {
78  // Interpolative addressing
79 
80  interpolationAddrPtr_ = new labelListList(pMesh_.size());
81  labelListList& addr = *interpolationAddrPtr_;
82 
83  weightsPtr_ = new scalarListList(pMesh_.size());
84  scalarListList& w = *weightsPtr_;
85 
86  // Points created from other points (i.e. points merged into it).
87  const List<objectMap>& cfc = map_.pointsFromPointsMap();
88 
89  forAll(cfc, cfcI)
90  {
91  // Get addressing
92  const labelList& mo = cfc[cfcI].masterObjects();
93 
94  label pointi = cfc[cfcI].index();
95 
96  if (addr[pointi].size())
97  {
99  << "Master point " << pointi
100  << " mapped from points " << mo
101  << " already destination of mapping." << abort(FatalError);
102  }
103 
104  // Map from masters, uniform weights
105  addr[pointi] = mo;
106  w[pointi] = scalarList(mo.size(), 1.0/mo.size());
107  }
108 
109 
110  // Do mapped points. Note that can already be set from pointsFromPoints
111  // so check if addressing size still zero.
112 
113  const labelList& cm = map_.pointMap();
114 
115  forAll(cm, pointi)
116  {
117  if (cm[pointi] > -1 && addr[pointi].empty())
118  {
119  // Mapped from a single point
120  addr[pointi] = labelList(1, cm[pointi]);
121  w[pointi] = scalarList(1, 1.0);
122  }
123  }
124 
125  // Grab inserted points (for them the size of addressing is still zero)
126 
127  insertedPointLabelsPtr_ = new labelList(pMesh_.size());
128  labelList& insertedPoints = *insertedPointLabelsPtr_;
129 
130  label nInsertedPoints = 0;
131 
132  forAll(addr, pointi)
133  {
134  if (addr[pointi].empty())
135  {
136  // Mapped from a dummy point. Take point 0 with weight 1.
137  addr[pointi] = labelList(1, label(0));
138  w[pointi] = scalarList(1, 1.0);
139 
140  insertedPoints[nInsertedPoints] = pointi;
141  nInsertedPoints++;
142  }
143  }
144 
145  insertedPoints.setSize(nInsertedPoints);
146  }
147 }
148 
149 
150 void Foam::pointMapper::clearOut()
151 {
152  deleteDemandDrivenData(directAddrPtr_);
153  deleteDemandDrivenData(interpolationAddrPtr_);
154  deleteDemandDrivenData(weightsPtr_);
155  deleteDemandDrivenData(insertedPointLabelsPtr_);
156 }
157 
158 
159 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
160 
162 (
163  const pointMesh& pMesh,
164  const polyTopoChangeMap& map
165 )
166 :
167  pMesh_(pMesh),
168  map_(map),
169  insertedPoints_(true),
170  direct_(false),
171  directAddrPtr_(nullptr),
172  interpolationAddrPtr_(nullptr),
173  weightsPtr_(nullptr),
174  insertedPointLabelsPtr_(nullptr)
175 {
176  // Check for possibility of direct mapping
177  if (map_.pointsFromPointsMap().empty())
178  {
179  direct_ = true;
180  }
181  else
182  {
183  direct_ = false;
184  }
185 
186  // Check for inserted points
187  if (direct_ && (map_.pointMap().empty() || min(map_.pointMap()) > -1))
188  {
189  insertedPoints_ = false;
190  }
191  else
192  {
193  // Check if there are inserted points with no owner
194 
195  // Make a copy of the point map, add the entries for points from points
196  // and check for left-overs
197  labelList cm(pMesh_.size(), -1);
198 
199  const List<objectMap>& cfc = map_.pointsFromPointsMap();
200 
201  forAll(cfc, cfcI)
202  {
203  cm[cfc[cfcI].index()] = 0;
204  }
205 
206  if (min(cm) < 0)
207  {
208  insertedPoints_ = true;
209  }
210  }
211 }
212 
213 
214 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
215 
217 {
218  clearOut();
219 }
220 
221 
222 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
223 
225 {
226  return map_.nOldPoints();
227 }
228 
229 
231 {
232  if (!direct())
233  {
235  << "Requested direct addressing for an interpolative mapper."
236  << abort(FatalError);
237  }
238 
239  if (!insertedObjects())
240  {
241  // No inserted points. Re-use pointMap
242  return map_.pointMap();
243  }
244  else
245  {
246  if (!directAddrPtr_)
247  {
248  calcAddressing();
249  }
250 
251  return *directAddrPtr_;
252  }
253 }
254 
255 
257 {
258  if (direct())
259  {
261  << "Requested interpolative addressing for a direct mapper."
262  << abort(FatalError);
263  }
264 
265  if (!interpolationAddrPtr_)
266  {
267  calcAddressing();
268  }
269 
270  return *interpolationAddrPtr_;
271 }
272 
273 
275 {
276  if (direct())
277  {
279  << "Requested interpolative weights for a direct mapper."
280  << abort(FatalError);
281  }
282 
283  if (!weightsPtr_)
284  {
285  calcAddressing();
286  }
287 
288  return *weightsPtr_;
289 }
290 
291 
293 {
294  if (!insertedPointLabelsPtr_)
295  {
296  if (!insertedObjects())
297  {
298  // There are no inserted points
299  insertedPointLabelsPtr_ = new labelList(0);
300  }
301  else
302  {
303  calcAddressing();
304  }
305  }
306 
307  return *insertedPointLabelsPtr_;
308 }
309 
310 
311 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
312 
313 
314 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
315 
316 
317 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
318 
319 
320 // ************************************************************************* //
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
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:325
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
virtual const labelListList & addressing() const
Return interpolated addressing.
Definition: pointMapper.C:256
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
virtual label sizeBeforeMapping() const
Return size before mapping.
Definition: pointMapper.C:224
const labelList & insertedObjectLabels() const
Return list of inserted points.
Definition: pointMapper.C:292
virtual ~pointMapper()
Destructor.
Definition: pointMapper.C:216
const labelList & pointMap() const
Old point map.
const List< objectMap > & pointsFromPointsMap() const
Points originating from points.
label nOldPoints() const
Number of old points.
bool insertedObjects() const
Are there any inserted points.
Definition: pointMapper.H:141
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
virtual const scalarListList & weights() const
Return interpolation weights.
Definition: pointMapper.C:274
virtual const labelUList & directAddressing() const
Return direct addressing.
Definition: pointMapper.C:230
List< scalarList > scalarListList
Definition: scalarList.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
virtual bool direct() const
Is the mapping direct.
Definition: pointMapper.H:119
List< label > labelList
A List of labels.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:131
label size() const
Return number of points.
Definition: pointMesh.H:85
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
Template functions to aid in the implementation of demand driven data.
pointMapper(const pointMesh &, const polyTopoChangeMap &map)
Construct from polyTopoChangeMap.
Definition: pointMapper.C:162
void deleteDemandDrivenData(DataPtr &dataPtr)