src/OpenFOAM/meshes/meshTools/matchPoints.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2009 OpenCFD Ltd.
   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 "matchPoints.H"
  28 #include "SortableList.H"
  29 #include "ListOps.H"
  30
  31 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  32
  33 bool Foam::matchPoints
  34 (
  35     const UList<point>& pts0,
  36     const UList<point>& pts1,
  37     const UList<scalar>& matchDistances,
  38     const bool verbose,
  39     labelList& from0To1,
  40     const point& origin
  41 )
  42 {
  43     from0To1.setSize(pts0.size());
  44     from0To1 = -1;
  45
  46     bool fullMatch = true;
  47
  48     point compareOrigin = origin;
  49
  50     if (origin == point(VGREAT, VGREAT, VGREAT))
  51     {
  52         if (pts1.size())
  53         {
  54             compareOrigin = sum(pts1)/pts1.size();
  55         }
  56     }
  57
  58     SortableList<scalar> pts0MagSqr(magSqr(pts0 - compareOrigin));
  59
  60     SortableList<scalar> pts1MagSqr(magSqr(pts1 - compareOrigin));
  61
  62     forAll(pts0MagSqr, i)
  63     {
  64         scalar dist0 = pts0MagSqr[i];
  65
  66         label face0I = pts0MagSqr.indices()[i];
  67
  68         scalar matchDist = matchDistances[face0I];
  69
  70         label startI = findLower(pts1MagSqr, 0.99999*dist0 - matchDist);
  71
  72         if (startI == -1)
  73         {
  74             startI = 0;
  75         }
  76
  77
  78         // Go through range of equal mag and find nearest vector.
  79         scalar minDistSqr = VGREAT;
  80         label minFaceI = -1;
  81
  82         for
  83         (
  84             label j = startI;
  85             (
  86                 (j < pts1MagSqr.size())
  87              && (pts1MagSqr[j] < 1.00001*dist0 + matchDist)
  88             );
  89             j++
  90         )
  91         {
  92             label faceI = pts1MagSqr.indices()[j];
  93             // Compare actual vectors
  94             scalar distSqr = magSqr(pts0[face0I] - pts1[faceI]);
  95
  96             if (distSqr <= sqr(matchDist) && distSqr < minDistSqr)
  97             {
  98                 minDistSqr = distSqr;
  99                 minFaceI = faceI;
 100             }
 101         }
 102
 103         if (minFaceI == -1)
 104         {
 105             fullMatch = false;
 106
 107             if (verbose)
 108             {
 109                 Pout<< "Cannot find point in pts1 matching point " << face0I
 110                     << " coord:" << pts0[face0I]
 111                     << " in pts0 when using tolerance " << matchDist << endl;
 112
 113                 // Go through range of equal mag and find equal vector.
 114                 Pout<< "Searching started from:" << startI << " in pts1"
 115                     << endl;
 116                 for
 117                 (
 118                     label j = startI;
 119                     (
 120                         (j < pts1MagSqr.size())
 121                      && (pts1MagSqr[j] < 1.00001*dist0 + matchDist)
 122                     );
 123                     j++
 124                 )
 125                 {
 126                     label faceI = pts1MagSqr.indices()[j];
 127
 128                     Pout<< "Compared coord:" << pts1[faceI]
 129                         << " with difference to point "
 130                         << mag(pts1[faceI] - pts0[face0I]) << endl;
 131                 }
 132             }
 133         }
 134
 135         from0To1[face0I] = minFaceI;
 136     }
 137
 138     return fullMatch;
 139 }
 140
 141 // ************************************************************************* //