2 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
4 // + This file is part of enGrid. +
6 // + Copyright 2008,2009 Oliver Gloth +
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. +
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. +
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/>. +
21 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
23 #include "guicreateboundarylayer.h"
24 #include "seedsimpleprismaticlayer.h"
25 #include "gridsmoother.h"
26 #include "createvolumemesh.h"
27 #include "swaptriangles.h"
28 #include "deletetetras.h"
29 #include "deletecells.h"
31 GuiCreateBoundaryLayer::GuiCreateBoundaryLayer()
33 getSet("boundary layer", "maximal relative error", 0.01, err_max
);
34 getSet("boundary layer", "maximal number of smoothing iterations", 10, max_iter
);
37 void GuiCreateBoundaryLayer::before()
39 populateBoundaryCodes(ui
.listWidget
);
42 void GuiCreateBoundaryLayer::operate()
44 getSelectedItems(ui
.listWidget
, boundary_codes
);
45 getSurfaceCells(boundary_codes
, layer_cells
, grid
);
47 cout
<< "\n\ncreating boundary layer mesh)" << endl
;
50 EG_VTKDCN(vtkIntArray
, node_status
, grid
, "node_status");
51 EG_VTKDCN(vtkIntArray
, node_layer
, grid
, "node_layer");
52 EG_VTKDCC(vtkIntArray
, bc
, grid
, "cell_code");
53 foreach(vtkIdType id_node
, nodes
) {
54 node_status
->SetValue(id_node
, 0);
56 foreach (int i_neigh_cells
, n2c
[_nodes
[id_node
]]) {
57 vtkIdType id_neigh_cell
= cells
[i_neigh_cells
];
58 if (isSurface(grid
, id_neigh_cell
)) {
59 if (boundary_codes
.contains(bc
->GetValue(id_neigh_cell
))) {
60 bcs
.insert(bc
->GetValue(id_neigh_cell
));
64 if (bcs
.size() >= 2) {
65 node_status
->SetValue(id_node
, 1);
67 node_layer
->SetValue(id_node
, -1);
69 foreach (vtkIdType id_cell
, layer_cells
) {
70 vtkIdType N_pts
, *pts
;
71 grid
->GetCellPoints(id_cell
, N_pts
, pts
);
72 for (int i_pts
= 0; i_pts
< N_pts
; ++i_pts
) {
73 node_layer
->SetValue(pts
[i_pts
], 0);
78 cout
<< "preparing prismatic layer" << endl
;
82 smooth
.setBoundaryCodes(boundary_codes
);
84 //smooth.setNumIterations(5);
86 SeedSimplePrismaticLayer seed_layer
;
92 swap
.setBoundaryCodes(boundary_codes
);
96 seed_layer
.setGrid(grid
);
99 seed_layer
.setAllCells();
100 seed_layer
.setLayerCells(layer_cells
);
101 seed_layer
.setBoundaryCodes(boundary_codes
);
103 seed_layer
.getLayerCells(layer_cells
);
105 for (int j
= 0; j
< max_iter
; ++j
) {
106 cout
<< "improving prismatic layer -> iteration " << j
+1 << "/" << max_iter
<< endl
;
107 smooth
.setAllCells();
112 vol
.setTraceCells(layer_cells
);
114 vol
.getTraceCells(layer_cells
);
115 if (smooth
.improvement() < err_max
) break;
117 //smooth.setAllCells();
123 cout << "finalising prismatic layer" << endl;
124 smooth.setAllCells();