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) 2011-2016 OpenFOAM Foundation
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/topoSet utilities to see how to select cells based on
32  various shapes.
33 
34  The mesh will subset all points, faces and cells needed to make a sub-mesh
35  but will not preserve attached boundary types.
36 
37 \*---------------------------------------------------------------------------*/
38 
39 #include "fvMeshSubset.H"
40 #include "argList.H"
41 #include "cellSet.H"
42 #include "IOobjectList.H"
43 #include "volFields.H"
44 
45 using namespace Foam;
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 
50 template<class Type>
51 void subsetVolFields
52 (
53  const fvMeshSubset& subsetter,
54  const wordList& fieldNames,
56 )
57 {
58  const fvMesh& baseMesh = subsetter.baseMesh();
59 
60  forAll(fieldNames, i)
61  {
62  const word& fieldName = fieldNames[i];
63 
64  Info<< "Subsetting field " << fieldName << endl;
65 
67  (
68  IOobject
69  (
70  fieldName,
71  baseMesh.time().timeName(),
72  baseMesh,
75  ),
76  baseMesh
77  );
78 
79  subFields.set(i, subsetter.interpolate(fld));
80  }
81 }
82 
83 
84 template<class Type>
85 void subsetSurfaceFields
86 (
87  const fvMeshSubset& subsetter,
88  const wordList& fieldNames,
90 )
91 {
92  const fvMesh& baseMesh = subsetter.baseMesh();
93 
94  forAll(fieldNames, i)
95  {
96  const word& fieldName = fieldNames[i];
97 
98  Info<< "Subsetting field " << fieldName << endl;
99 
101  (
102  IOobject
103  (
104  fieldName,
105  baseMesh.time().timeName(),
106  baseMesh,
109  ),
110  baseMesh
111  );
112 
113  subFields.set(i, subsetter.interpolate(fld));
114  }
115 }
116 
117 
118 template<class Type>
119 void subsetPointFields
120 (
121  const fvMeshSubset& subsetter,
122  const pointMesh& pMesh,
123  const wordList& fieldNames,
125 )
126 {
127  const fvMesh& baseMesh = subsetter.baseMesh();
128 
129  forAll(fieldNames, i)
130  {
131  const word& fieldName = fieldNames[i];
132 
133  Info<< "Subsetting field " << fieldName << endl;
134 
136  (
137  IOobject
138  (
139  fieldName,
140  baseMesh.time().timeName(),
141  baseMesh,
144  ),
145  pMesh
146  );
147 
148  subFields.set(i, subsetter.interpolate(fld));
149  }
150 }
151 
152 
153 template<class Type>
154 void subsetDimensionedFields
155 (
156  const fvMeshSubset& subsetter,
157  const wordList& fieldNames,
159 )
160 {
161  const fvMesh& baseMesh = subsetter.baseMesh();
162 
163  forAll(fieldNames, i)
164  {
165  const word& fieldName = fieldNames[i];
166 
167  Info<< "Subsetting field " << fieldName << endl;
168 
170  (
171  IOobject
172  (
173  fieldName,
174  baseMesh.time().timeName(),
175  baseMesh,
178  ),
179  baseMesh
180  );
181 
182  subFields.set(i, subsetter.interpolate(fld));
183  }
184 }
185 
186 
187 
188 int main(int argc, char *argv[])
189 {
191  (
192  "select a mesh subset based on a cellSet"
193  );
194 
195  #include "addOverwriteOption.H"
196  #include "addRegionOption.H"
197  argList::validArgs.append("cellSet");
199  (
200  "patch",
201  "name",
202  "add exposed internal faces to specified patch instead of to "
203  "'oldInternalFaces'"
204  );
206  (
207  "resultTime",
208  "time",
209  "specify a time for the resulting mesh"
210  );
211  #include "setRootCase.H"
212  #include "createTime.H"
213  runTime.functionObjects().off();
214 
215  Foam::word meshRegionName = polyMesh::defaultRegion;
216  args.optionReadIfPresent("region", meshRegionName);
217 
218  #include "createNamedMesh.H"
219 
220 
221  const word setName = args[1];
222 
223  word meshInstance = mesh.pointsInstance();
224  word fieldsInstance = runTime.timeName();
225 
226  const bool overwrite = args.optionFound("overwrite");
227  const bool specifiedInstance = args.optionReadIfPresent
228  (
229  "resultTime",
230  fieldsInstance
231  );
232  if (specifiedInstance)
233  {
234  // Set both mesh and field to this time
235  meshInstance = fieldsInstance;
236  }
237 
238 
239  Info<< "Reading cell set from " << setName << endl << endl;
240 
241  // Create mesh subsetting engine
242  fvMeshSubset subsetter(mesh);
243 
244  label patchi = -1;
245 
246  if (args.optionFound("patch"))
247  {
248  const word patchName = args["patch"];
249 
250  patchi = mesh.boundaryMesh().findPatchID(patchName);
251 
252  if (patchi == -1)
253  {
255  << nl << "Valid patches are " << mesh.boundaryMesh().names()
256  << exit(FatalError);
257  }
258 
259  Info<< "Adding exposed internal faces to patch " << patchName << endl
260  << endl;
261  }
262  else
263  {
264  Info<< "Adding exposed internal faces to a patch called"
265  << " \"oldInternalFaces\" (created if necessary)" << endl
266  << endl;
267  }
268 
269 
270  cellSet currentSet(mesh, setName);
271 
272  subsetter.setLargeCellSubset(currentSet, patchi, true);
273 
274  IOobjectList objects(mesh, runTime.timeName());
275 
276 
277  // Read vol fields and subset
278  // ~~~~~~~~~~~~~~~~~~~~~~~~~~
279 
280  wordList scalarNames(objects.names(volScalarField::typeName));
281  PtrList<volScalarField> scalarFlds(scalarNames.size());
282  subsetVolFields(subsetter, scalarNames, scalarFlds);
283 
284  wordList vectorNames(objects.names(volVectorField::typeName));
285  PtrList<volVectorField> vectorFlds(vectorNames.size());
286  subsetVolFields(subsetter, vectorNames, vectorFlds);
287 
288  wordList sphericalTensorNames
289  (
290  objects.names(volSphericalTensorField::typeName)
291  );
292  PtrList<volSphericalTensorField> sphericalTensorFlds
293  (
294  sphericalTensorNames.size()
295  );
296  subsetVolFields(subsetter, sphericalTensorNames, sphericalTensorFlds);
297 
298  wordList symmTensorNames(objects.names(volSymmTensorField::typeName));
299  PtrList<volSymmTensorField> symmTensorFlds(symmTensorNames.size());
300  subsetVolFields(subsetter, symmTensorNames, symmTensorFlds);
301 
302  wordList tensorNames(objects.names(volTensorField::typeName));
303  PtrList<volTensorField> tensorFlds(tensorNames.size());
304  subsetVolFields(subsetter, tensorNames, tensorFlds);
305 
306 
307  // Read surface fields and subset
308  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
309 
310  wordList surfScalarNames(objects.names(surfaceScalarField::typeName));
311  PtrList<surfaceScalarField> surfScalarFlds(surfScalarNames.size());
312  subsetSurfaceFields(subsetter, surfScalarNames, surfScalarFlds);
313 
314  wordList surfVectorNames(objects.names(surfaceVectorField::typeName));
315  PtrList<surfaceVectorField> surfVectorFlds(surfVectorNames.size());
316  subsetSurfaceFields(subsetter, surfVectorNames, surfVectorFlds);
317 
318  wordList surfSphericalTensorNames
319  (
321  );
322  PtrList<surfaceSphericalTensorField> surfSphericalTensorFlds
323  (
324  surfSphericalTensorNames.size()
325  );
326  subsetSurfaceFields
327  (
328  subsetter,
329  surfSphericalTensorNames,
330  surfSphericalTensorFlds
331  );
332 
333  wordList surfSymmTensorNames
334  (
335  objects.names(surfaceSymmTensorField::typeName)
336  );
337  PtrList<surfaceSymmTensorField> surfSymmTensorFlds
338  (
339  surfSymmTensorNames.size()
340  );
341  subsetSurfaceFields(subsetter, surfSymmTensorNames, surfSymmTensorFlds);
342 
343  wordList surfTensorNames(objects.names(surfaceTensorField::typeName));
344  PtrList<surfaceTensorField> surfTensorFlds(surfTensorNames.size());
345  subsetSurfaceFields(subsetter, surfTensorNames, surfTensorFlds);
346 
347 
348  // Read point fields and subset
349  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
350 
351  const pointMesh& pMesh = pointMesh::New(mesh);
352 
353  wordList pointScalarNames(objects.names(pointScalarField::typeName));
354  PtrList<pointScalarField> pointScalarFlds(pointScalarNames.size());
355  subsetPointFields(subsetter, pMesh, pointScalarNames, pointScalarFlds);
356 
357  wordList pointVectorNames(objects.names(pointVectorField::typeName));
358  PtrList<pointVectorField> pointVectorFlds(pointVectorNames.size());
359  subsetPointFields(subsetter, pMesh, pointVectorNames, pointVectorFlds);
360 
361  wordList pointSphericalTensorNames
362  (
364  );
365  PtrList<pointSphericalTensorField> pointSphericalTensorFlds
366  (
367  pointSphericalTensorNames.size()
368  );
369  subsetPointFields
370  (
371  subsetter,
372  pMesh,
373  pointSphericalTensorNames,
374  pointSphericalTensorFlds
375  );
376 
377  wordList pointSymmTensorNames
378  (
379  objects.names(pointSymmTensorField::typeName)
380  );
381  PtrList<pointSymmTensorField> pointSymmTensorFlds
382  (
383  pointSymmTensorNames.size()
384  );
385  subsetPointFields
386  (
387  subsetter,
388  pMesh,
389  pointSymmTensorNames,
390  pointSymmTensorFlds
391  );
392 
393  wordList pointTensorNames(objects.names(pointTensorField::typeName));
394  PtrList<pointTensorField> pointTensorFlds(pointTensorNames.size());
395  subsetPointFields(subsetter, pMesh, pointTensorNames, pointTensorFlds);
396 
397 
398  // Read dimensioned fields and subset
399  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
400 
401  typedef volScalarField::Internal dimScalType;
402  wordList scalarDimNames(objects.names(dimScalType::typeName));
403  PtrList<dimScalType> scalarDimFlds(scalarDimNames.size());
404  subsetDimensionedFields(subsetter, scalarDimNames, scalarDimFlds);
405 
406  typedef volVectorField::Internal dimVecType;
407  wordList vectorDimNames(objects.names(dimVecType::typeName));
408  PtrList<dimVecType> vectorDimFlds(vectorDimNames.size());
409  subsetDimensionedFields(subsetter, vectorDimNames, vectorDimFlds);
410 
411  typedef volSphericalTensorField::Internal dimSphereType;
412  wordList sphericalTensorDimNames(objects.names(dimSphereType::typeName));
413  PtrList<dimSphereType> sphericalTensorDimFlds
414  (
415  sphericalTensorDimNames.size()
416  );
417  subsetDimensionedFields
418  (
419  subsetter,
420  sphericalTensorDimNames,
421  sphericalTensorDimFlds
422  );
423 
424  typedef volSymmTensorField::Internal dimSymmTensorType;
425  wordList symmTensorDimNames(objects.names(dimSymmTensorType::typeName));
426  PtrList<dimSymmTensorType> symmTensorDimFlds(symmTensorDimNames.size());
427  subsetDimensionedFields(subsetter, symmTensorDimNames, symmTensorDimFlds);
428 
429  typedef volTensorField::Internal dimTensorType;
430  wordList tensorDimNames(objects.names(dimTensorType::typeName));
431  PtrList<dimTensorType> tensorDimFlds(tensorDimNames.size());
432  subsetDimensionedFields(subsetter, tensorDimNames, tensorDimFlds);
433 
434 
435  // Write mesh and fields to new time
436  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
437 
438  if (overwrite || specifiedInstance)
439  {
440  runTime.setTime(instant(fieldsInstance), 0);
441  subsetter.subMesh().setInstance(meshInstance);
442  }
443  else
444  {
445  runTime++;
446  }
447 
448  Info<< "Writing subsetted mesh and fields to time " << runTime.timeName()
449  << endl;
450  subsetter.subMesh().write();
451 
452 
453  // Subsetting adds 'subset' prefix. Rename field to be like original.
454  forAll(scalarFlds, i)
455  {
456  scalarFlds[i].rename(scalarNames[i]);
457  scalarFlds[i].write();
458  }
459  forAll(vectorFlds, i)
460  {
461  vectorFlds[i].rename(vectorNames[i]);
462  vectorFlds[i].write();
463  }
464  forAll(sphericalTensorFlds, i)
465  {
466  sphericalTensorFlds[i].rename(sphericalTensorNames[i]);
467  sphericalTensorFlds[i].write();
468  }
469  forAll(symmTensorFlds, i)
470  {
471  symmTensorFlds[i].rename(symmTensorNames[i]);
472  symmTensorFlds[i].write();
473  }
474  forAll(tensorFlds, i)
475  {
476  tensorFlds[i].rename(tensorNames[i]);
477  tensorFlds[i].write();
478  }
479 
480  // Surface ones.
481  forAll(surfScalarFlds, i)
482  {
483  surfScalarFlds[i].rename(surfScalarNames[i]);
484  surfScalarFlds[i].write();
485  }
486  forAll(surfVectorFlds, i)
487  {
488  surfVectorFlds[i].rename(surfVectorNames[i]);
489  surfVectorFlds[i].write();
490  }
491  forAll(surfSphericalTensorFlds, i)
492  {
493  surfSphericalTensorFlds[i].rename(surfSphericalTensorNames[i]);
494  surfSphericalTensorFlds[i].write();
495  }
496  forAll(surfSymmTensorFlds, i)
497  {
498  surfSymmTensorFlds[i].rename(surfSymmTensorNames[i]);
499  surfSymmTensorFlds[i].write();
500  }
501  forAll(surfTensorNames, i)
502  {
503  surfTensorFlds[i].rename(surfTensorNames[i]);
504  surfTensorFlds[i].write();
505  }
506 
507  // Point ones
508  forAll(pointScalarFlds, i)
509  {
510  pointScalarFlds[i].rename(pointScalarNames[i]);
511  pointScalarFlds[i].write();
512  }
513  forAll(pointVectorFlds, i)
514  {
515  pointVectorFlds[i].rename(pointVectorNames[i]);
516  pointVectorFlds[i].write();
517  }
518  forAll(pointSphericalTensorFlds, i)
519  {
520  pointSphericalTensorFlds[i].rename(pointSphericalTensorNames[i]);
521  pointSphericalTensorFlds[i].write();
522  }
523  forAll(pointSymmTensorFlds, i)
524  {
525  pointSymmTensorFlds[i].rename(pointSymmTensorNames[i]);
526  pointSymmTensorFlds[i].write();
527  }
528  forAll(pointTensorNames, i)
529  {
530  pointTensorFlds[i].rename(pointTensorNames[i]);
531  pointTensorFlds[i].write();
532  }
533 
534  // DimensionedFields
535  forAll(scalarDimFlds, i)
536  {
537  scalarDimFlds[i].rename(scalarDimNames[i]);
538  scalarDimFlds[i].write();
539  }
540  forAll(vectorDimFlds, i)
541  {
542  vectorDimFlds[i].rename(vectorDimNames[i]);
543  vectorDimFlds[i].write();
544  }
545  forAll(sphericalTensorDimFlds, i)
546  {
547  sphericalTensorDimFlds[i].rename(sphericalTensorDimNames[i]);
548  sphericalTensorDimFlds[i].write();
549  }
550  forAll(symmTensorDimFlds, i)
551  {
552  symmTensorDimFlds[i].rename(symmTensorDimNames[i]);
553  symmTensorDimFlds[i].write();
554  }
555  forAll(tensorDimFlds, i)
556  {
557  tensorDimFlds[i].rename(tensorDimNames[i]);
558  tensorDimFlds[i].write();
559  }
560 
561 
562  Info<< "End\n" << endl;
563 
564  return 0;
565 }
566 
567 
568 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
wordList names() const
Return a list of patch names.
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static const char *const typeName
Definition: Field.H:94
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
Definition: argListI.H:198
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
static SLList< string > validArgs
A list of valid (mandatory) arguments.
Definition: argList.H:154
Generic GeometricField class.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
const fvMesh & subMesh() const
Return reference to subset mesh.
static const pointMesh & New(const polyMesh &mesh)
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
dynamicFvMesh & mesh
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
A class for handling words, derived from string.
Definition: word.H:59
static void addOption(const word &opt, const string &param="", const string &usage="")
Add to an option to validOptions with usage information.
Definition: argList.C:93
static const char nl
Definition: Ostream.H:262
Post-processing mesh subset tool. Given the original mesh and the list of selected cells...
Definition: fvMeshSubset.H:73
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
An instant of time. Contains the time value and name.
Definition: instant.H:64
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:772
A collection of cell labels.
Definition: cellSet.H:48
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const fvMesh & baseMesh() const
Original mesh.
Definition: fvMeshSubset.H:217
messageStream Info
static void addNote(const string &)
Add extra notes for the usage information.
Definition: argList.C:124
Foam::argList args(argc, argv)
label findPatchID(const word &patchName) const
Find patch index given a name.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
static tmp< GeometricField< Type, fvPatchField, volMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &, const fvMesh &sMesh, const labelList &patchMap, const labelList &cellMap, const labelList &faceMap)
Map volume field.
void setLargeCellSubset(const labelList &region, const label currentRegion, const label patchID=-1, const bool syncCouples=true)
Set the subset from all cells with region == currentRegion.
Definition: fvMeshSubset.C:795
virtual bool write() const
Write mesh using IO settings from time.
Definition: fvMesh.C:870
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243