calcPointCells.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Description
28 calculate point cells - ie, the cells attached to each point
29
30 - remove unused points, adjust pointCells and cellFaces accordingly
31\*---------------------------------------------------------------------------*/
32
33#include "meshReader.H"
34
35// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36
37void Foam::meshReader::calcPointCells() const
38{
39 static const label UNIT_POINT_CELLS = 12;
40
41 if (pointCellsPtr_)
42 {
44 << "pointCells already calculated"
45 << abort(FatalError);
46 }
47
48 label nPoints = points_.size();
49
50 pointCellsPtr_.reset(new labelListList(nPoints));
51 auto& ptCells = *pointCellsPtr_;
52
53 forAll(ptCells, i)
54 {
55 ptCells[i].setSize(UNIT_POINT_CELLS);
56 }
57
58 // Initialize the list of labels which will hold the count of the
59 // actual number of cells per point during the analysis
60 labelList cellCount(nPoints, Zero);
61
62 // Note. Unlike the standard point-cell algorithm, which asks the cell for
63 // the supporting point labels, we need to work based on the cell faces.
64 // This is because some of the faces do not come from the cell shape.
65 // It is also advantageous to remove duplicates from the point-cell
66 // addressing, because this removes a lot of waste later.
67
68 faceListList& cFaces = cellFaces();
69
70 // For each cell
71 forAll(cFaces, celli)
72 {
73 const faceList& faces = cFaces[celli];
74
75 forAll(faces, i)
76 {
77 // For each vertex
78 const labelList& labels = faces[i];
79
80 forAll(labels, j)
81 {
82 // Set working point label
83 label curPoint = labels[j];
84 labelList& curPointCells = ptCells[curPoint];
85 label curCount = cellCount[curPoint];
86
87 // check if the cell has been added before
88 bool found = false;
89
90 for (label f = 0; f < curCount; f++)
91 {
92 if (curPointCells[f] == celli)
93 {
94 found = true;
95 break;
96 }
97 }
98
99 if (!found)
100 {
101 // If the list of pointCells is not big enough, double it
102 if (curPointCells.size() <= curCount)
103 {
104 curPointCells.setSize(curPointCells.size()*2);
105 }
106
107 // Enter the cell label in the point's cell list
108 curPointCells[curCount] = celli;
109
110 // Increment the cell count for the point addressed
111 cellCount[curPoint]++;
112 }
113 }
114 }
115 }
116
117 // report and remove unused points
118 // - adjust points, pointCells, and cellFaces accordingly
119 label pointi = 0;
120 labelList oldToNew(nPoints, -1);
121
122 forAll(ptCells, i)
123 {
124 ptCells[i].setSize(cellCount[i]);
125 if (cellCount[i] > 0)
126 {
127 oldToNew[i] = pointi++;
128 }
129 }
130
131 // report unused points
132 if (nPoints > pointi)
133 {
134 Info<< "removing " << (nPoints - pointi) << " unused points" << endl;
135
136 nPoints = pointi;
137
138 // adjust points and truncate - bend const-ness
139 pointField& adjustedPoints = const_cast<pointField&>(points_);
140
141 inplaceReorder(oldToNew, adjustedPoints);
142 adjustedPoints.setSize(nPoints);
143
144 // adjust pointCells and truncate
145 inplaceReorder(oldToNew, ptCells);
146 ptCells.setSize(nPoints);
147
148 // adjust cellFaces - this could be faster
149 // For each cell
150 forAll(cFaces, celli)
151 {
152 faceList& faces = cFaces[celli];
153
154 // For each face
155 forAll(faces, i)
156 {
157 inplaceRenumber(oldToNew, faces[i]);
158 }
159 }
160 }
161}
162
163
165{
166 if (!pointCellsPtr_)
167 {
168 calcPointCells();
169 }
170
171 return *pointCellsPtr_;
172}
173
174
175// ************************************************************************* //
bool found
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
pointField points_
Points supporting the mesh.
Definition: meshReader.H:216
Smooth ATC in cells having a point to a set of patches supplied by type.
Definition: pointCells.H:59
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
label nPoints
void inplaceRenumber(const labelUList &oldToNew, IntListType &input)
Inplace renumber the values (not the indices) of a list.
List< label > labelList
A List of labels.
Definition: List.H:66
messageStream Info
Information stream (stdout output on master, null elsewhere)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
void inplaceReorder(const labelUList &oldToNew, ListType &input, const bool prune=false)
Inplace reorder the elements of a list.
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:56
errorManip< error > abort(error &err)
Definition: errorManip.H:144
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
List< faceList > faceListList
A List of faceList.
Definition: faceListFwd.H:49
labelList f(nPoints)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333