src/libengrid/checkforoverlap.cpp

   1 //
   2 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   3 // +                                                                      +
   4 // + This file is part of enGrid.                                         +
   5 // +                                                                      +
   6 // + Copyright 2008-2012 enGits GmbH                                     +
   7 // +                                                                      +
   8 // + enGrid is free software: you can redistribute it and/or modify       +
   9 // + it under the terms of the GNU General Public License as published by +
  10 // + the Free Software Foundation, either version 3 of the License, or    +
  11 // + (at your option) any later version.                                  +
  12 // +                                                                      +
  13 // + enGrid is distributed in the hope that it will be useful,            +
  14 // + but WITHOUT ANY WARRANTY; without even the implied warranty of       +
  15 // + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the        +
  16 // + GNU General Public License for more details.                         +
  17 // +                                                                      +
  18 // + You should have received a copy of the GNU General Public License    +
  19 // + along with enGrid. If not, see <http://www.gnu.org/licenses/>.       +
  20 // +                                                                      +
  21 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  22 //
  23
  24 #include "checkforoverlap.h"
  25 #include "guimainwindow.h"
  26 #include "timer.h"
  27 #include "facefinder.h"
  28
  29 CheckForOverlap::CheckForOverlap()
  30 {
  31   EG_TYPENAME;
  32 }
  33
  34 void CheckForOverlap::operate()
  35 {
  36   cout << "testing mesh for overlapping faces ..." << endl;
  37   Timer timer;
  38   m_BoundaryCodes = GuiMainWindow::pointer()->getAllBoundaryCodes();
  39   QVector<bool> is_overlap(m_Grid->GetNumberOfCells(), false);
  40   QVector<vtkIdType> faces;
  41   getAllSurfaceCells(faces, m_Grid);
  42   setAllSurfaceCells();
  43   int bc_max = 0;
  44   EG_VTKDCC(vtkIntArray, cell_code, m_Grid, "cell_code");
  45
  46   FaceFinder find;
  47   find.setMaxNumFaces(100);
  48   find.setGrid(m_Grid);
  49
  50   int N_searches = 0;
  51   int N_buckets = 0;
  52   int N_faces = 0;
  53   for (int i_faces = 0; i_faces < faces.size(); ++i_faces) {
  54     vtkIdType id_face1 = faces[i_faces];
  55     bc_max = max(bc_max, cell_code->GetValue(id_face1));
  56     QVector<vec3_t> x1;
  57     vec3_t xc(0,0,0);
  58     QSet<vtkIdType> neighbours;
  59     {
  60       vtkIdType N_pts, *pts;
  61       m_Grid->GetCellPoints(id_face1, N_pts, pts);
  62       x1.resize(N_pts + 1);
  63       for (int i = 0; i < N_pts; ++i) {
  64         m_Grid->GetPoint(pts[i], x1[i].data());
  65         xc += x1[i];
  66         for (int j = 0; j < m_Part.n2cGSize(pts[i]); ++j) {
  67           neighbours.insert(m_Part.n2cGG(pts[i], j));
  68         }
  69       }
  70       xc *= 1.0/N_pts;
  71       x1[N_pts] = x1[0];
  72     }
  73     QVector<vtkIdType> close_faces;
  74     find.getCloseFaces(xc, close_faces);
  75     if (close_faces.size() == 0) {
  76       EG_BUG;
  77     }
  78     ++N_buckets;
  79     N_faces += close_faces.size();
  80     foreach (vtkIdType id_face2, close_faces) {
  81       if (!neighbours.contains(id_face2)) {
  82         ++N_searches;
  83         vtkIdType N_pts, *pts;
  84         m_Grid->GetCellPoints(id_face2, N_pts, pts);
  85         QVector<vec3_t> x2(N_pts);
  86         for (int i = 0; i < N_pts; ++i) {
  87           m_Grid->GetPoint(pts[i], x2[i].data());
  88         }
  89         for (int i = 0; i < x1.size() - 1; ++i) {
  90           vec3_t x, r;
  91           if (intersectEdgeAndTriangle(x2[0], x2[1], x2[2], x1[i], x1[i+1], x, r, 1e-6)) {
  92             is_overlap[id_face1] = true;
  93             is_overlap[id_face2] = true;
  94           }
  95         }
  96       }
  97     }
  98     if (timer()) {
  99       cout << "  " << i_faces + 1 << " of " << faces.size() << "  " << N_searches << " searches" << endl;
 100       cout << N_faces/N_buckets << " faces/bucket" << endl;
 101       N_searches = 0;
 102     }
 103   }
 104   int N = 0;
 105   foreach (vtkIdType id_face, faces) {
 106     if (is_overlap[id_face]) {
 107       cell_code->SetValue(id_face, bc_max + 1);
 108       ++N;
 109     }
 110   }
 111   if (N > 0) {
 112     GuiMainWindow::pointer()->addBC(bc_max + 1, BoundaryCondition("overlapping_faces", "patch"));
 113     GuiMainWindow::pointer()->updateBoundaryCodes(true);
 114   }
 115   cout << N << " overlapping or close faces found" << endl;
 116   cout << N_faces/N_buckets << " faces/bucket" << endl;
 117 }