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

surfaceIntersectionFuncs.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 Description
00025 
00026 \*---------------------------------------------------------------------------*/
00027 
00028 #include "surfaceIntersection.H"
00029 #include <meshTools/triSurfaceSearch.H>
00030 #include <triSurface/labelPairLookup.H>
00031 #include <OpenFOAM/OFstream.H>
00032 #include <OpenFOAM/HashSet.H>
00033 #include <triSurface/triSurface.H>
00034 #include <meshTools/pointIndexHit.H>
00035 #include <meshTools/octreeDataTriSurface.H>
00036 #include <meshTools/octree.H>
00037 #include <meshTools/meshTools.H>
00038 
00039 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00040 
00041 
00042 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00043 
00044 void Foam::surfaceIntersection::writeOBJ(const point& pt, Ostream& os)
00045 {
00046     os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << endl;
00047 }
00048 
00049 
00050 void Foam::surfaceIntersection::writeOBJ
00051 (
00052     const List<point>& pts,
00053     const List<edge>& edges,
00054     Ostream& os
00055 )
00056 {
00057     forAll(pts, i)
00058     {
00059         writeOBJ(pts[i], os);
00060     }
00061     forAll(edges, i)
00062     {
00063         const edge& e = edges[i];
00064 
00065         os << "l " << e.start()+1 << ' ' << e.end()+1 << endl;
00066     }
00067 }
00068 
00069 
00070 // Get minimum length of all edges connected to point
00071 Foam::scalar Foam::surfaceIntersection::minEdgeLen
00072 (
00073     const triSurface& surf,
00074     const label pointI
00075 )
00076 {
00077     const labelList& pEdges = surf.pointEdges()[pointI];
00078 
00079     scalar minLen = GREAT;
00080 
00081     forAll(pEdges, pEdgeI)
00082     {
00083         const edge& e = surf.edges()[pEdges[pEdgeI]];
00084 
00085         minLen = min(minLen, e.mag(surf.localPoints()));
00086     }
00087 
00088     return minLen;
00089 }
00090 
00091 
00092 // Get edge between fp and fp+1 on faceI.
00093 Foam::label Foam::surfaceIntersection::getEdge
00094 (
00095     const triSurface& surf,
00096     const label faceI,
00097     const label fp
00098 )
00099 {
00100     const labelledTri& f = surf.localFaces()[faceI];
00101 
00102     edge faceEdge(f[fp], f[(fp+1) % 3]);
00103 
00104     const labelList& eLabels = surf.faceEdges()[faceI];
00105 
00106     forAll(eLabels, eI)
00107     {
00108         const label edgeI = eLabels[eI];
00109 
00110         if (surf.edges()[edgeI] == faceEdge)
00111         {
00112             return edgeI;
00113         }
00114     }
00115 
00116     FatalErrorIn
00117     (
00118         "surfaceIntersection::getEdge(const triSurface&"
00119         ", const label, const label"
00120     )   << "Problem:: Cannot find edge with vertices " << faceEdge
00121         << " in face " << faceI
00122         << abort(FatalError);
00123 
00124     return -1;
00125 }
00126 
00127 
00128 // Given a map remove all consecutive duplicate elements.
00129 void Foam::surfaceIntersection::removeDuplicates
00130 (
00131     const labelList& map,
00132     labelList& elems
00133 )
00134 {
00135     bool hasDuplicate = false;
00136 
00137     label prevVertI = -1;
00138 
00139     forAll(elems, elemI)
00140     {
00141         label newVertI = map[elems[elemI]];
00142 
00143         if (newVertI == prevVertI)
00144         {
00145             hasDuplicate = true;
00146 
00147             break;
00148         }
00149         prevVertI = newVertI;
00150     }
00151 
00152     if (hasDuplicate)
00153     {
00154         // Create copy
00155         labelList oldElems(elems);
00156 
00157         label elemI = 0;
00158 
00159         // Insert first
00160         elems[elemI++] = map[oldElems[0]];
00161 
00162         for(label vertI = 1; vertI < oldElems.size(); vertI++)
00163         {
00164             // Insert others only if they differ from one before
00165             label newVertI = map[oldElems[vertI]];
00166 
00167             if (newVertI != elems[elems.size()-1])
00168             {
00169                 elems[elemI++] = newVertI;
00170             }
00171         }
00172         elems.setSize(elemI);
00173     }
00174 }
00175 
00176 
00177 // Remap.
00178 void Foam::surfaceIntersection::inlineRemap
00179 (
00180     const labelList& map,
00181     labelList& elems
00182 )
00183 {
00184     forAll(elems, elemI)
00185     {
00186         elems[elemI] = map[elems[elemI]];
00187     }
00188 }
00189 
00190 
00191 // Remove all duplicate and degenerate elements. Return unique elements and
00192 // map from old to new.
00193 Foam::edgeList Foam::surfaceIntersection::filterEdges
00194 (
00195     const edgeList& edges,
00196     labelList& map
00197 )
00198 {
00199     HashSet<edge, Hash<edge> > uniqueEdges(10*edges.size());
00200 
00201     edgeList newEdges(edges.size());
00202 
00203     map.setSize(edges.size());
00204     map = -1;
00205 
00206     label newEdgeI = 0;
00207 
00208     forAll(edges, edgeI)
00209     {
00210         const edge& e = edges[edgeI];
00211 
00212         if
00213         (
00214             (e.start() != e.end())
00215          && (uniqueEdges.find(e) == uniqueEdges.end())
00216         )
00217         {
00218             // Edge is -non degenerate and -not yet seen.
00219             uniqueEdges.insert(e);
00220 
00221             map[edgeI] = newEdgeI;
00222 
00223             newEdges[newEdgeI++] = e;
00224         }
00225     }
00226 
00227     newEdges.setSize(newEdgeI);
00228 
00229     return newEdges;
00230 }
00231 
00232 
00233 // Remove all duplicate elements.
00234 Foam::labelList Foam::surfaceIntersection::filterLabels
00235 (
00236     const labelList& elems,
00237     labelList& map
00238 )
00239 {
00240     labelHashSet uniqueElems(10*elems.size());
00241 
00242     labelList newElems(elems.size());
00243 
00244     map.setSize(elems.size());
00245     map = -1;
00246 
00247     label newElemI = 0;
00248 
00249     forAll(elems, elemI)
00250     {
00251         label elem = elems[elemI];
00252 
00253         if (uniqueElems.find(elem) == uniqueElems.end())
00254         {
00255             // First time elem is seen
00256             uniqueElems.insert(elem);
00257 
00258             map[elemI] = newElemI;
00259 
00260             newElems[newElemI++] = elem;
00261         }
00262     }
00263 
00264     newElems.setSize(newElemI);
00265 
00266     return newElems;
00267 }
00268 
00269 
00270 void Foam::surfaceIntersection::writeIntersectedEdges
00271 (
00272     const triSurface& surf,
00273     const labelListList& edgeCutVerts,
00274     Ostream& os
00275 ) const
00276 {
00277     // Dump all points (surface followed by cutPoints)
00278     const pointField& pts = surf.localPoints();
00279 
00280     forAll(pts, pointI)
00281     {
00282         writeOBJ(pts[pointI], os);
00283     }
00284     forAll(cutPoints(), cutPointI)
00285     {
00286         writeOBJ(cutPoints()[cutPointI], os);
00287     }
00288 
00289     forAll(edgeCutVerts, edgeI)
00290     {
00291         const labelList& extraVerts = edgeCutVerts[edgeI];
00292 
00293         if (extraVerts.size())
00294         {
00295             const edge& e = surf.edges()[edgeI];
00296 
00297             // Start of original edge to first extra point
00298             os  << "l " << e.start()+1 << ' '
00299                 << extraVerts[0] + surf.nPoints() + 1 << endl;
00300 
00301             for(label i = 1; i < extraVerts.size(); i++)
00302             {
00303                 os  << "l " << extraVerts[i-1] + surf.nPoints() + 1  << ' '
00304                     << extraVerts[i] + surf.nPoints() + 1 << endl;
00305             }
00306 
00307             os  << "l " << extraVerts[extraVerts.size()-1] + surf.nPoints() + 1
00308                 << ' ' << e.end()+1 << endl;
00309         }
00310     }
00311 }
00312 
00313 
00314 // Return 0 (p close to start), 1(close to end) or -1.
00315 Foam::label Foam::surfaceIntersection::classify
00316 (
00317     const scalar startTol,
00318     const scalar endTol,
00319     const point& p,
00320     const edge& e,
00321     const pointField& points
00322 )
00323 {
00324     if (mag(p - points[e.start()]) < startTol)
00325     {
00326         return 0;
00327     }
00328     else if (mag(p - points[e.end()]) < endTol)
00329     {
00330         return 1;
00331     }
00332     else
00333     {
00334         return -1;
00335     }
00336 }
00337 
00338 
00339 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines