OpenFOAM logo
Open Source CFD Toolkit

edgeVertex.H

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-2005 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 the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     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, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     edgeVertex
00027 
00028 Description
00029     Combines edge or vertex in single label. Used to specify cuts across
00030     cell circumference.
00031 
00032 SourceFiles
00033 
00034 \*---------------------------------------------------------------------------*/
00035 
00036 #ifndef edgeVertex_H
00037 #define edgeVertex_H
00038 
00039 #include "label.H"
00040 #include "polyMesh.H"
00041 
00042 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00043 
00044 namespace Foam
00045 {
00046 
00047 // Class forward declarations
00048 class refineCell;
00049 
00050 /*---------------------------------------------------------------------------*\
00051                            Class edgeVertex Declaration
00052 \*---------------------------------------------------------------------------*/
00053 
00054 class edgeVertex
00055 {
00056     // Private data
00057 
00058         //- Reference to mesh. (could be primitive mesh but keeping polyMesh
00059         //  here saves storing reference at higher levels where we do need it)
00060         const polyMesh& mesh_;
00061 
00062     // Private Member Functions
00063 
00064         //- Disallow default bitwise copy construct
00065         edgeVertex(const edgeVertex&);
00066 
00067         //- Disallow default bitwise assignment
00068         void operator=(const edgeVertex&);
00069 
00070 
00071 public:
00072 
00073     // Static Functions
00074 
00075         //- Update refine list from map. Used to update cell/face labels
00076         //  after morphing
00077         static void updateLabels(const labelList& map, List<refineCell>&);
00078 
00079         //- Update map from map. Used to update cell/face labels
00080         //  after morphing
00081         static void updateLabels(const labelList& map, Map<label>&);
00082 
00083         //- Update map from map. Used to update cell/face labels
00084         //  after morphing
00085         static void updateLabels(const labelList& map, labelHashSet&);
00086 
00087 
00088 
00089     // Constructors
00090 
00091         //- Construct from mesh
00092         edgeVertex(const polyMesh& mesh)
00093         :
00094             mesh_(mesh)
00095         {}
00096 
00097 
00098     // Member Functions
00099 
00100         const polyMesh& mesh() const
00101         {
00102             return mesh_;
00103         }
00104 
00105 
00106     // EdgeVertex handling
00107 
00108         //- is eVert an edge?
00109         bool isEdge(const label eVert) const
00110         {
00111             if (eVert < 0 || eVert >= (mesh_.nPoints() + mesh_.nEdges()))
00112             {
00113                 FatalErrorIn("edgeVertex::isEdge(const label)")
00114                     << "EdgeVertex " << eVert << " out of range "
00115                     << mesh_.nPoints() << " to "
00116                     << (mesh_.nPoints() + mesh_.nEdges() - 1)
00117                     << abort(FatalError);
00118             }
00119             
00120             return eVert >= mesh_.nPoints();
00121         }
00122 
00123         //- convert eVert to edge label
00124         label getEdge(const label eVert) const
00125         {
00126             if (!isEdge(eVert))
00127             {
00128                 FatalErrorIn("edgeVertex::getEdge(const label)")
00129                     << "EdgeVertex " << eVert << " not an edge"
00130                     << abort(FatalError);
00131             }
00132             return eVert - mesh_.nPoints();
00133         }
00134 
00135         //- convert eVert to vertex label
00136         label getVertex(const label eVert) const
00137         {
00138             if (isEdge(eVert) || (eVert < 0))
00139             {
00140                 FatalErrorIn("edgeVertex::getVertex(const label)")
00141                     << "EdgeVertex " << eVert << " not a vertex"
00142                     << abort(FatalError);
00143             }            
00144             return eVert;
00145         }
00146 
00147         //- Convert pointI to eVert
00148         label vertToEVert(const label vertI) const
00149         {
00150             if ((vertI < 0) || (vertI >= mesh_.nPoints()))
00151             {
00152                 FatalErrorIn("edgeVertex::vertToEVert(const label)")
00153                     << "Illegal vertex number " << vertI
00154                     << abort(FatalError);
00155             }
00156             return vertI;
00157         }
00158 
00159         //- Convert edgeI to eVert
00160         label edgeToEVert(const label edgeI) const
00161         {
00162             if ((edgeI < 0) || (edgeI >= mesh_.nEdges()))
00163             {
00164                 FatalErrorIn("edgeVertex::edgeToEVert(const label)")
00165                     << "Illegal edge number " << edgeI
00166                     << abort(FatalError);
00167             }
00168             return mesh_.nPoints() + edgeI;
00169         }
00170 
00171         //- Return coordinate of cut (uses weight if edgeCut)
00172         point coord(const label cut, const scalar weight) const;
00173 
00174         //- Find mesh edge (or -1) between two cuts. 
00175         label cutPairToEdge(const label cut0, const label cut1) const;
00176 
00177         //- Write cut description to Ostream
00178         Ostream& writeCut(Ostream& os, const label cut, const scalar) const;
00179 
00180         //- Write cut descriptions to Ostream
00181         Ostream& writeCuts(Ostream& os, const labelList&, const scalarField&)
00182          const;
00183 };
00184 
00185 
00186 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00187 
00188 } // End namespace Foam
00189 
00190 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00191 
00192 #endif
00193 
00194 // ************************************************************************* //
For further information go to www.openfoam.org