OpenFOAM logo
Open Source CFD Toolkit

pointEdgePointI.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 \*---------------------------------------------------------------------------*/
00026 
00027 #include "polyMesh.H"
00028 #include "transform.H"
00029 
00030 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00031 
00032 // Update this with w2 if w2 nearer to pt.
00033 inline bool Foam::pointEdgePoint::update
00034 (
00035     const point& pt,
00036     const pointEdgePoint& w2,
00037     const scalar tol
00038 )
00039 {
00040     scalar dist2 = magSqr(pt - w2.origin());
00041 
00042     if (!valid())
00043     {
00044         // current not yet set so use any value
00045         distSqr_ = dist2;
00046         origin_ = w2.origin();
00047 
00048         return true;
00049     }        
00050 
00051     scalar diff = distSqr_ - dist2;
00052 
00053     if (diff < 0)
00054     {
00055         // already nearer to pt
00056         return false;
00057     }
00058 
00059     if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
00060     {
00061         // don't propagate small changes
00062         return false;
00063     }
00064     else
00065     {
00066         // update with new values
00067         distSqr_ = dist2;
00068         origin_ = w2.origin();
00069 
00070         return true;
00071     }
00072 }
00073 
00074     
00075 // Update this with w2 (information on same point)
00076 inline bool Foam::pointEdgePoint::update
00077 (
00078     const pointEdgePoint& w2,
00079     const scalar tol
00080 )
00081 {
00082     if (!valid())
00083     {
00084         // current not yet set so use any value
00085         distSqr_ = w2.distSqr();
00086         origin_ = w2.origin();
00087 
00088         return true;
00089     }        
00090 
00091     scalar diff = distSqr_ - w2.distSqr();
00092 
00093     if (diff < 0)
00094     {
00095         // already nearer to pt
00096         return false;
00097     }
00098 
00099     if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
00100     {
00101         // don't propagate small changes
00102         return false;
00103     }
00104     else
00105     {
00106         // update with new values
00107         distSqr_ =  w2.distSqr();
00108         origin_ = w2.origin();
00109 
00110         return true;
00111     }
00112 }
00113 
00114 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00115 
00116 // Null constructor
00117 inline Foam::pointEdgePoint::pointEdgePoint()
00118 :
00119     origin_(greatPoint),
00120     distSqr_(GREAT)
00121 {}
00122 
00123 
00124 // Construct from origin, distance
00125 inline Foam::pointEdgePoint::pointEdgePoint
00126 (
00127     const point& origin, 
00128     const scalar distSqr
00129 )
00130 :
00131     origin_(origin),
00132     distSqr_(distSqr)
00133 {}
00134 
00135 
00136 // Construct as copy
00137 inline Foam::pointEdgePoint::pointEdgePoint(const pointEdgePoint& wpt)
00138 :
00139     origin_(wpt.origin()),
00140     distSqr_(wpt.distSqr())
00141 {}
00142 
00143 
00144 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00145 
00146 inline const Foam::point& Foam::pointEdgePoint::origin() const
00147 {
00148     return origin_;
00149 }
00150 
00151 
00152 inline Foam::scalar Foam::pointEdgePoint::distSqr() const
00153 {
00154     return distSqr_;
00155 }
00156 
00157 
00158 inline bool Foam::pointEdgePoint::valid() const
00159 {
00160     return origin_ != greatPoint;
00161 }
00162 
00163 
00164 // Checks for cyclic points
00165 inline bool Foam::pointEdgePoint::sameGeometry
00166 (
00167     const pointEdgePoint& w2,
00168     const scalar tol
00169 ) const
00170 {
00171     scalar diff = Foam::mag(distSqr() - w2.distSqr());
00172 
00173     if (diff < SMALL)
00174     {
00175         return true;
00176     }
00177     else
00178     {
00179         if ((distSqr() > SMALL) && ((diff/distSqr()) < tol))
00180         {
00181             return true;
00182         }
00183         else
00184         {
00185             return false;
00186         }
00187     }
00188 }
00189 
00190 
00191 inline void Foam::pointEdgePoint::leaveDomain
00192 (
00193     const polyPatch& patch,
00194     const label patchPointI,
00195     const point& coord
00196 )
00197 {
00198     origin_ -= coord;
00199 }
00200 
00201 
00202 inline void Foam::pointEdgePoint::transform(const tensor& rotTensor)
00203 {
00204     origin_ = Foam::transform(rotTensor, origin_);
00205 }
00206 
00207 
00208 // Update absolute geometric quantities. Note that distance (distSqr_)
00209 // is not affected by leaving/entering domain.
00210 inline void Foam::pointEdgePoint::enterDomain
00211 (
00212     const polyPatch& patch,
00213     const label patchPointI,
00214     const point& coord
00215 )
00216 {
00217     // back to absolute form
00218     origin_ += coord;
00219 }
00220 
00221 
00222 // Update this with information from connected edge
00223 inline bool Foam::pointEdgePoint::updatePoint
00224 (
00225     const polyMesh& mesh,
00226     const label pointI,
00227     const label edgeI,
00228     const pointEdgePoint& edgeInfo,
00229     const scalar tol
00230 )
00231 {
00232     return 
00233         update
00234         (
00235             mesh.points()[pointI],
00236             edgeInfo,
00237             tol
00238         );
00239     }    
00240 
00241 
00242 // Update this with new information on same point
00243 inline bool Foam::pointEdgePoint::updatePoint
00244 (
00245     const polyMesh& mesh,
00246     const label pointI,
00247     const pointEdgePoint& newPointInfo,
00248     const scalar tol
00249 )
00250 {
00251     return
00252         update
00253         (
00254             mesh.points()[pointI],
00255             newPointInfo,
00256             tol
00257         );
00258 }    
00259 
00260 
00261 // Update this with new information on same point. No extra information.
00262 inline bool Foam::pointEdgePoint::updatePoint
00263 (
00264     const pointEdgePoint& newPointInfo,
00265     const scalar tol
00266 )
00267 {
00268     return update(newPointInfo, tol);
00269 }    
00270 
00271 
00272 // Update this with information from connected point
00273 inline bool Foam::pointEdgePoint::updateEdge
00274 (
00275     const polyMesh& mesh,
00276     const label edgeI,
00277     const label pointI,
00278     const pointEdgePoint& pointInfo,
00279     const scalar tol
00280 )
00281 {
00282     const pointField& points = mesh.points();
00283 
00284     const edge& e = mesh.edges()[edgeI];
00285 
00286     const point edgeMid(0.5*(points[e[0]] + points[e[1]]));
00287 
00288     return
00289         update
00290         (
00291             edgeMid,
00292             pointInfo,
00293             tol
00294         );
00295 }    
00296 
00297 
00298 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
00299 
00300 inline bool Foam::pointEdgePoint::operator==(const Foam::pointEdgePoint& rhs)
00301  const
00302 {
00303     return origin() == rhs.origin();
00304 }
00305 
00306 
00307 inline bool Foam::pointEdgePoint::operator!=(const Foam::pointEdgePoint& rhs)
00308  const
00309 {
00310     return !(*this == rhs);
00311 }
00312 
00313 
00314 // ************************************************************************* //
For further information go to www.openfoam.org