src/meshTools/searchableSurface/searchablePlate.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original author
   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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "searchablePlate.H"
  28 #include "addToRunTimeSelectionTable.H"
  29 #include "SortableList.H"
  30
  31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  32
  33 namespace Foam
  34 {
  35
  36 defineTypeNameAndDebug(searchablePlate, 0);
  37 addToRunTimeSelectionTable(searchableSurface, searchablePlate, dict);
  38
  39 }
  40
  41
  42 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  43
  44 Foam::direction Foam::searchablePlate::calcNormal(const point& span)
  45 {
  46     direction normalDir = 3;
  47
  48     for (direction dir = 0; dir < vector::nComponents; dir++)
  49     {
  50         if (span[dir] < 0)
  51         {
  52             FatalErrorIn("searchablePlate::calcNormal()")
  53                 << "Span should have two positive and one zero entry. Now:"
  54                 << span << exit(FatalError);
  55         }
  56         else if (span[dir] < VSMALL)
  57         {
  58             if (normalDir == 3)
  59             {
  60                 normalDir = dir;
  61             }
  62             else
  63             {
  64                 // Multiple zero entries. Flag and exit.
  65                 normalDir = 3;
  66                 break;
  67             }
  68         }
  69     }
  70
  71     if (normalDir == 3)
  72     {
  73         FatalErrorIn("searchablePlate::calcNormal()")
  74             << "Span should have one and only zero entry. Now:" << span
  75             << exit(FatalError);
  76     }
  77
  78     return normalDir;
  79 }
  80
  81
  82 // Returns miss or hit with face (always 0)
  83 Foam::pointIndexHit Foam::searchablePlate::findNearest
  84 (
  85     const point& sample,
  86     const scalar nearestDistSqr
  87 ) const
  88 {
  89     // For every component direction can be
  90     // left of min, right of max or inbetween.
  91     // - outside points: project first one x plane (either min().x()
  92     // or max().x()), then onto y plane and finally z. You should be left
  93     // with intersection point
  94     // - inside point: find nearest side (compare to mid point). Project onto
  95     //   that.
  96
  97     // Project point on plane.
  98     pointIndexHit info(true, sample, 0);
  99     info.rawPoint()[normalDir_] = origin_[normalDir_];
 100
 101     // Clip to edges if outside
 102     for (direction dir = 0; dir < vector::nComponents; dir++)
 103     {
 104         if (dir != normalDir_)
 105         {
 106             if (info.rawPoint()[dir] < origin_[dir])
 107             {
 108                 info.rawPoint()[dir] = origin_[dir];
 109             }
 110             else if (info.rawPoint()[dir] > origin_[dir]+span_[dir])
 111             {
 112                 info.rawPoint()[dir] = origin_[dir]+span_[dir];
 113             }
 114         }
 115     }
 116
 117     // Check if outside. Optimisation: could do some checks on distance already
 118     // on components above
 119     if (magSqr(info.rawPoint() - sample) > nearestDistSqr)
 120     {
 121         info.setMiss();
 122         info.setIndex(-1);
 123     }
 124
 125     return info;
 126 }
 127
 128
 129 Foam::pointIndexHit Foam::searchablePlate::findLine
 130 (
 131     const point& start,
 132     const point& end
 133 ) const
 134 {
 135     pointIndexHit info
 136     (
 137         true,
 138         vector::zero,
 139         0
 140     );
 141
 142     const vector dir(end-start);
 143
 144     if (mag(dir[normalDir_]) < VSMALL)
 145     {
 146         info.setMiss();
 147         info.setIndex(-1);
 148     }
 149     else
 150     {
 151         scalar t = (origin_[normalDir_]-start[normalDir_]) / dir[normalDir_];
 152
 153         if (t < 0 || t > 1)
 154         {
 155             info.setMiss();
 156             info.setIndex(-1);
 157         }
 158         else
 159         {
 160             info.rawPoint() = start+t*dir;
 161             info.rawPoint()[normalDir_] = origin_[normalDir_];
 162
 163             // Clip to edges
 164             for (direction dir = 0; dir < vector::nComponents; dir++)
 165             {
 166                 if (dir != normalDir_)
 167                 {
 168                     if (info.rawPoint()[dir] < origin_[dir])
 169                     {
 170                         info.setMiss();
 171                         info.setIndex(-1);
 172                         break;
 173                     }
 174                     else if (info.rawPoint()[dir] > origin_[dir]+span_[dir])
 175                     {
 176                         info.setMiss();
 177                         info.setIndex(-1);
 178                         break;
 179                     }
 180                 }
 181             }
 182         }
 183     }
 184
 185     // Debug
 186     if (info.hit())
 187     {
 188         treeBoundBox bb(origin_, origin_+span_);
 189         bb.min()[normalDir_] -= 1E-6;
 190         bb.max()[normalDir_] += 1E-6;
 191
 192         if (!bb.contains(info.hitPoint()))
 193         {
 194             FatalErrorIn("searchablePlate::findLine(..)")
 195                 << "bb:" << bb << endl
 196                 << "origin_:" << origin_ << endl
 197                 << "span_:" << span_ << endl
 198                 << "normalDir_:" << normalDir_ << endl
 199                 << "hitPoint:" << info.hitPoint()
 200                 << abort(FatalError);
 201         }
 202     }
 203
 204     return info;
 205 }
 206
 207
 208 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 209
 210 Foam::searchablePlate::searchablePlate
 211 (
 212     const IOobject& io,
 213     const point& origin,
 214     const vector& span
 215 )
 216 :
 217     searchableSurface(io),
 218     origin_(origin),
 219     span_(span),
 220     normalDir_(calcNormal(span_))
 221 {
 222     if (debug)
 223     {
 224         Info<< "searchablePlate::searchablePlate :"
 225             << " origin:" << origin_
 226             << " origin+span:" << origin_+span_
 227             << " normal:" << vector::componentNames[normalDir_]
 228             << endl;
 229     }
 230 }
 231
 232
 233 Foam::searchablePlate::searchablePlate
 234 (
 235     const IOobject& io,
 236     const dictionary& dict
 237 )
 238 :
 239     searchableSurface(io),
 240     origin_(dict.lookup("origin")),
 241     span_(dict.lookup("span")),
 242     normalDir_(calcNormal(span_))
 243 {
 244     if (debug)
 245     {
 246         Info<< "searchablePlate::searchablePlate :"
 247             << " origin:" << origin_
 248             << " origin+span:" << origin_+span_
 249             << " normal:" << vector::componentNames[normalDir_]
 250             << endl;
 251     }
 252 }
 253
 254
 255 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 256
 257 Foam::searchablePlate::~searchablePlate()
 258 {}
 259
 260
 261 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 262
 263 const Foam::wordList& Foam::searchablePlate::regions() const
 264 {
 265     if (regions_.empty())
 266     {
 267         regions_.setSize(1);
 268         regions_[0] = "region0";
 269     }
 270     return regions_;
 271 }
 272
 273
 274 void Foam::searchablePlate::findNearest
 275 (
 276     const pointField& samples,
 277     const scalarField& nearestDistSqr,
 278     List<pointIndexHit>& info
 279 ) const
 280 {
 281     info.setSize(samples.size());
 282
 283     forAll(samples, i)
 284     {
 285         info[i] = findNearest(samples[i], nearestDistSqr[i]);
 286     }
 287 }
 288
 289
 290 void Foam::searchablePlate::findLine
 291 (
 292     const pointField& start,
 293     const pointField& end,
 294     List<pointIndexHit>& info
 295 ) const
 296 {
 297     info.setSize(start.size());
 298
 299     forAll(start, i)
 300     {
 301         info[i] = findLine(start[i], end[i]);
 302     }
 303 }
 304
 305
 306 void Foam::searchablePlate::findLineAny
 307 (
 308     const pointField& start,
 309     const pointField& end,
 310     List<pointIndexHit>& info
 311 ) const
 312 {
 313     findLine(start, end, info);
 314 }
 315
 316
 317 void Foam::searchablePlate::findLineAll
 318 (
 319     const pointField& start,
 320     const pointField& end,
 321     List<List<pointIndexHit> >& info
 322 ) const
 323 {
 324     List<pointIndexHit> nearestInfo;
 325     findLine(start, end, nearestInfo);
 326
 327     info.setSize(start.size());
 328     forAll(info, pointI)
 329     {
 330         if (nearestInfo[pointI].hit())
 331         {
 332             info[pointI].setSize(1);
 333             info[pointI][0] = nearestInfo[pointI];
 334         }
 335         else
 336         {
 337             info[pointI].clear();
 338         }
 339     }
 340 }
 341
 342
 343 void Foam::searchablePlate::getRegion
 344 (
 345     const List<pointIndexHit>& info,
 346     labelList& region
 347 ) const
 348 {
 349     region.setSize(info.size());
 350     region = 0;
 351 }
 352
 353
 354 void Foam::searchablePlate::getNormal
 355 (
 356     const List<pointIndexHit>& info,
 357     vectorField& normal
 358 ) const
 359 {
 360     normal.setSize(info.size());
 361     normal = vector::zero;
 362     forAll(normal, i)
 363     {
 364         normal[i][normalDir_] = 1.0;
 365     }
 366 }
 367
 368
 369 void Foam::searchablePlate::getVolumeType
 370 (
 371     const pointField& points,
 372     List<volumeType>& volType
 373 ) const
 374 {
 375     FatalErrorIn
 376     (
 377         "searchableCollection::getVolumeType(const pointField&"
 378         ", List<volumeType>&) const"
 379     )   << "Volume type not supported for plate."
 380         << exit(FatalError);
 381 }
 382
 383
 384 // ************************************************************************* //