edgeVertex.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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  {
109  << "master cell:" << iter.key()
110  << " has disappeared" << endl;
111  }
112  else
113  {
114  newCellPairs.insert(newMaster, newSlave);
115  }
116  }
117 
118  cellPairs = newCellPairs;
119  }
120 }
121 
122 
123 // Update stored cell numbers using map.
124 // Do in two passes to prevent allocation if nothing changed.
126 (
127  const labelList& map,
128  labelHashSet& cells
129 )
130 {
131  // Iterate over map to see if anything changed
132  bool changed = false;
133 
134  forAllConstIter(labelHashSet, cells, iter)
135  {
136  const label newCelli = map[iter.key()];
137 
138  if (newCelli != iter.key())
139  {
140  changed = true;
141 
142  break;
143  }
144  }
145 
146  // Relabel (use second Map to prevent overlapping)
147  if (changed)
148  {
149  labelHashSet newCells(2*cells.size());
150 
151  forAllConstIter(labelHashSet, cells, iter)
152  {
153  const label newCelli = map[iter.key()];
154 
155  if (newCelli != -1)
156  {
157  newCells.insert(newCelli);
158  }
159  }
160 
161  cells = newCells;
162  }
163 }
164 
165 
166 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
167 
169 (
170  const primitiveMesh& mesh,
171  const label cut,
172  const scalar weight
173 )
174 {
175  const pointField& pts = mesh.points();
176 
177  if (isEdge(mesh, cut))
178  {
179  const edge& e = mesh.edges()[getEdge(mesh, cut)];
180 
181  return weight*pts[e.end()] + (1-weight)*pts[e.start()];
182  }
183  else
184  {
185  return pts[getVertex(mesh, cut)];
186  }
187 }
188 
189 
191 (
192  const primitiveMesh& mesh,
193  const label cut0,
194  const label cut1
195 )
196 {
197  if (!isEdge(mesh, cut0) && !isEdge(mesh, cut1))
198  {
199  return meshTools::findEdge
200  (
201  mesh,
202  getVertex(mesh, cut0),
203  getVertex(mesh, cut1)
204  );
205  }
206  else
207  {
208  return -1;
209  }
210 }
211 
212 
214 (
215  Ostream& os,
216  const label cut,
217  const scalar weight
218 ) const
219 {
220  if (isEdge(cut))
221  {
222  label edgeI = getEdge(cut);
223 
224  const edge& e = mesh().edges()[edgeI];
225 
226  os << "edge:" << edgeI << e << ' ' << coord(cut, weight);
227  }
228  else
229  {
230  label vertI = getVertex(cut);
231 
232  os << "vertex:" << vertI << ' ' << coord(cut, weight);
233  }
234  return os;
235 }
236 
237 
239 (
240  Ostream& os,
241  const labelList& cuts,
242  const scalarField& weights
243 ) const
244 {
245  forAll(cuts, cutI)
246  {
247  if (cutI > 0)
248  {
249  os << ' ';
250  }
251  writeCut(os, cuts[cutI], weights[cutI]);
252  }
253  return os;
254 }
255 
256 
257 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
static label cutPairToEdge(const primitiveMesh &, const label cut0, const label cut1)
Find mesh edge (or -1) between two cuts.
Definition: edgeVertex.C:191
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:74
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Ostream & writeCut(Ostream &os, const label cut, const scalar) const
Write cut description to Ostream.
Definition: edgeVertex.C:214
Ostream & writeCuts(Ostream &os, const labelList &, const scalarField &) const
Write cut descriptions to Ostream.
Definition: edgeVertex.C:239
virtual const pointField & points() const =0
Return mesh points.
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:337
fvMesh & mesh
label size() const
Return number of elements in table.
Definition: HashTableI.H:65
const vector & direction() const
Definition: refineCell.H:87
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
static void updateLabels(const labelList &map, List< refineCell > &)
Update refine list from map. Used to update cell/face labels.
Definition: edgeVertex.C:36
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:56
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
label cellNo() const
Definition: refineCell.H:82
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
void setSize(const label)
Reset size of List.
Definition: List.C:281
#define WarningInFunction
Report a warning using Foam::Warning.
label end() const
Return end vertex label.
Definition: edgeI.H:92
static point coord(const primitiveMesh &, const label cut, const scalar weight)
Return coordinate of cut (uses weight if edgeCut)
Definition: edgeVertex.C:169
const dimensionedScalar e
Elementary charge.
Definition: doubleScalar.H:105
label start() const
Return start vertex label.
Definition: edgeI.H:81