FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 Description
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "edgeVertex.H"
29 #include <meshTools/meshTools.H>
30 #include <dynamicMesh/refineCell.H>
31 
32 
33 // * * * * * * * * * * * * * * * Static Functions * * * * * * * * * * * * * //
34 
35 
36 // Update stored refine list using map
38 (
39  const labelList& map,
40  List<refineCell>& refCells
41 )
42 {
43  label newRefI = 0;
44 
45  forAll(refCells, refI)
46  {
47  const refineCell& refCell = refCells[refI];
48 
49  label oldCellI = refCell.cellNo();
50 
51  label newCellI = map[oldCellI];
52 
53  if (newCellI != -1)
54  {
55  refCells[newRefI++] = refineCell(newCellI, refCell.direction());
56  }
57  }
58  refCells.setSize(newRefI);
59 }
60 
61 
62 // Update stored cell numbers using map.
63 // Do in two passes to prevent allocation if nothing changed.
65 (
66  const labelList& map,
67  Map<label>& cellPairs
68 )
69 {
70  // Iterate over map to see if anything changed
71  bool changed = false;
72 
73  for
74  (
75  Map<label>::const_iterator iter = cellPairs.begin();
76  iter != cellPairs.end();
77  ++iter
78  )
79  {
80  label newMaster = map[iter.key()];
81 
82  label newSlave = -1;
83 
84  if (iter() != -1)
85  {
86  newSlave = map[iter()];
87  }
88 
89  if ((newMaster != iter.key()) || (newSlave != iter()))
90  {
91  changed = true;
92 
93  break;
94  }
95  }
96 
97  // Relabel (use second Map to prevent overlapping)
98  if (changed)
99  {
100  Map<label> newCellPairs(2*cellPairs.size());
101 
102  for
103  (
104  Map<label>::const_iterator iter = cellPairs.begin();
105  iter != cellPairs.end();
106  ++iter
107  )
108  {
109  label newMaster = map[iter.key()];
110 
111  label newSlave = -1;
112 
113  if (iter() != -1)
114  {
115  newSlave = map[iter()];
116  }
117 
118  if (newMaster == -1)
119  {
120  WarningIn
121  (
122  "edgeVertex::updateLabels(const labelList&, "
123  "Map<label>&)"
124  ) << "master cell:" << iter.key()
125  << " has disappeared" << endl;
126  }
127  else
128  {
129  newCellPairs.insert(newMaster, newSlave);
130  }
131  }
132 
133  cellPairs = newCellPairs;
134  }
135 }
136 
137 
138 // Update stored cell numbers using map.
139 // Do in two passes to prevent allocation if nothing changed.
141 (
142  const labelList& map,
144 )
145 {
146  // Iterate over map to see if anything changed
147  bool changed = false;
148 
149  for
150  (
151  labelHashSet::const_iterator iter = cells.begin();
152  iter != cells.end();
153  ++iter
154  )
155  {
156  label newCellI = map[iter.key()];
157 
158  if (newCellI != iter.key())
159  {
160  changed = true;
161 
162  break;
163  }
164  }
165 
166  // Relabel (use second Map to prevent overlapping)
167  if (changed)
168  {
169  labelHashSet newCells(2*cells.size());
170 
171  for
172  (
173  labelHashSet::const_iterator iter = cells.begin();
174  iter != cells.end();
175  ++iter
176  )
177  {
178  label newCellI = map[iter.key()];
179 
180  if (newCellI != -1)
181  {
182  newCells.insert(newCellI);
183  }
184  }
185 
186  cells = newCells;
187  }
188 }
189 
190 
191 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
192 
194 (
195  const primitiveMesh& mesh,
196  const label cut,
197  const scalar weight
198 )
199 {
200  const pointField& pts = mesh.points();
201 
202  if (isEdge(mesh, cut))
203  {
204  const edge& e = mesh.edges()[getEdge(mesh, cut)];
205 
206  return weight*pts[e.end()] + (1-weight)*pts[e.start()];
207  }
208  else
209  {
210  return pts[getVertex(mesh, cut)];
211  }
212 }
213 
214 
216 (
217  const primitiveMesh& mesh,
218  const label cut0,
219  const label cut1
220 )
221 {
222  if (!isEdge(mesh, cut0) && !isEdge(mesh, cut1))
223  {
224  return meshTools::findEdge
225  (
226  mesh,
227  getVertex(mesh, cut0),
228  getVertex(mesh, cut1)
229  );
230  }
231  else
232  {
233  return -1;
234  }
235 }
236 
237 
239 (
240  Ostream& os,
241  const label cut,
242  const scalar weight
243 ) const
244 {
245  if (isEdge(cut))
246  {
247  label edgeI = getEdge(cut);
248 
249  const edge& e = mesh().edges()[edgeI];
250 
251  os << "edge:" << edgeI << e << ' ' << coord(cut, weight);
252  }
253  else
254  {
255  label vertI = getVertex(cut);
256 
257  os << "vertex:" << vertI << ' ' << coord(cut, weight);
258  }
259  return os;
260 }
261 
262 
264 (
265  Ostream& os,
266  const labelList& cuts,
267  const scalarField& weights
268 ) const
269 {
270  forAll(cuts, cutI)
271  {
272  if (cutI > 0)
273  {
274  os << ' ';
275  }
276  writeCut(os, cuts[cutI], weights[cutI]);
277  }
278  return os;
279 }
280 
281 
282 // ************************ vim: set sw=4 sts=4 et: ************************ //