FreeFOAM The Cross-Platform CFD Toolkit
subsetMesh.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  subsetMesh
26 
27 Description
28  Selects a section of mesh based on a cellSet.
29 
30  The utility sub-sets the mesh to choose only a part of interest. Check
31  the setSet/cellSet utilities to see how to select cells based on various.
32 
33  The mesh will subset all points, faces and cells needed to make a sub-mesh
34  but will not preserve attached boundary types.
35 
36 Usage
37 
38  - subsetMesh [OPTIONS] <cellSet name> <patchName>
39 
40  @param <cellSet name> \n
41  @todo Detailed description of argument.
42 
43  @param <patchName> \n
44  @todo Detailed description of argument.
45 
46  @param -overwrite \n
47  Overwrite existing data.
48 
49  @param -case <dir>\n
50  Case directory.
51 
52  @param -parallel \n
53  Run in parallel.
54 
55  @param -help \n
56  Display help message.
57 
58  @param -doc \n
59  Display Doxygen API documentation page for this application.
60 
61  @param -srcDoc \n
62  Display Doxygen source documentation page for this application.
63 
64 \*---------------------------------------------------------------------------*/
65 
67 #include <OpenFOAM/argList.H>
68 #include <meshTools/cellSet.H>
69 #include <OpenFOAM/IOobjectList.H>
70 #include <finiteVolume/volFields.H>
71 
72 using namespace Foam;
73 
74 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
75 
76 
77 template<class Type>
78 void subsetVolFields
79 (
80  const fvMeshSubset& subsetter,
81  const wordList& fieldNames,
83 )
84 {
85  const fvMesh& baseMesh = subsetter.baseMesh();
86 
87  forAll(fieldNames, i)
88  {
89  const word& fieldName = fieldNames[i];
90 
91  Info<< "Subsetting field " << fieldName << endl;
92 
94  (
95  IOobject
96  (
97  fieldName,
98  baseMesh.time().timeName(),
99  baseMesh,
102  ),
103  baseMesh
104  );
105 
106  subFields.set(i, subsetter.interpolate(fld));
107  }
108 }
109 
110 
111 template<class Type>
112 void subsetSurfaceFields
113 (
114  const fvMeshSubset& subsetter,
115  const wordList& fieldNames,
117 )
118 {
119  const fvMesh& baseMesh = subsetter.baseMesh();
120 
121  forAll(fieldNames, i)
122  {
123  const word& fieldName = fieldNames[i];
124 
125  Info<< "Subsetting field " << fieldName << endl;
126 
128  (
129  IOobject
130  (
131  fieldName,
132  baseMesh.time().timeName(),
133  baseMesh,
136  ),
137  baseMesh
138  );
139 
140  subFields.set(i, subsetter.interpolate(fld));
141  }
142 }
143 
144 
145 template<class Type>
146 void subsetPointFields
147 (
148  const fvMeshSubset& subsetter,
149  const pointMesh& pMesh,
150  const wordList& fieldNames,
152 )
153 {
154  const fvMesh& baseMesh = subsetter.baseMesh();
155 
156  forAll(fieldNames, i)
157  {
158  const word& fieldName = fieldNames[i];
159 
160  Info<< "Subsetting field " << fieldName << endl;
161 
163  (
164  IOobject
165  (
166  fieldName,
167  baseMesh.time().timeName(),
168  baseMesh,
171  ),
172  pMesh
173  );
174 
175  subFields.set(i, subsetter.interpolate(fld));
176  }
177 }
178 
179 
180 // Main program:
181 
182 int main(int argc, char *argv[])
183 {
184  argList::validArgs.append("set");
185  argList::validOptions.insert("patch", "patch name");
186  argList::validOptions.insert("overwrite", "");
187 
188 # include <OpenFOAM/setRootCase.H>
189 # include <OpenFOAM/createTime.H>
190  runTime.functionObjects().off();
191 # include <OpenFOAM/createMesh.H>
192  const word oldInstance = mesh.pointsInstance();
193 
194  word setName(args.additionalArgs()[0]);
195  bool overwrite = args.optionFound("overwrite");
196 
197 
198  Info<< "Reading cell set from " << setName << endl << endl;
199 
200  // Create mesh subsetting engine
201  fvMeshSubset subsetter(mesh);
202 
203  label patchI = -1;
204 
205  if (args.optionFound("patch"))
206  {
207  word patchName(args.option("patch"));
208 
209  patchI = mesh.boundaryMesh().findPatchID(patchName);
210 
211  if (patchI == -1)
212  {
213  FatalErrorIn(args.executable()) << "Illegal patch " << patchName
214  << nl << "Valid patches are " << mesh.boundaryMesh().names()
215  << exit(FatalError);
216  }
217 
218  Info<< "Adding exposed internal faces to patch " << patchName << endl
219  << endl;
220  }
221  else
222  {
223  Info<< "Adding exposed internal faces to a patch called"
224  << " \"oldInternalFaces\" (created if nessecary)" << endl
225  << endl;
226  }
227 
228 
229  cellSet currentSet(mesh, setName);
230 
231  subsetter.setLargeCellSubset(currentSet, patchI, true);
232 
233  IOobjectList objects(mesh, runTime.timeName());
234 
235 
236  // Read vol fields and subset
237  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
238 
239  wordList scalarNames(objects.names(volScalarField::typeName));
240  PtrList<volScalarField> scalarFlds(scalarNames.size());
241  subsetVolFields(subsetter, scalarNames, scalarFlds);
242 
243  wordList vectorNames(objects.names(volVectorField::typeName));
244  PtrList<volVectorField> vectorFlds(vectorNames.size());
245  subsetVolFields(subsetter, vectorNames, vectorFlds);
246 
247  wordList sphericalTensorNames
248  (
249  objects.names(volSphericalTensorField::typeName)
250  );
251  PtrList<volSphericalTensorField> sphericalTensorFlds
252  (
253  sphericalTensorNames.size()
254  );
255  subsetVolFields(subsetter, sphericalTensorNames, sphericalTensorFlds);
256 
257  wordList symmTensorNames(objects.names(volSymmTensorField::typeName));
258  PtrList<volSymmTensorField> symmTensorFlds(symmTensorNames.size());
259  subsetVolFields(subsetter, symmTensorNames, symmTensorFlds);
260 
261  wordList tensorNames(objects.names(volTensorField::typeName));
262  PtrList<volTensorField> tensorFlds(tensorNames.size());
263  subsetVolFields(subsetter, tensorNames, tensorFlds);
264 
265 
266  // Read surface fields and subset
267  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
268 
271  subsetSurfaceFields(subsetter, surfScalarNames, surfScalarFlds);
272 
275  subsetSurfaceFields(subsetter, surfVectorNames, surfVectorFlds);
276 
277  wordList surfSphericalTensorNames
278  (
280  );
281  PtrList<surfaceSphericalTensorField> surfSphericalTensorFlds
282  (
283  surfSphericalTensorNames.size()
284  );
285  subsetSurfaceFields
286  (
287  subsetter,
288  surfSphericalTensorNames,
289  surfSphericalTensorFlds
290  );
291 
292  wordList surfSymmTensorNames
293  (
294  objects.names(surfaceSymmTensorField::typeName)
295  );
296  PtrList<surfaceSymmTensorField> surfSymmTensorFlds
297  (
298  surfSymmTensorNames.size()
299  );
300  subsetSurfaceFields(subsetter, surfSymmTensorNames, surfSymmTensorFlds);
301 
302  wordList surfTensorNames(objects.names(surfaceTensorField::typeName));
303  PtrList<surfaceTensorField> surfTensorFlds(surfTensorNames.size());
304  subsetSurfaceFields(subsetter, surfTensorNames, surfTensorFlds);
305 
306 
307  // Read point fields and subset
308  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
309 
310  const pointMesh& pMesh = pointMesh::New(mesh);
311 
312  wordList pointScalarNames(objects.names(pointScalarField::typeName));
313  PtrList<pointScalarField> pointScalarFlds(pointScalarNames.size());
314  subsetPointFields(subsetter, pMesh, pointScalarNames, pointScalarFlds);
315 
316  wordList pointVectorNames(objects.names(pointVectorField::typeName));
317  PtrList<pointVectorField> pointVectorFlds(pointVectorNames.size());
318  subsetPointFields(subsetter, pMesh, pointVectorNames, pointVectorFlds);
319 
320  wordList pointSphericalTensorNames
321  (
323  );
324  PtrList<pointSphericalTensorField> pointSphericalTensorFlds
325  (
326  pointSphericalTensorNames.size()
327  );
328  subsetPointFields
329  (
330  subsetter,
331  pMesh,
332  pointSphericalTensorNames,
333  pointSphericalTensorFlds
334  );
335 
336  wordList pointSymmTensorNames
337  (
338  objects.names(pointSymmTensorField::typeName)
339  );
340  PtrList<pointSymmTensorField> pointSymmTensorFlds
341  (
342  pointSymmTensorNames.size()
343  );
344  subsetPointFields
345  (
346  subsetter,
347  pMesh,
348  pointSymmTensorNames,
349  pointSymmTensorFlds
350  );
351 
352  wordList pointTensorNames(objects.names(pointTensorField::typeName));
353  PtrList<pointTensorField> pointTensorFlds(pointTensorNames.size());
354  subsetPointFields(subsetter, pMesh, pointTensorNames, pointTensorFlds);
355 
356 
357 
358  // Write mesh and fields to new time
359  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
360 
361  if (!overwrite)
362  {
363  runTime++;
364  }
365  else
366  {
367  subsetter.subMesh().setInstance(oldInstance);
368  }
369 
370  Info<< "Writing subsetted mesh and fields to time " << runTime.timeName()
371  << endl;
372  subsetter.subMesh().write();
373 
374 
375  // Subsetting adds 'subset' prefix. Rename field to be like original.
376  forAll(scalarFlds, i)
377  {
378  scalarFlds[i].rename(scalarNames[i]);
379 
380  scalarFlds[i].write();
381  }
382  forAll(vectorFlds, i)
383  {
384  vectorFlds[i].rename(vectorNames[i]);
385 
386  vectorFlds[i].write();
387  }
388  forAll(sphericalTensorFlds, i)
389  {
390  sphericalTensorFlds[i].rename(sphericalTensorNames[i]);
391 
392  sphericalTensorFlds[i].write();
393  }
394  forAll(symmTensorFlds, i)
395  {
396  symmTensorFlds[i].rename(symmTensorNames[i]);
397 
398  symmTensorFlds[i].write();
399  }
400  forAll(tensorFlds, i)
401  {
402  tensorFlds[i].rename(tensorNames[i]);
403 
404  tensorFlds[i].write();
405  }
406 
407  // Surface ones.
408  forAll(surfScalarFlds, i)
409  {
410  surfScalarFlds[i].rename(surfScalarNames[i]);
411 
412  surfScalarFlds[i].write();
413  }
414  forAll(surfVectorFlds, i)
415  {
416  surfVectorFlds[i].rename(surfVectorNames[i]);
417 
418  surfVectorFlds[i].write();
419  }
420  forAll(surfSphericalTensorFlds, i)
421  {
422  surfSphericalTensorFlds[i].rename(surfSphericalTensorNames[i]);
423 
424  surfSphericalTensorFlds[i].write();
425  }
426  forAll(surfSymmTensorFlds, i)
427  {
428  surfSymmTensorFlds[i].rename(surfSymmTensorNames[i]);
429 
430  surfSymmTensorFlds[i].write();
431  }
432  forAll(surfTensorNames, i)
433  {
434  surfTensorFlds[i].rename(surfTensorNames[i]);
435 
436  surfTensorFlds[i].write();
437  }
438 
439  // Point ones
440  forAll(pointScalarFlds, i)
441  {
442  pointScalarFlds[i].rename(pointScalarNames[i]);
443 
444  pointScalarFlds[i].write();
445  }
446  forAll(pointVectorFlds, i)
447  {
448  pointVectorFlds[i].rename(pointVectorNames[i]);
449 
450  pointVectorFlds[i].write();
451  }
452  forAll(pointSphericalTensorFlds, i)
453  {
454  pointSphericalTensorFlds[i].rename(pointSphericalTensorNames[i]);
455 
456  pointSphericalTensorFlds[i].write();
457  }
458  forAll(pointSymmTensorFlds, i)
459  {
460  pointSymmTensorFlds[i].rename(pointSymmTensorNames[i]);
461 
462  pointSymmTensorFlds[i].write();
463  }
464  forAll(pointTensorNames, i)
465  {
466  pointTensorFlds[i].rename(pointTensorNames[i]);
467 
468  pointTensorFlds[i].write();
469  }
470 
471 
472  Info << nl << "End" << endl;
473 
474  return 0;
475 }
476 
477 
478 // ************************ vim: set sw=4 sts=4 et: ************************ //