patchSummary.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  patchSummary
26 
27 Description
28  Writes fields and boundary condition info for each patch at each requested
29  time instance.
30 
31  Default action is to write a single entry for patches/patchGroups with the
32  same boundary conditions. Use the -expand option to print every patch
33  separately. In case of multiple groups matching it will print only the
34  first one.
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #include "fvCFD.H"
39 #include "volFields.H"
40 #include "pointFields.H"
41 #include "IOobjectList.H"
42 #include "patchSummaryTemplates.H"
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 int main(int argc, char *argv[])
47 {
48  timeSelector::addOptions();
49 
50  #include "addRegionOption.H"
51  argList::addBoolOption
52  (
53  "expand",
54  "Do not combine patches"
55  );
56  #include "setRootCase.H"
57  #include "createTime.H"
58 
59  instantList timeDirs = timeSelector::select0(runTime, args);
60 
61  const bool expand = args.optionFound("expand");
62 
63 
64  #include "createNamedMesh.H"
65  const polyBoundaryMesh& bm = mesh.boundaryMesh();
66 
67 
68  forAll(timeDirs, timeI)
69  {
70  runTime.setTime(timeDirs[timeI], timeI);
71 
72  Info<< "Time = " << runTime.timeName() << nl << endl;
73 
74  // Update the mesh if changed
75  if (mesh.readUpdate() == polyMesh::TOPO_PATCH_CHANGE)
76  {
77  Info<< "Detected changed patches. Recreating patch group table."
78  << endl;
79  }
80 
81 
82  const IOobjectList fieldObjs(mesh, runTime.timeName());
83  const wordList objNames = fieldObjs.names();
84 
85  PtrList<volScalarField> vsf(objNames.size());
86  PtrList<volVectorField> vvf(objNames.size());
87  PtrList<volSphericalTensorField> vsptf(objNames.size());
88  PtrList<volSymmTensorField> vsytf(objNames.size());
89  PtrList<volTensorField> vtf(objNames.size());
90 
91  PtrList<pointScalarField> psf(objNames.size());
92  PtrList<pointVectorField> pvf(objNames.size());
93  PtrList<pointSphericalTensorField> psptf(objNames.size());
94  PtrList<pointSymmTensorField> psytf(objNames.size());
95  PtrList<pointTensorField> ptf(objNames.size());
96 
97  Info<< "Valid fields:" << endl;
98 
99  forAll(objNames, objI)
100  {
101  IOobject obj
102  (
103  objNames[objI],
104  runTime.timeName(),
105  mesh,
106  IOobject::MUST_READ
107  );
108 
109  if (obj.headerOk())
110  {
111  addToFieldList(vsf, obj, objI, mesh);
112  addToFieldList(vvf, obj, objI, mesh);
113  addToFieldList(vsptf, obj, objI, mesh);
114  addToFieldList(vsytf, obj, objI, mesh);
115  addToFieldList(vtf, obj, objI, mesh);
116 
117  addToFieldList(psf, obj, objI, pointMesh::New(mesh));
118  addToFieldList(pvf, obj, objI, pointMesh::New(mesh));
119  addToFieldList(psptf, obj, objI, pointMesh::New(mesh));
120  addToFieldList(psytf, obj, objI, pointMesh::New(mesh));
121  addToFieldList(ptf, obj, objI, pointMesh::New(mesh));
122  }
123  }
124 
125  Info<< endl;
126 
127 
128  if (expand)
129  {
130  // Print each patch separately
131 
132  forAll(bm, patchi)
133  {
134  Info<< bm[patchi].type() << "\t: " << bm[patchi].name() << nl;
135  outputFieldList(vsf, patchi);
136  outputFieldList(vvf, patchi);
137  outputFieldList(vsptf, patchi);
138  outputFieldList(vsytf, patchi);
139  outputFieldList(vtf, patchi);
140 
141  outputFieldList(psf, patchi);
142  outputFieldList(pvf, patchi);
143  outputFieldList(psptf, patchi);
144  outputFieldList(psytf, patchi);
145  outputFieldList(ptf, patchi);
146  Info<< endl;
147  }
148  }
149  else
150  {
151  // Collect for each patch the bc type per field. Merge similar
152  // patches.
153 
154  // Per 'group', the map from fieldname to patchfield type
155  DynamicList<HashTable<word>> fieldToTypes(bm.size());
156  // Per 'group' the patches
157  DynamicList<DynamicList<label>> groupToPatches(bm.size());
158  forAll(bm, patchi)
159  {
160  HashTable<word> fieldToType;
161  collectFieldList(vsf, patchi, fieldToType);
162  collectFieldList(vvf, patchi, fieldToType);
163  collectFieldList(vsptf, patchi, fieldToType);
164  collectFieldList(vsytf, patchi, fieldToType);
165  collectFieldList(vtf, patchi, fieldToType);
166 
167  collectFieldList(psf, patchi, fieldToType);
168  collectFieldList(pvf, patchi, fieldToType);
169  collectFieldList(psptf, patchi, fieldToType);
170  collectFieldList(psytf, patchi, fieldToType);
171  collectFieldList(ptf, patchi, fieldToType);
172 
173  label groupI = findIndex(fieldToTypes, fieldToType);
174  if (groupI == -1)
175  {
176  DynamicList<label> group(1);
177  group.append(patchi);
178  groupToPatches.append(group);
179  fieldToTypes.append(fieldToType);
180  }
181  else
182  {
183  groupToPatches[groupI].append(patchi);
184  }
185  }
186 
187 
188  forAll(groupToPatches, groupI)
189  {
190  const DynamicList<label>& patchIDs = groupToPatches[groupI];
191 
192  if (patchIDs.size() > 1)
193  {
194  // Check if part of a group
195  wordList groups;
196  labelHashSet nonGroupPatches;
197  bm.matchGroups(patchIDs, groups, nonGroupPatches);
198 
199  const labelList sortedPatches(nonGroupPatches.sortedToc());
200  forAll(sortedPatches, i)
201  {
202  Info<< bm[sortedPatches[i]].type()
203  << "\t: " << bm[sortedPatches[i]].name() << nl;
204  }
205  if (groups.size())
206  {
207  forAll(groups, i)
208  {
209  Info<< "group\t: " << groups[i] << nl;
210  }
211  }
212  outputFieldList(vsf, patchIDs[0]);
213  outputFieldList(vvf, patchIDs[0]);
214  outputFieldList(vsptf, patchIDs[0]);
215  outputFieldList(vsytf, patchIDs[0]);
216  outputFieldList(vtf, patchIDs[0]);
217 
218  outputFieldList(psf, patchIDs[0]);
219  outputFieldList(pvf, patchIDs[0]);
220  outputFieldList(psptf, patchIDs[0]);
221  outputFieldList(psytf, patchIDs[0]);
222  outputFieldList(ptf, patchIDs[0]);
223  Info<< endl;
224  }
225  else
226  {
227  // No group.
228  forAll(patchIDs, i)
229  {
230  label patchi = patchIDs[i];
231  Info<< bm[patchi].type()
232  << "\t: " << bm[patchi].name() << nl;
233  outputFieldList(vsf, patchi);
234  outputFieldList(vvf, patchi);
235  outputFieldList(vsptf, patchi);
236  outputFieldList(vsytf, patchi);
237  outputFieldList(vtf, patchi);
238 
239  outputFieldList(psf, patchi);
240  outputFieldList(pvf, patchi);
241  outputFieldList(psptf, patchi);
242  outputFieldList(psytf, patchi);
243  outputFieldList(ptf, patchi);
244  Info<< endl;
245  }
246  }
247  }
248  }
249  }
250 
251  Info<< "End\n" << endl;
252 
253  return 0;
254 }
255 
256 
257 // ************************************************************************* //
void addToFieldList(PtrList< GeoField > &fieldList, const IOobject &obj, const label fieldi, const typename GeoField::Mesh &mesh)
const char *const group
Group name for atomic constants.
List< instant > instantList
List of instants.
Definition: instantList.H:42
#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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:210
dynamicFvMesh & mesh
void collectFieldList(const PtrList< GeoField > &fieldList, const label patchi, HashTable< word > &fieldToType)
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const char nl
Definition: Ostream.H:262
label findIndex(const ListType &, typename ListType::const_reference, const label start=0)
Find first occurence of given element and return index,.
List< word > wordList
A List of words.
Definition: fileName.H:54
label patchi
messageStream Info
Foam::argList args(argc, argv)
void outputFieldList(const PtrList< GeoField > &fieldList, const label patchi)