FreeFOAM The Cross-Platform CFD Toolkit
refineHexMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 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
13  the Free Software Foundation, either version 3 of the License, or
14  (at your 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, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  refineHexMesh
26 
27 Description
28  Refines a hex mesh by 2x2x2 cell splitting.
29 
30 Usage
31 
32  - refineHexMesh [OPTIONS] <cellSet name>
33 
34  @param <cellSet name> \n
35  @todo Detailed description of argument.
36 
37  @param -overwrite \n
38  Overwrite existing data.
39 
40  @param -case <dir>\n
41  Case directory.
42 
43  @param -parallel \n
44  Run in parallel.
45 
46  @param -help \n
47  Display help message.
48 
49  @param -doc \n
50  Display Doxygen API documentation page for this application.
51 
52  @param -srcDoc \n
53  Display Doxygen source documentation page for this application.
54 
55 \*---------------------------------------------------------------------------*/
56 
57 #include <finiteVolume/fvMesh.H>
58 #include <OpenFOAM/pointMesh.H>
59 #include <OpenFOAM/argList.H>
60 #include <OpenFOAM/Time.H>
61 #include <dynamicMesh/hexRef8.H>
62 #include <meshTools/cellSet.H>
63 #include <OpenFOAM/OFstream.H>
64 #include <meshTools/meshTools.H>
65 #include <OpenFOAM/IFstream.H>
67 #include <OpenFOAM/mapPolyMesh.H>
68 #include <finiteVolume/volMesh.H>
70 #include <finiteVolume/volFields.H>
72 #include <OpenFOAM/pointFields.H>
73 #include <OpenFOAM/ReadFields.H>
74 
75 using namespace Foam;
76 
77 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
78 
79 // Main program:
80 int main(int argc, char *argv[])
81 {
82  argList::validOptions.insert("overwrite", "");
83  argList::validArgs.append("cellSet");
84 # include <OpenFOAM/setRootCase.H>
85 # include <OpenFOAM/createTime.H>
86  runTime.functionObjects().off();
87 # include <OpenFOAM/createMesh.H>
88  const word oldInstance = mesh.pointsInstance();
89 
90  pointMesh pMesh(mesh);
91 
92  word cellSetName(args.args()[1]);
93  bool overwrite = args.optionFound("overwrite");
94 
95  Info<< "Reading cells to refine from cellSet " << cellSetName
96  << nl << endl;
97 
98  cellSet cellsToRefine(mesh, cellSetName);
99 
100  Info<< "Read " << returnReduce(cellsToRefine.size(), sumOp<label>())
101  << " cells to refine from cellSet " << cellSetName << nl
102  << endl;
103 
104 
105  // Read objects in time directory
106  IOobjectList objects(mesh, runTime.timeName());
107 
108  // Read vol fields.
109 
111  ReadFields(mesh, objects, vsFlds);
112 
114  ReadFields(mesh, objects, vvFlds);
115 
117  ReadFields(mesh, objects, vstFlds);
118 
119  PtrList<volSymmTensorField> vsymtFlds;
120  ReadFields(mesh, objects, vsymtFlds);
121 
123  ReadFields(mesh, objects, vtFlds);
124 
125  // Read surface fields.
126 
128  ReadFields(mesh, objects, ssFlds);
129 
131  ReadFields(mesh, objects, svFlds);
132 
134  ReadFields(mesh, objects, sstFlds);
135 
137  ReadFields(mesh, objects, ssymtFlds);
138 
140  ReadFields(mesh, objects, stFlds);
141 
142  // Read point fields
144  ReadFields(pMesh, objects, psFlds);
145 
147  ReadFields(pMesh, objects, pvFlds);
148 
149 
150 
151  // Construct refiner without unrefinement. Read existing point/cell level.
153 
154  // Some stats
155  Info<< "Read mesh:" << nl
156  << " cells:" << mesh.globalData().nTotalCells() << nl
157  << " faces:" << mesh.globalData().nTotalFaces() << nl
158  << " points:" << mesh.globalData().nTotalPoints() << nl
159  << " cellLevel :"
160  << " min:" << gMin(meshCutter.cellLevel())
161  << " max:" << gMax(meshCutter.cellLevel()) << nl
162  << " pointLevel :"
163  << " min:" << gMin(meshCutter.pointLevel())
164  << " max:" << gMax(meshCutter.pointLevel()) << nl
165  << endl;
166 
167 
168  // Maintain 2:1 ratio
169  labelList newCellsToRefine
170  (
171  meshCutter.consistentRefinement
172  (
173  cellsToRefine.toc(),
174  true // extend set
175  )
176  );
177 
178  // Mesh changing engine.
179  polyTopoChange meshMod(mesh);
180 
181  // Play refinement commands into mesh changer.
182  meshCutter.setRefinement(newCellsToRefine, meshMod);
183 
184  if (!overwrite)
185  {
186  runTime++;
187  }
188 
189  // Create mesh, return map from old to new mesh.
190  autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
191 
192  // Update fields
193  mesh.updateMesh(map);
194  pMesh.updateMesh(map);
195 
196  // Update numbering of cells/vertices.
197  meshCutter.updateMesh(map);
198 
199  // Optionally inflate mesh
200  if (map().hasMotionPoints())
201  {
202  mesh.movePoints(map().preMotionPoints());
203  pMesh.movePoints(map().preMotionPoints());
204  }
205 
206  Pout<< "Refined from " << returnReduce(map().nOldCells(), sumOp<label>())
207  << " to " << mesh.globalData().nTotalCells() << " cells." << nl << endl;
208 
209  if (overwrite)
210  {
211  mesh.setInstance(oldInstance);
212  meshCutter.setInstance(oldInstance);
213  }
214  Info<< "Writing mesh to " << runTime.timeName() << endl;
215 
216  mesh.write();
217  meshCutter.write();
218 
219  Info<< "End\n" << endl;
220 
221  return 0;
222 }
223 
224 
225 // ************************ vim: set sw=4 sts=4 et: ************************ //