src/dynamicMesh/meshCut/directions/directionInfo/directionInfo.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 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
  27 #include "directionInfo.H"
  28 #include "polyMesh.H"
  29
  30 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  31
  32 // Find edge among edgeLabels that uses v0 and v1
  33 Foam::label Foam::directionInfo::findEdge
  34 (
  35     const primitiveMesh& mesh,
  36     const labelList& edgeLabels,
  37     const label v1,
  38     const label v0
  39 )
  40 {
  41     forAll(edgeLabels, edgeLabelI)
  42     {
  43         label edgeI = edgeLabels[edgeLabelI];
  44
  45         if (mesh.edges()[edgeI] == edge(v0, v1))
  46         {
  47             return edgeI;
  48         }
  49     }
  50
  51     FatalErrorIn("directionInfo::findEdge")
  52         << "Cannot find an edge among " << edgeLabels << endl
  53         << "that uses vertices " << v0
  54         << " and " << v1
  55         << abort(FatalError);
  56
  57     return -1;
  58 }
  59
  60
  61 Foam::label Foam::directionInfo::lowest
  62 (
  63     const label size,
  64     const label a,
  65     const label b
  66 )
  67 {
  68     // Get next point
  69     label a1 = (a + 1) % size;
  70
  71     if (a1 == b)
  72     {
  73         return a;
  74     }
  75     else
  76     {
  77         label b1 = (b + 1) % size;
  78
  79         if (b1 != a)
  80         {
  81             FatalErrorIn("directionInfo::lowest")
  82                 << "Problem : a:" << a << " b:" << b << " size:" << size
  83                 << abort(FatalError);
  84         }
  85
  86         return b;
  87     }
  88 }
  89
  90
  91 // Have edge on hex cell. Find corresponding edge on face. Return -1 if none
  92 // found.
  93 Foam::label Foam::directionInfo::edgeToFaceIndex
  94 (
  95     const primitiveMesh& mesh,
  96     const label cellI,
  97     const label faceI,
  98     const label edgeI
  99 )
 100 {
 101     if ((edgeI < 0) || (edgeI >= mesh.nEdges()))
 102     {
 103         FatalErrorIn("directionInfo::edgeToFaceIndex")
 104             << "Illegal edge label:" << edgeI
 105             << " when projecting cut edge from cell " << cellI
 106             << " to face " << faceI
 107             << abort(FatalError);
 108     }
 109
 110     const edge& e = mesh.edges()[edgeI];
 111
 112     const face& f = mesh.faces()[faceI];
 113
 114     // edgeI is either
 115     // - in faceI. Convert into index in face.
 116     // - connected (but not in) to face. Return -1.
 117     // - in face opposite faceI. Convert into index in face.
 118
 119     label fpA = findIndex(f, e.start());
 120     label fpB = findIndex(f, e.end());
 121
 122     if (fpA != -1)
 123     {
 124         if (fpB != -1)
 125         {
 126             return lowest(f.size(), fpA, fpB);
 127         }
 128         else
 129         {
 130             // e.start() in face, e.end() not
 131             return -1;
 132         }
 133     }
 134     else
 135     {
 136         if (fpB != -1)
 137         {
 138             // e.end() in face, e.start() not
 139             return -1;
 140         }
 141         else
 142         {
 143             // Both not in face.
 144             // e is on opposite face. Determine corresponding edge on this face:
 145             // - determine two faces using edge (one is the opposite face,
 146             //   one is 'side' face
 147             // - walk on both these faces to opposite edge
 148             // - check if this opposite edge is on faceI
 149
 150             label f0I, f1I;
 151
 152             meshTools::getEdgeFaces(mesh, cellI, edgeI, f0I, f1I);
 153
 154             // Walk to opposite edge on face f0
 155             label edge0I =
 156                 meshTools::walkFace(mesh, f0I, edgeI, e.start(), 2);
 157
 158             // Check if edge on faceI.
 159
 160             const edge& e0 = mesh.edges()[edge0I];
 161
 162             fpA = findIndex(f, e0.start());
 163             fpB = findIndex(f, e0.end());
 164
 165             if ((fpA != -1) && (fpB != -1))
 166             {
 167                 return lowest(f.size(), fpA, fpB);
 168             }
 169
 170             // Face0 is doesn't have an edge on faceI (so must be the opposite
 171             // face) so try face1.
 172
 173             // Walk to opposite edge on face f1
 174             label edge1I =
 175                 meshTools::walkFace(mesh, f1I, edgeI, e.start(), 2);
 176
 177             // Check if edge on faceI.
 178             const edge& e1 = mesh.edges()[edge1I];
 179
 180             fpA = findIndex(f, e1.start());
 181             fpB = findIndex(f, e1.end());
 182
 183             if ((fpA != -1) && (fpB != -1))
 184             {
 185                 return lowest(f.size(), fpA, fpB);
 186             }
 187
 188             FatalErrorIn("directionInfo::edgeToFaceIndex")
 189                 << "Found connected faces " << mesh.faces()[f0I] << " and "
 190                 << mesh.faces()[f1I] << " sharing edge " << edgeI << endl
 191                 << "But none seems to be connected to face " << faceI
 192                 << " vertices:" << f
 193                 << abort(FatalError);
 194
 195             return -1;
 196         }
 197     }
 198 }
 199
 200
 201 // * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
 202
 203 Foam::Ostream& Foam::operator<<
 204 (
 205     Foam::Ostream& os,
 206     const Foam::directionInfo& wDist
 207 )
 208 {
 209     if (os.format() == IOstream::ASCII)
 210     {
 211         os << wDist.index_ << wDist.n_;
 212     }
 213     else
 214     {
 215         os.write
 216         (
 217             reinterpret_cast<const char*>(&wDist.index_),
 218             sizeof(directionInfo)
 219         );
 220     }
 221
 222     // Check state of Ostream
 223     os.check("Ostream& operator<<(Ostream&, const directionInfo&)");
 224     return os;
 225
 226 }
 227
 228
 229 Foam::Istream& Foam::operator>>(Foam::Istream& is, Foam::directionInfo& wDist)
 230 {
 231     if (is.format() == IOstream::ASCII)
 232     {
 233         is >> wDist.index_ >> wDist.n_;
 234     }
 235     else
 236     {
 237         is.read
 238         (
 239             reinterpret_cast<char*>(&wDist.index_),
 240             sizeof(directionInfo)
 241         );
 242     }
 243
 244     // Check state of Istream
 245     is.check("Istream& operator>>(Istream&, directionInfo&)");
 246     return is;
 247 }
 248
 249 // ************************************************************************* //