FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

PatchToolsSortEdges.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "PatchTools.H"
00027 #include <OpenFOAM/SortableList.H>
00028 #include <OpenFOAM/transform.H>
00029 
00030 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00031 
00032 template
00033 <
00034     class Face,
00035     template<class> class FaceList,
00036     class PointField,
00037     class PointType
00038 >
00039 
00040 Foam::labelListList
00041 Foam::PatchTools::sortedEdgeFaces
00042 (
00043     const PrimitivePatch<Face, FaceList, PointField, PointType>& p
00044 )
00045 {
00046     const edgeList& edges = p.edges();
00047     const labelListList& edgeFaces = p.edgeFaces();
00048     const List<Face>& localFaces = p.localFaces();
00049     const Field<PointType>& localPoints = p.localPoints();
00050 
00051     // create the lists for the various results. (resized on completion)
00052     labelListList sortedEdgeFaces = labelListList(edgeFaces.size());
00053 
00054     forAll(edgeFaces, edgeI)
00055     {
00056         const labelList& faceNbs = edgeFaces[edgeI];
00057 
00058         if (faceNbs.size() > 2)
00059         {
00060             // Get point on edge and normalized direction of edge (= e2 base
00061             // of our coordinate system)
00062             const edge& e = edges[edgeI];
00063 
00064             const point& edgePt = localPoints[e.start()];
00065 
00066             vector e2 = e.vec(localPoints);
00067             e2 /= mag(e2) + VSMALL;
00068 
00069             // Get opposite vertex for 0th face
00070             const Face& f = localFaces[faceNbs[0]];
00071 
00072             label fp0 = findIndex(f, e[0]);
00073             label fp1 = f.fcIndex(fp0);
00074             label vertI = (f[fp1] != e[1] ? f[fp1] : f.fcIndex(fp1));
00075 
00076             // Get vector normal both to e2 and to edge from opposite vertex
00077             // to edge (will be x-axis of our coordinate system)
00078             vector e0 = e2 ^ (localPoints[vertI] - edgePt);
00079             e0 /= mag(e0) + VSMALL;
00080 
00081             // Get y-axis of coordinate system
00082             vector e1 = e2 ^ e0;
00083 
00084             SortableList<scalar> faceAngles(faceNbs.size());
00085 
00086             // e0 is reference so angle is 0
00087             faceAngles[0] = 0;
00088 
00089             for (label nbI = 1; nbI < faceNbs.size(); nbI++)
00090             {
00091                 // Get opposite vertex
00092                 const Face& f = localFaces[faceNbs[nbI]];
00093                 label fp0 = findIndex(f, e[0]);
00094                 label fp1 = f.fcIndex(fp0);
00095                 label vertI = (f[fp1] != e[1] ? f[fp1] : f.fcIndex(fp1));
00096 
00097                 vector vec = e2 ^ (localPoints[vertI] - edgePt);
00098                 vec /= mag(vec) + VSMALL;
00099 
00100                 faceAngles[nbI] = pseudoAngle
00101                 (
00102                     e0,
00103                     e1,
00104                     vec
00105                 );
00106             }
00107 
00108             faceAngles.sort();
00109 
00110             sortedEdgeFaces[edgeI] = UIndirectList<label>
00111             (
00112                 faceNbs,
00113                 faceAngles.indices()
00114             );
00115         }
00116         else
00117         {
00118             // No need to sort. Just copy.
00119             sortedEdgeFaces[edgeI] = faceNbs;
00120         }
00121     }
00122 
00123     return sortedEdgeFaces;
00124 }
00125 
00126 
00127 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines