Merge branch 'master' of ssh://git.code.sf.net/p/foam-extend/foam-extend-3.2
[foam-extend-3.2.git] / src / mesh / cfMesh / utilities / meshes / polyMeshGen2DEngine / polyMeshGen2DEngine.C
blobc510e626ee5f4db8c8ff6d197235081f2e11e49a
1 /*---------------------------------------------------------------------------*\
2   =========                 |
3   \\      /  F ield         | cfMesh: A library for mesh generation
4    \\    /   O peration     |
5     \\  /    A nd           | Author: Franjo Juretic (franjo.juretic@c-fields.com)
6      \\/     M anipulation  | Copyright (C) Creative Fields, Ltd.
7 -------------------------------------------------------------------------------
8 License
9     This file is part of cfMesh.
11     cfMesh is free software; you can redistribute it and/or modify it
12     under the terms of the GNU General Public License as published by the
13     Free Software Foundation; either version 3 of the License, or (at your
14     option) any later version.
16     cfMesh 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.
21     You should have received a copy of the GNU General Public License
22     along with cfMesh.  If not, see <http://www.gnu.org/licenses/>.
24 Description
26 \*---------------------------------------------------------------------------*/
28 #include "polyMeshGen2DEngine.H"
29 #include "polyMeshGenAddressing.H"
30 #include "demandDrivenData.H"
32 # ifdef USE_OMP
33 #include <omp.h>
34 # endif
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 namespace Foam
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 void polyMeshGen2DEngine::findActiveFaces() const
45     const faceListPMG& faces = mesh_.faces();
46     const boolList& zMinPoints = this->zMinPoints();
47     const boolList& zMaxPoints = this->zMaxPoints();
49     activeFacePtr_ = new boolList(faces.size());
51     # ifdef USE_OMP
52     # pragma omp parallel for schedule(dynamic, 50)
53     # endif
54     forAll(faces, faceI)
55     {
56         bool hasZMin(false), hasZMax(false);
58         const face& f = faces[faceI];
59         forAll(f, pI)
60         {
61             hasZMin |= zMinPoints[f[pI]];
62             hasZMax |= zMaxPoints[f[pI]];
63         }
65         activeFacePtr_->operator[](faceI) = (hasZMin && hasZMax);
66     }
69 void polyMeshGen2DEngine::findActiveFaceLabels() const
71     const boolList& activeFace = this->activeFace();
73     label counter(0);
75     forAll(activeFace, faceI)
76         if( activeFace[faceI] )
77             ++counter;
79     activeFaceLabelsPtr_ = new labelList(counter);
81     counter = 0;
82     forAll(activeFace, faceI)
83         if( activeFace[faceI] )
84             activeFaceLabelsPtr_->operator[](counter++) = faceI;
87 void polyMeshGen2DEngine::findZMinPoints() const
89     const pointFieldPMG& points = mesh_.points();
91     zMinPointPtr_ = new boolList(points.size());
93     const scalar tZ = 0.05 * (bb_.max().z() - bb_.min().z());;
95     # ifdef USE_OMP
96     # pragma omp parallel for schedule(dynamic, 50)
97     # endif
98     forAll(points, pointI)
99     {
100         if( Foam::mag(points[pointI].z() - bb_.min().z()) < tZ )
101         {
102             zMinPointPtr_->operator[](pointI) = true;
103         }
104         else
105         {
106             zMinPointPtr_->operator[](pointI) = false;
107         }
108     }
111 void polyMeshGen2DEngine::findZMinPointLabels() const
113     const boolList& zMinPoints = this->zMinPoints();
115     label counter(0);
117     forAll(zMinPoints, pointI)
118         if( zMinPoints[pointI] )
119             ++counter;
121     if( 2 * counter != zMinPoints.size() )
122     {
123         FatalErrorIn
124         (
125             "void polyMeshGen2DEngine::findZMinPointLabels()"
126         ) << "The number of points at smallest z coordinate is"
127           << " not half of the total number of points."
128           << " This is not a 2D mesh or is not aligned with the z axis"
129           << exit(FatalError);
130     }
132     zMinPointLabelsPtr_ = new labelList(counter);
134     counter = 0;
135     forAll(zMinPoints, pointI)
136         if( zMinPoints[pointI] )
137             zMinPointLabelsPtr_->operator[](counter++) = pointI;
140 void polyMeshGen2DEngine::findZMinOffsetPoints() const
142     const boolList& zMinPoints = this->zMinPoints();
143     const labelList& zMinPointLabels = this->zMinPointLabels();
145     zMinToZMaxPtr_ = new labelList(zMinPointLabels.size());
147     const VRWGraph& pointPoints = mesh_.addressingData().pointPoints();
148     # ifdef USE_OMP
149     # pragma omp parallel for schedule(dynamic, 50)
150     # endif
151     forAll(zMinPointLabels, pI)
152     {
153         const label pointI = zMinPointLabels[pI];
155         label nInactive(0), offsetPoint(-1);
156         forAllRow(pointPoints, pointI, ppI)
157         {
158             if( !zMinPoints[pointPoints(pointI, ppI)] )
159             {
160                 ++nInactive;
161                 offsetPoint = pointPoints(pointI, ppI);
162             }
163         }
165         if( nInactive == 1 )
166         {
167             zMinToZMaxPtr_->operator[](pI) = offsetPoint;
168         }
169         else
170         {
171             FatalErrorIn
172             (
173                 "void polyMeshGen2DEngine::findZMinOffsetPoints()"
174             ) << "This cannot be a 2D mesh" << exit(FatalError);
175         }
176     }
179 void polyMeshGen2DEngine::findZMaxPoints() const
181     const pointFieldPMG& points = mesh_.points();
183     zMaxPointPtr_ = new boolList(points.size());
185     const scalar tZ = 0.05 * (bb_.max().z() - bb_.min().z());
187     # ifdef USE_OMP
188     # pragma omp parallel for schedule(dynamic, 50)
189     # endif
190     forAll(points, pointI)
191     {
192         if( Foam::mag(points[pointI].z() - bb_.max().z()) < tZ )
193         {
194             zMaxPointPtr_->operator[](pointI) = true;
195         }
196         else
197         {
198             zMaxPointPtr_->operator[](pointI) = false;
199         }
200     }
203 void polyMeshGen2DEngine::findZMaxPointLabels() const
205     const boolList& zMaxPoints = this->zMaxPoints();
207     label counter(0);
209     forAll(zMaxPoints, pointI)
210         if( zMaxPoints[pointI] )
211             ++counter;
213     if( 2 * counter != zMaxPoints.size() )
214     {
215         FatalErrorIn
216         (
217             "void polyMeshGen2DEngine::findZMaxPointLabels()"
218         ) << "The number of points at largest z coordinate is"
219           << " not half of the total number of points."
220           << " This is not a 2D mesh or is not aligned with the z axis"
221           << exit(FatalError);
222     }
224     zMaxPointLabelsPtr_ = new labelList(counter);
226     counter = 0;
227     forAll(zMaxPoints, pointI)
228         if( zMaxPoints[pointI] )
229             zMaxPointLabelsPtr_->operator[](counter++) = pointI;
232 void polyMeshGen2DEngine::findZMaxOffsetPoints() const
234     const boolList& zMaxPoints = this->zMaxPoints();
235     const labelList& zMaxPointLabels = this->zMaxPointLabels();
237     zMaxToZMinPtr_ = new labelList(zMaxPointLabels.size());
239     const VRWGraph& pointPoints = mesh_.addressingData().pointPoints();
240     # ifdef USE_OMP
241     # pragma omp parallel for schedule(dynamic, 50)
242     # endif
243     forAll(zMaxPointLabels, pI)
244     {
245         const label pointI = zMaxPointLabels[pI];
247         label nInactive(0), offsetPoint(-1);
248         forAllRow(pointPoints, pointI, ppI)
249         {
250             if( !zMaxPoints[pointPoints(pointI, ppI)] )
251             {
252                 ++nInactive;
253                 offsetPoint = pointPoints(pointI, ppI);
254             }
255         }
257         if( nInactive == 1 )
258         {
259             zMaxToZMinPtr_->operator[](pI) = offsetPoint;
260         }
261         else
262         {
263             FatalErrorIn
264             (
265                 "void polyMeshGen2DEngine::findZMaxOffsetPoints()"
266             ) << "This cannot be a 2D mesh" << exit(FatalError);
267         }
268     }
271 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
273 polyMeshGen2DEngine::polyMeshGen2DEngine(const polyMeshGen& mesh)
275     mesh_(mesh),
276     bb_(),
277     activeFacePtr_(NULL),
278     activeFaceLabelsPtr_(NULL),
279     zMinPointPtr_(NULL),
280     zMinPointLabelsPtr_(NULL),
281     zMinToZMaxPtr_(NULL),
282     zMaxPointPtr_(NULL),
283     zMaxPointLabelsPtr_(NULL),
284     zMaxToZMinPtr_(NULL)
286     const pointFieldPMG& points = mesh_.points();
288     bb_.min() = point(VGREAT, VGREAT, VGREAT);
289     bb_.max() = point(-VGREAT, -VGREAT, -VGREAT);
291     # ifdef USE_OMP
292     # pragma omp parallel
293     # endif
294     {
295         point localMin(VGREAT, VGREAT, VGREAT);
296         point localMax(-VGREAT, -VGREAT, -VGREAT);
298         # ifdef USE_OMP
299         # pragma omp for schedule(dynamic, 50)
300         # endif
301         forAll(points, pointI)
302         {
303             localMin = Foam::min(localMin, points[pointI]);
304             localMax = Foam::max(localMax, points[pointI]);
305         }
307         # ifdef USE_OMP
308         # pragma omp critical
309         # endif
310         {
311             bb_.min() = Foam::min(bb_.min(), localMin);
312             bb_.max() = Foam::max(bb_.max(), localMax);
313         }
314     }
316     if( Pstream::parRun() )
317     {
318         reduce(bb_.min(), minOp<point>());
319         reduce(bb_.max(), maxOp<point>());
320     }
323 polyMeshGen2DEngine::~polyMeshGen2DEngine()
325     clearOut();
328 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
330 void polyMeshGen2DEngine::correctPoints()
332     pointFieldPMG& points = const_cast<pointFieldPMG&>(mesh_.points());
334     const labelList& zMinPointLabels = this->zMinPointLabels();
335     const labelList& zMinOffset = this->zMinToZMax();
337     # ifdef USE_OMP
338     # pragma omp parallel for schedule(dynamic, 50)
339     # endif
340     forAll(zMinPointLabels, apI)
341     {
342         point& p = points[zMinPointLabels[apI]];
343         point& op = points[zMinOffset[apI]];
345         op.x() = p.x();
346         op.y() = p.y();
347         p.z() = bb_.min().z();
348         op.z() = bb_.max().z();
349     }
352 void polyMeshGen2DEngine::clearOut()
354     deleteDemandDrivenData(activeFacePtr_);
355     deleteDemandDrivenData(activeFaceLabelsPtr_);
356     deleteDemandDrivenData(zMinPointPtr_);
357     deleteDemandDrivenData(zMinPointLabelsPtr_);
358     deleteDemandDrivenData(zMinToZMaxPtr_);
359     deleteDemandDrivenData(zMaxPointPtr_);
360     deleteDemandDrivenData(zMaxPointLabelsPtr_);
361     deleteDemandDrivenData(zMaxToZMinPtr_);
364 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
366 } // End namespace Foam
368 // ************************************************************************* //