pointMapper.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 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 "mapPolyMesh.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(mpm_.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 = mpm_.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 = mpm_.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 
161 Foam::pointMapper::pointMapper(const pointMesh& pMesh, const mapPolyMesh& mpm)
162 :
163  pMesh_(pMesh),
164  mpm_(mpm),
165  insertedPoints_(true),
166  direct_(false),
167  directAddrPtr_(nullptr),
168  interpolationAddrPtr_(nullptr),
169  weightsPtr_(nullptr),
170  insertedPointLabelsPtr_(nullptr)
171 {
172  // Check for possibility of direct mapping
173  if (mpm_.pointsFromPointsMap().empty())
174  {
175  direct_ = true;
176  }
177  else
178  {
179  direct_ = false;
180  }
181 
182  // Check for inserted points
183  if (direct_ && (mpm_.pointMap().empty() || min(mpm_.pointMap()) > -1))
184  {
185  insertedPoints_ = false;
186  }
187  else
188  {
189  //Check if there are inserted points with no owner
190 
191  // Make a copy of the point map, add the entries for points from points
192  // and check for left-overs
193  labelList cm(pMesh_.size(), -1);
194 
195  const List<objectMap>& cfc = mpm_.pointsFromPointsMap();
196 
197  forAll(cfc, cfcI)
198  {
199  cm[cfc[cfcI].index()] = 0;
200  }
201 
202  if (min(cm) < 0)
203  {
204  insertedPoints_ = true;
205  }
206  }
207 }
208 
209 
210 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
211 
213 {
214  clearOut();
215 }
216 
217 
218 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
219 
221 {
222  return mpm_.pointMap().size();
223 }
224 
225 
227 {
228  return mpm_.nOldPoints();
229 }
230 
231 
233 {
234  if (!direct())
235  {
237  << "Requested direct addressing for an interpolative mapper."
238  << abort(FatalError);
239  }
240 
241  if (!insertedObjects())
242  {
243  // No inserted points. Re-use pointMap
244  return mpm_.pointMap();
245  }
246  else
247  {
248  if (!directAddrPtr_)
249  {
250  calcAddressing();
251  }
252 
253  return *directAddrPtr_;
254  }
255 }
256 
257 
259 {
260  if (direct())
261  {
263  << "Requested interpolative addressing for a direct mapper."
264  << abort(FatalError);
265  }
266 
267  if (!interpolationAddrPtr_)
268  {
269  calcAddressing();
270  }
271 
272  return *interpolationAddrPtr_;
273 }
274 
275 
277 {
278  if (direct())
279  {
281  << "Requested interpolative weights for a direct mapper."
282  << abort(FatalError);
283  }
284 
285  if (!weightsPtr_)
286  {
287  calcAddressing();
288  }
289 
290  return *weightsPtr_;
291 }
292 
293 
295 {
296  if (!insertedPointLabelsPtr_)
297  {
298  if (!insertedObjects())
299  {
300  // There are no inserted points
301  insertedPointLabelsPtr_ = new labelList(0);
302  }
303  else
304  {
305  calcAddressing();
306  }
307  }
308 
309  return *insertedPointLabelsPtr_;
310 }
311 
312 
313 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
314 
315 
316 // * * * * * * * * * * * * * * * Friend Functions * * * * * * * * * * * * * //
317 
318 
319 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
320 
321 
322 // ************************************************************************* //
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
bool empty() const
Return true if the UList is empty (ie, size() is zero)
Definition: UListI.H:313
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 List< objectMap > & pointsFromPointsMap() const
Points originating from points.
Definition: mapPolyMesh.H:396
virtual const labelListList & addressing() const
Return interpolated addressing.
Definition: pointMapper.C:258
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
virtual label sizeBeforeMapping() const
Return size before mapping.
Definition: pointMapper.C:226
const labelList & insertedObjectLabels() const
Return list of inserted points.
Definition: pointMapper.C:294
virtual ~pointMapper()
Destructor.
Definition: pointMapper.C:212
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
bool insertedObjects() const
Are there any inserted points.
Definition: pointMapper.H:148
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
virtual const scalarListList & weights() const
Return interpolaion weights.
Definition: pointMapper.C:276
virtual const labelUList & directAddressing() const
Return direct addressing.
Definition: pointMapper.C:232
virtual label size() const
Return size.
Definition: pointMapper.C:220
List< scalarList > scalarListList
Definition: scalarList.H:51
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:50
virtual bool direct() const
Is the mapping direct.
Definition: pointMapper.H:126
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:94
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:61
const labelList & pointMap() const
Old point map.
Definition: mapPolyMesh.H:390
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Template functions to aid in the implementation of demand driven data.
void deleteDemandDrivenData(DataPtr &dataPtr)
label nOldPoints() const
Number of old points.
Definition: mapPolyMesh.H:363