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 | Copyright (C) 2011 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 "edgeVertex.H"
27 #include "meshTools.H"
28 #include "refineCell.H"
29 
30 
31 // * * * * * * * * * * * * * * * Static Functions * * * * * * * * * * * * * //
32 
33 
34 // Update stored refine list using map
36 (
37  const labelList& map,
38  List<refineCell>& refCells
39 )
40 {
41  label newRefI = 0;
42 
43  forAll(refCells, refI)
44  {
45  const refineCell& refCell = refCells[refI];
46 
47  label oldCellI = refCell.cellNo();
48 
49  label newCellI = map[oldCellI];
50 
51  if (newCellI != -1)
52  {
53  refCells[newRefI++] = refineCell(newCellI, refCell.direction());
54  }
55  }
56  refCells.setSize(newRefI);
57 }
58 
59 
60 // Update stored cell numbers using map.
61 // Do in two passes to prevent allocation if nothing changed.
63 (
64  const labelList& map,
65  Map<label>& cellPairs
66 )
67 {
68  // Iterate over map to see if anything changed
69  bool changed = false;
70 
71  forAllConstIter(Map<label>, cellPairs, iter)
72  {
73  label newMaster = map[iter.key()];
74 
75  label newSlave = -1;
76 
77  if (iter() != -1)
78  {
79  newSlave = map[iter()];
80  }
81 
82  if ((newMaster != iter.key()) || (newSlave != iter()))
83  {
84  changed = true;
85 
86  break;
87  }
88  }
89 
90  // Relabel (use second Map to prevent overlapping)
91  if (changed)
92  {
93  Map<label> newCellPairs(2*cellPairs.size());
94 
95  forAllConstIter(Map<label>, cellPairs, iter)
96  {
97  label newMaster = map[iter.key()];
98 
99  label newSlave = -1;
100 
101  if (iter() != -1)
102  {
103  newSlave = map[iter()];
104  }
105 
106  if (newMaster == -1)
107  {
108  WarningIn
109  (
110  "edgeVertex::updateLabels(const labelList&, "
111  "Map<label>&)"
112  ) << "master cell:" << iter.key()
113  << " has disappeared" << endl;
114  }
115  else
116  {
117  newCellPairs.insert(newMaster, newSlave);
118  }
119  }
120 
121  cellPairs = newCellPairs;
122  }
123 }
124 
125 
126 // Update stored cell numbers using map.
127 // Do in two passes to prevent allocation if nothing changed.
129 (
130  const labelList& map,
131  labelHashSet& cells
132 )
133 {
134  // Iterate over map to see if anything changed
135  bool changed = false;
136 
137  forAllConstIter(labelHashSet, cells, iter)
138  {
139  const label newCellI = map[iter.key()];
140 
141  if (newCellI != iter.key())
142  {
143  changed = true;
144 
145  break;
146  }
147  }
148 
149  // Relabel (use second Map to prevent overlapping)
150  if (changed)
151  {
152  labelHashSet newCells(2*cells.size());
153 
154  forAllConstIter(labelHashSet, cells, iter)
155  {
156  const label newCellI = map[iter.key()];
157 
158  if (newCellI != -1)
159  {
160  newCells.insert(newCellI);
161  }
162  }
163 
164  cells = newCells;
165  }
166 }
167 
168 
169 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
170 
172 (
173  const primitiveMesh& mesh,
174  const label cut,
175  const scalar weight
176 )
177 {
178  const pointField& pts = mesh.points();
179 
180  if (isEdge(mesh, cut))
181  {
182  const edge& e = mesh.edges()[getEdge(mesh, cut)];
183 
184  return weight*pts[e.end()] + (1-weight)*pts[e.start()];
185  }
186  else
187  {
188  return pts[getVertex(mesh, cut)];
189  }
190 }
191 
192 
194 (
195  const primitiveMesh& mesh,
196  const label cut0,
197  const label cut1
198 )
199 {
200  if (!isEdge(mesh, cut0) && !isEdge(mesh, cut1))
201  {
202  return meshTools::findEdge
203  (
204  mesh,
205  getVertex(mesh, cut0),
206  getVertex(mesh, cut1)
207  );
208  }
209  else
210  {
211  return -1;
212  }
213 }
214 
215 
217 (
218  Ostream& os,
219  const label cut,
220  const scalar weight
221 ) const
222 {
223  if (isEdge(cut))
224  {
225  label edgeI = getEdge(cut);
226 
227  const edge& e = mesh().edges()[edgeI];
228 
229  os << "edge:" << edgeI << e << ' ' << coord(cut, weight);
230  }
231  else
232  {
233  label vertI = getVertex(cut);
234 
235  os << "vertex:" << vertI << ' ' << coord(cut, weight);
236  }
237  return os;
238 }
239 
240 
242 (
243  Ostream& os,
244  const labelList& cuts,
245  const scalarField& weights
246 ) const
247 {
248  forAll(cuts, cutI)
249  {
250  if (cutI > 0)
251  {
252  os << ' ';
253  }
254  writeCut(os, cuts[cutI], weights[cutI]);
255  }
256  return os;
257 }
258 
259 
260 // ************************************************************************* //
virtual const pointField & points() const =0
Return mesh points.
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:52
label cellNo() const
Definition: refineCell.H:78
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
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
dynamicFvMesh & mesh
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:74
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
const dimensionedScalar e
Elementary charge.
Definition: doubleFloat.H:78
#define forAll(list, i)
Definition: UList.H:421
Ostream & writeCuts(Ostream &os, const labelList &, const scalarField &) const
Write cut descriptions to Ostream.
Definition: edgeVertex.C:242
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
static void updateLabels(const labelList &map, List< refineCell > &)
Update refine list from map. Used to update cell/face labels.
Definition: edgeVertex.C:36
static point coord(const primitiveMesh &, const label cut, const scalar weight)
Return coordinate of cut (uses weight if edgeCut)
Definition: edgeVertex.C:172
const vector & direction() const
Definition: refineCell.H:83
Ostream & writeCut(Ostream &os, const label cut, const scalar) const
Write cut description to Ostream.
Definition: edgeVertex.C:217
label end() const
Return end vertex label.
Definition: edgeI.H:92
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
label start() const
Return start vertex label.
Definition: edgeI.H:81
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:343
static label cutPairToEdge(const primitiveMesh &, const label cut0, const label cut1)
Find mesh edge (or -1) between two cuts.
Definition: edgeVertex.C:194