edgeVertex.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) 2018-2019 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
27\*---------------------------------------------------------------------------*/
28
29#include "edgeVertex.H"
30#include "meshTools.H"
31#include "refineCell.H"
32
33// * * * * * * * * * * * * * * * Static Functions * * * * * * * * * * * * * //
34
35// Update stored refine list using map
37(
38 const labelList& map,
39 List<refineCell>& refCells
40)
41{
42 label newRefI = 0;
43
44 forAll(refCells, refI)
45 {
46 const refineCell& refCell = refCells[refI];
47
48 label oldCelli = refCell.cellNo();
49
50 label newCelli = map[oldCelli];
51
52 if (newCelli != -1)
53 {
54 refCells[newRefI++] = refineCell(newCelli, refCell.direction());
55 }
56 }
57 refCells.setSize(newRefI);
58}
59
60
61// Update stored cell numbers using map.
62// Do in two passes to prevent allocation if nothing changed.
64(
65 const labelList& map,
66 Map<label>& cellPairs
67)
68{
69 // Iterate over map to see if anything changed
70 bool changed = false;
71
72 forAllConstIters(cellPairs, iter)
73 {
74 label newMaster = map[iter.key()];
75
76 label newSlave = -1;
77
78 if (iter.val() != -1)
79 {
80 newSlave = map[iter.val()];
81 }
82
83 if ((newMaster != iter.key()) || (newSlave != iter.val()))
84 {
85 changed = true;
86
87 break;
88 }
89 }
90
91 // Relabel (use second Map to prevent overlapping)
92 if (changed)
93 {
94 Map<label> newCellPairs(2*cellPairs.size());
95
96 forAllConstIters(cellPairs, iter)
97 {
98 label newMaster = map[iter.key()];
99
100 label newSlave = -1;
101
102 if (iter.val() != -1)
103 {
104 newSlave = map[iter.val()];
105 }
106
107 if (newMaster == -1)
108 {
110 << "master cell:" << iter.key()
111 << " has disappeared" << endl;
112 }
113 else
114 {
115 newCellPairs.insert(newMaster, newSlave);
116 }
117 }
118
119 cellPairs = newCellPairs;
120 }
121}
122
123
124// Update stored cell numbers using map.
125// Do in two passes to prevent allocation if nothing changed.
127(
128 const labelList& map,
130)
131{
132 // Iterate over map to see if anything changed
133 bool changed = false;
134
135 for (const label celli : cells)
136 {
137 const label newCelli = map[celli];
138
139 if (newCelli != celli)
140 {
141 changed = true;
142
143 break;
144 }
145 }
146
147 // Relabel (use second Map to prevent overlapping)
148 if (changed)
149 {
150 labelHashSet newCells(2*cells.size());
151
152 for (const label celli : cells)
153 {
154 const label newCelli = map[celli];
155
156 if (newCelli != -1)
157 {
158 newCells.insert(newCelli);
159 }
160 }
161
162 cells = newCells;
163 }
164}
165
166
167// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
168
170(
171 const primitiveMesh& mesh,
172 const label cut,
173 const scalar weight
174)
175{
176 const pointField& pts = mesh.points();
177
178 if (isEdge(mesh, cut))
179 {
180 const edge& e = mesh.edges()[getEdge(mesh, cut)];
181
182 return weight*pts[e.end()] + (1-weight)*pts[e.start()];
183 }
184 else
185 {
186 return pts[getVertex(mesh, cut)];
187 }
188}
189
190
192(
193 const primitiveMesh& mesh,
194 const label cut0,
195 const label cut1
196)
197{
198 if (!isEdge(mesh, cut0) && !isEdge(mesh, cut1))
199 {
201 (
202 mesh,
203 getVertex(mesh, cut0),
204 getVertex(mesh, cut1)
205 );
206 }
207 else
208 {
209 return -1;
210 }
211}
212
213
215(
216 Ostream& os,
217 const label cut,
218 const scalar weight
219) const
220{
221 if (isEdge(cut))
222 {
223 label edgeI = getEdge(cut);
224
225 const edge& e = mesh().edges()[edgeI];
226
227 os << "edge:" << edgeI << e << ' ' << coord(cut, weight);
228 }
229 else
230 {
231 label vertI = getVertex(cut);
232
233 os << "vertex:" << vertI << ' ' << coord(cut, weight);
234 }
235 return os;
236}
237
238
240(
241 Ostream& os,
242 const labelList& cuts,
243 const scalarField& weights
244) const
245{
246 forAll(cuts, cutI)
247 {
248 if (cutI > 0)
249 {
250 os << ' ';
251 }
252 writeCut(os, cuts[cutI], weights[cutI]);
253 }
254 return os;
255}
256
257
258// ************************************************************************* //
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Definition: HashTableI.H:180
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
void setSize(const label n)
Alias for resize()
Definition: List.H:218
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
iterator end() noexcept
Return an iterator to end traversing the UList.
Definition: UListI.H:350
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Ostream & writeCuts(Ostream &os, const labelList &, const scalarField &) const
Write cut descriptions to Ostream.
Definition: edgeVertex.C:240
static void updateLabels(const labelList &map, List< refineCell > &)
Update refine list from map. Used to update cell/face labels.
Definition: edgeVertex.C:37
Ostream & writeCut(Ostream &os, const label cut, const scalar) const
Write cut description to Ostream.
Definition: edgeVertex.C:215
static label cutPairToEdge(const primitiveMesh &, const label cut0, const label cut1)
Find mesh edge (or -1) between two cuts.
Definition: edgeVertex.C:192
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1083
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:79
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:57
const vector & direction() const
Definition: refineCell.H:88
label cellNo() const
Definition: refineCell.H:83
Patchify triangles based on orientation w.r.t other (triangulated or triangulatable) surfaces.
const coordinateSystem & coord() const noexcept
Return const access to the coordinate system.
dynamicFvMesh & mesh
OBJstream os(runTime.globalPath()/outputName)
const cellShapeList & cells
#define WarningInFunction
Report a warning using Foam::Warning.
label findEdge(const edgeList &edges, const labelList &candidates, const label v0, const label v1)
Return edge among candidates that uses the two vertices.
Definition: meshTools.C:359
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278