hexRef8Data.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) 2015-2017 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 "IOobject.H"
27 #include "UList.H"
28 
29 #include "hexRef8Data.H"
30 #include "mapPolyMesh.H"
31 #include "mapDistributePolyMesh.H"
32 #include "polyMesh.H"
33 #include "syncTools.H"
34 #include "refinementHistory.H"
35 #include "fvMesh.H"
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
39 Foam::hexRef8Data::hexRef8Data(const IOobject& io)
40 {
41  {
42  IOobject rio(io);
43  rio.rename("cellLevel");
44  bool haveFile = returnReduce
45  (
46  rio.typeHeaderOk<labelIOList>(true),
47  orOp<bool>()
48  );
49  if (haveFile)
50  {
51  Info<< "Reading hexRef8 data : " << rio.name() << endl;
52  cellLevelPtr_.reset(new labelIOList(rio));
53  }
54  }
55  {
56  IOobject rio(io);
57  rio.rename("pointLevel");
58  bool haveFile = returnReduce
59  (
60  rio.typeHeaderOk<labelIOList>(true),
61  orOp<bool>()
62  );
63  if (haveFile)
64  {
65  Info<< "Reading hexRef8 data : " << rio.name() << endl;
66  pointLevelPtr_.reset(new labelIOList(rio));
67  }
68  }
69  {
70  IOobject rio(io);
71  rio.rename("level0Edge");
72  bool haveFile = returnReduce
73  (
75  orOp<bool>()
76  );
77  if (haveFile)
78  {
79  Info<< "Reading hexRef8 data : " << rio.name() << endl;
80  level0EdgePtr_.reset(new uniformDimensionedScalarField(rio));
81  }
82  }
83  {
84  IOobject rio(io);
85  rio.rename("refinementHistory");
86  bool haveFile = returnReduce
87  (
89  orOp<bool>()
90  );
91  if (haveFile)
92  {
93  Info<< "Reading hexRef8 data : " << rio.name() << endl;
94  refHistoryPtr_.reset(new refinementHistory(rio));
95  }
96  }
97 }
98 
99 
100 Foam::hexRef8Data::hexRef8Data
101 (
102  const IOobject& io,
103  const hexRef8Data& data,
104  const labelList& cellMap,
105  const labelList& pointMap
106 )
107 {
108  if (data.cellLevelPtr_.valid())
109  {
110  IOobject rio(io);
111  rio.rename(data.cellLevelPtr_().name());
112 
113  cellLevelPtr_.reset
114  (
115  new labelIOList
116  (
117  rio,
118  UIndirectList<label>(data.cellLevelPtr_(), cellMap)()
119  )
120  );
121  }
122  if (data.pointLevelPtr_.valid())
123  {
124  IOobject rio(io);
125  rio.rename(data.pointLevelPtr_().name());
126 
127  pointLevelPtr_.reset
128  (
129  new labelIOList
130  (
131  rio,
132  UIndirectList<label>(data.pointLevelPtr_(), pointMap)()
133  )
134  );
135  }
136  if (data.level0EdgePtr_.valid())
137  {
138  IOobject rio(io);
139  rio.rename(data.level0EdgePtr_().name());
140 
141  level0EdgePtr_.reset
142  (
143  new uniformDimensionedScalarField(rio, data.level0EdgePtr_())
144  );
145  }
146  if (data.refHistoryPtr_.valid())
147  {
148  IOobject rio(io);
149  rio.rename(data.refHistoryPtr_().name());
150 
151  refHistoryPtr_ = data.refHistoryPtr_().clone(rio, cellMap);
152  }
153 }
154 
155 
156 Foam::hexRef8Data::hexRef8Data
157 (
158  const IOobject& io,
159  const UPtrList<const labelList>& cellMaps,
160  const UPtrList<const labelList>& pointMaps,
161  const UPtrList<const hexRef8Data>& procDatas
162 )
163 {
164  const polyMesh& mesh = dynamic_cast<const polyMesh&>(io.db());
165 
166  // cellLevel
167 
168  if (procDatas[0].cellLevelPtr_.valid())
169  {
170  IOobject rio(io);
171  rio.rename(procDatas[0].cellLevelPtr_().name());
172 
173  cellLevelPtr_.reset(new labelIOList(rio, mesh.nCells()));
174  labelList& cellLevel = cellLevelPtr_();
175 
176  forAll(procDatas, procI)
177  {
178  const labelList& procCellLevel = procDatas[procI].cellLevelPtr_();
179  UIndirectList<label>(cellLevel, cellMaps[procI]) = procCellLevel;
180  }
181  }
182 
183 
184  // pointLevel
185 
186  if (procDatas[0].pointLevelPtr_.valid())
187  {
188  IOobject rio(io);
189  rio.rename(procDatas[0].pointLevelPtr_().name());
190 
191  pointLevelPtr_.reset(new labelIOList(rio, mesh.nPoints()));
192  labelList& pointLevel = pointLevelPtr_();
193 
194  forAll(procDatas, procI)
195  {
196  const labelList& procPointLevel = procDatas[procI].pointLevelPtr_();
197  UIndirectList<label>(pointLevel, pointMaps[procI]) = procPointLevel;
198  }
199  }
200 
201 
202  // level0Edge
203 
204  if (procDatas[0].level0EdgePtr_.valid())
205  {
206  IOobject rio(io);
207  rio.rename(procDatas[0].level0EdgePtr_().name());
208 
209  level0EdgePtr_.reset
210  (
212  (
213  rio,
214  procDatas[0].level0EdgePtr_()
215  )
216  );
217  }
218 
219 
220  // refinementHistory
221 
222  if (procDatas[0].refHistoryPtr_.valid())
223  {
224  IOobject rio(io);
225  rio.rename(procDatas[0].refHistoryPtr_().name());
226 
227  UPtrList<const refinementHistory> procRefs(procDatas.size());
228  forAll(procDatas, i)
229  {
230  procRefs.set(i, &procDatas[i].refHistoryPtr_());
231  }
232 
233  refHistoryPtr_.reset
234  (
236  (
237  rio,
238  cellMaps,
239  procRefs
240  )
241  );
242  }
243 }
244 
245 
246 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
247 
249 {}
250 
251 
252 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
253 
255 {
256  const polyMesh& mesh = dynamic_cast<const polyMesh&>(io.db());
257 
258  bool hasCellLevel = returnReduce(cellLevelPtr_.valid(), orOp<bool>());
259  if (hasCellLevel && !cellLevelPtr_.valid())
260  {
261  IOobject rio(io);
262  rio.rename("cellLevel");
263  rio.readOpt() = IOobject::NO_READ;
264  cellLevelPtr_.reset(new labelIOList(rio, labelList(mesh.nCells(), 0)));
265  }
266 
267  bool hasPointLevel = returnReduce(pointLevelPtr_.valid(), orOp<bool>());
268  if (hasPointLevel && !pointLevelPtr_.valid())
269  {
270  IOobject rio(io);
271  rio.rename("pointLevel");
272  rio.readOpt() = IOobject::NO_READ;
273  pointLevelPtr_.reset
274  (
275  new labelIOList(rio, labelList(mesh.nPoints(), 0))
276  );
277  }
278 
279  bool hasLevel0Edge = returnReduce(level0EdgePtr_.valid(), orOp<bool>());
280  if (hasLevel0Edge)
281  {
282  // Get master length
283  scalar masterLen = level0EdgePtr_().value();
284  Pstream::scatter(masterLen);
285  if (!level0EdgePtr_.valid())
286  {
287  IOobject rio(io);
288  rio.rename("level0Edge");
289  rio.readOpt() = IOobject::NO_READ;
290  level0EdgePtr_.reset
291  (
293  (
294  rio,
295  dimensionedScalar("zero", dimLength, masterLen)
296  )
297  );
298  }
299  }
300 
301  bool hasHistory = returnReduce(refHistoryPtr_.valid(), orOp<bool>());
302  if (hasHistory && !refHistoryPtr_.valid())
303  {
304  IOobject rio(io);
305  rio.rename("refinementHistory");
306  rio.readOpt() = IOobject::NO_READ;
307  refHistoryPtr_.reset(new refinementHistory(rio, mesh.nCells(), true));
308  }
309 }
310 
311 
313 {
314  if (cellLevelPtr_.valid())
315  {
316  map.cellMap().distribute(cellLevelPtr_());
317  }
318  if (pointLevelPtr_.valid())
319  {
320  map.pointMap().distribute(pointLevelPtr_());
321  }
322 
323  // No need to distribute the level0Edge
324 
325  if (refHistoryPtr_.valid() && refHistoryPtr_().active())
326  {
327  refHistoryPtr_().distribute(map);
328  }
329 }
330 
331 
333 {
334  bool ok = true;
335  if (cellLevelPtr_.valid())
336  {
337  ok = ok && cellLevelPtr_().write();
338  }
339  if (pointLevelPtr_.valid())
340  {
341  ok = ok && pointLevelPtr_().write();
342  }
343  if (level0EdgePtr_.valid())
344  {
345  ok = ok && level0EdgePtr_().write();
346  }
347  if (refHistoryPtr_.valid())
348  {
349  ok = ok && refHistoryPtr_().write();
350  }
351  return ok;
352 }
353 
354 
355 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
const mapDistribute & pointMap() const
Point distribute map.
const word & name() const
Return name.
Definition: IOobject.H:291
bool typeHeaderOk(const bool checkType=true)
Read header (uses typeFilePath to find file) and check header.
void sync(const IOobject &io)
Parallel synchronise. This enforces valid objects on all processors.
Definition: hexRef8Data.C:254
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
const mapDistribute & cellMap() const
Cell distribute map.
virtual void rename(const word &newName)
Rename.
Definition: IOobject.H:321
label nCells() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
UniformDimensionedField< scalar > uniformDimensionedScalarField
dynamicFvMesh & mesh
List< label > labelList
A List of labels.
Definition: labelList.H:56
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:54
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
bool write() const
Write.
Definition: hexRef8Data.C:332
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A List with indirect addressing.
Definition: fvMatrix.H:106
messageStream Info
label nPoints() const
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
void distribute(const mapDistributePolyMesh &)
In-place distribute.
Definition: hexRef8Data.C:312
Various for reading/decomposing/reconstructing/distributing refinement data.
Definition: hexRef8Data.H:57
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
All refinement history. Used in unrefinement.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:331
readOption readOpt() const
Definition: IOobject.H:353
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
~hexRef8Data()
Destructor.
Definition: hexRef8Data.C:248
IOList< label > labelIOList
Label container classes.
Definition: labelIOList.H:42