setFields.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 Description
25  Set values on a selected set of cells/patchfaces through a dictionary.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "argList.H"
30 #include "Time.H"
31 #include "fvMesh.H"
32 #include "topoSetSource.H"
33 #include "cellSet.H"
34 #include "faceSet.H"
35 #include "volFields.H"
36 
37 using namespace Foam;
38 
39 template<class Type>
40 bool setCellFieldType
41 (
42  const word& fieldTypeDesc,
43  const fvMesh& mesh,
44  const labelList& selectedCells,
45  Istream& fieldValueStream
46 )
47 {
49 
50  if (fieldTypeDesc != fieldType::typeName + "Value")
51  {
52  return false;
53  }
54 
55  word fieldName(fieldValueStream);
56 
57  // Check the current time directory
58  IOobject fieldHeader
59  (
60  fieldName,
61  mesh.time().timeName(),
62  mesh,
64  );
65 
66  // Check the "constant" directory
67  if (!fieldHeader.headerOk())
68  {
69  fieldHeader = IOobject
70  (
71  fieldName,
72  mesh.time().constant(),
73  mesh,
75  );
76  }
77 
78  // Check field exists
79  if (fieldHeader.headerOk())
80  {
81  Info<< " Setting internal values of "
82  << fieldHeader.headerClassName()
83  << " " << fieldName << endl;
84 
85  fieldType field(fieldHeader, mesh);
86 
87  const Type& value = pTraits<Type>(fieldValueStream);
88 
89  if (selectedCells.size() == field.size())
90  {
91  field.primitiveFieldRef() = value;
92  }
93  else
94  {
95  forAll(selectedCells, celli)
96  {
97  field[selectedCells[celli]] = value;
98  }
99  }
100 
102  Boundary& fieldBf = field.boundaryFieldRef();
103 
104  forAll(field.boundaryField(), patchi)
105  {
106  fieldBf[patchi] = fieldBf[patchi].patchInternalField();
107  }
108 
109  if (!field.write())
110  {
112  << "Failed writing field " << fieldName << endl;
113  }
114  }
115  else
116  {
118  << "Field " << fieldName << " not found" << endl;
119 
120  // Consume value
121  (void)pTraits<Type>(fieldValueStream);
122  }
123 
124  return true;
125 }
126 
127 
128 class setCellField
129 {
130 
131 public:
132 
133  setCellField()
134  {}
135 
136  autoPtr<setCellField> clone() const
137  {
138  return autoPtr<setCellField>(new setCellField());
139  }
140 
141  class iNew
142  {
143  const fvMesh& mesh_;
144  const labelList& selectedCells_;
145 
146  public:
147 
148  iNew(const fvMesh& mesh, const labelList& selectedCells)
149  :
150  mesh_(mesh),
151  selectedCells_(selectedCells)
152  {}
153 
154  autoPtr<setCellField> operator()(Istream& fieldValues) const
155  {
156  word fieldType(fieldValues);
157 
158  if
159  (
160  !(
161  setCellFieldType<scalar>
162  (fieldType, mesh_, selectedCells_, fieldValues)
163  || setCellFieldType<vector>
164  (fieldType, mesh_, selectedCells_, fieldValues)
165  || setCellFieldType<sphericalTensor>
166  (fieldType, mesh_, selectedCells_, fieldValues)
167  || setCellFieldType<symmTensor>
168  (fieldType, mesh_, selectedCells_, fieldValues)
169  || setCellFieldType<tensor>
170  (fieldType, mesh_, selectedCells_, fieldValues)
171  )
172  )
173  {
175  << "field type " << fieldType << " not currently supported"
176  << endl;
177  }
178 
179  return autoPtr<setCellField>(new setCellField());
180  }
181  };
182 };
183 
184 
185 template<class Type>
186 bool setFaceFieldType
187 (
188  const word& fieldTypeDesc,
189  const fvMesh& mesh,
190  const labelList& selectedFaces,
191  Istream& fieldValueStream
192 )
193 {
195 
196  if (fieldTypeDesc != fieldType::typeName + "Value")
197  {
198  return false;
199  }
200 
201  word fieldName(fieldValueStream);
202 
203  // Check the current time directory
204  IOobject fieldHeader
205  (
206  fieldName,
207  mesh.time().timeName(),
208  mesh,
210  );
211 
212  // Check the "constant" directory
213  if (!fieldHeader.headerOk())
214  {
215  fieldHeader = IOobject
216  (
217  fieldName,
218  mesh.time().constant(),
219  mesh,
221  );
222  }
223 
224  // Check field exists
225  if (fieldHeader.headerOk())
226  {
227  Info<< " Setting patchField values of "
228  << fieldHeader.headerClassName()
229  << " " << fieldName << endl;
230 
231  fieldType field(fieldHeader, mesh);
232 
233  const Type& value = pTraits<Type>(fieldValueStream);
234 
235  // Create flat list of selected faces and their value.
236  Field<Type> allBoundaryValues(mesh.nFaces()-mesh.nInternalFaces());
237  forAll(field.boundaryField(), patchi)
238  {
240  (
241  allBoundaryValues,
242  field.boundaryField()[patchi].size(),
243  field.boundaryField()[patchi].patch().start()
244  - mesh.nInternalFaces()
245  ) = field.boundaryField()[patchi];
246  }
247 
248  // Override
249  bool hasWarned = false;
250  labelList nChanged
251  (
252  returnReduce(field.boundaryField().size(), maxOp<label>()),
253  0
254  );
255  forAll(selectedFaces, i)
256  {
257  label facei = selectedFaces[i];
258  if (mesh.isInternalFace(facei))
259  {
260  if (!hasWarned)
261  {
262  hasWarned = true;
264  << "Ignoring internal face " << facei
265  << ". Suppressing further warnings." << endl;
266  }
267  }
268  else
269  {
270  label bFacei = facei-mesh.nInternalFaces();
271  allBoundaryValues[bFacei] = value;
272  nChanged[mesh.boundaryMesh().patchID()[bFacei]]++;
273  }
274  }
275 
277  Pstream::listCombineScatter(nChanged);
278 
280  Boundary& fieldBf = field.boundaryFieldRef();
281 
282  // Reassign.
283  forAll(field.boundaryField(), patchi)
284  {
285  if (nChanged[patchi] > 0)
286  {
287  Info<< " On patch "
288  << field.boundaryField()[patchi].patch().name()
289  << " set " << nChanged[patchi] << " values" << endl;
290  fieldBf[patchi] == SubField<Type>
291  (
292  allBoundaryValues,
293  fieldBf[patchi].size(),
294  fieldBf[patchi].patch().start()
295  - mesh.nInternalFaces()
296  );
297  }
298  }
299 
300  if (!field.write())
301  {
303  << "Failed writing field " << field.name() << exit(FatalError);
304  }
305  }
306  else
307  {
309  << "Field " << fieldName << " not found" << endl;
310 
311  // Consume value
312  (void)pTraits<Type>(fieldValueStream);
313  }
314 
315  return true;
316 }
317 
318 
319 class setFaceField
320 {
321 
322 public:
323 
324  setFaceField()
325  {}
326 
327  autoPtr<setFaceField> clone() const
328  {
329  return autoPtr<setFaceField>(new setFaceField());
330  }
331 
332  class iNew
333  {
334  const fvMesh& mesh_;
335  const labelList& selectedFaces_;
336 
337  public:
338 
339  iNew(const fvMesh& mesh, const labelList& selectedFaces)
340  :
341  mesh_(mesh),
342  selectedFaces_(selectedFaces)
343  {}
344 
345  autoPtr<setFaceField> operator()(Istream& fieldValues) const
346  {
347  word fieldType(fieldValues);
348 
349  if
350  (
351  !(
352  setFaceFieldType<scalar>
353  (fieldType, mesh_, selectedFaces_, fieldValues)
354  || setFaceFieldType<vector>
355  (fieldType, mesh_, selectedFaces_, fieldValues)
356  || setFaceFieldType<sphericalTensor>
357  (fieldType, mesh_, selectedFaces_, fieldValues)
358  || setFaceFieldType<symmTensor>
359  (fieldType, mesh_, selectedFaces_, fieldValues)
360  || setFaceFieldType<tensor>
361  (fieldType, mesh_, selectedFaces_, fieldValues)
362  )
363  )
364  {
366  << "field type " << fieldType << " not currently supported"
367  << endl;
368  }
369 
370  return autoPtr<setFaceField>(new setFaceField());
371  }
372  };
373 };
374 
375 
376 
377 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
378 
379 int main(int argc, char *argv[])
380 {
381  #include "addRegionOption.H"
382  #include "setRootCase.H"
383  #include "createTime.H"
384  #include "createNamedMesh.H"
385 
386  Info<< "Reading setFieldsDict\n" << endl;
387 
388  IOdictionary setFieldsDict
389  (
390  IOobject
391  (
392  "setFieldsDict",
393  runTime.system(),
394  mesh,
397  )
398  );
399 
400  if (setFieldsDict.found("defaultFieldValues"))
401  {
402  Info<< "Setting field default values" << endl;
403  PtrList<setCellField> defaultFieldValues
404  (
405  setFieldsDict.lookup("defaultFieldValues"),
406  setCellField::iNew(mesh, labelList(mesh.nCells()))
407  );
408  Info<< endl;
409  }
410 
411 
412  Info<< "Setting field region values" << endl;
413 
414  PtrList<entry> regions(setFieldsDict.lookup("regions"));
415 
416  forAll(regions, regionI)
417  {
418  const entry& region = regions[regionI];
419 
420  autoPtr<topoSetSource> source =
421  topoSetSource::New(region.keyword(), mesh, region.dict());
422 
423  if (source().setType() == topoSetSource::CELLSETSOURCE)
424  {
425  cellSet selectedCellSet
426  (
427  mesh,
428  "cellSet",
429  mesh.nCells()/10+1 // Reasonable size estimate.
430  );
431 
432  source->applyToSet
433  (
435  selectedCellSet
436  );
437 
438  PtrList<setCellField> fieldValues
439  (
440  region.dict().lookup("fieldValues"),
441  setCellField::iNew(mesh, selectedCellSet.toc())
442  );
443  }
444  else if (source().setType() == topoSetSource::FACESETSOURCE)
445  {
446  faceSet selectedFaceSet
447  (
448  mesh,
449  "faceSet",
450  (mesh.nFaces()-mesh.nInternalFaces())/10+1
451  );
452 
453  source->applyToSet
454  (
456  selectedFaceSet
457  );
458 
459  PtrList<setFaceField> fieldValues
460  (
461  region.dict().lookup("fieldValues"),
462  setFaceField::iNew(mesh, selectedFaceSet.toc())
463  );
464  }
465  }
466 
467 
468  Info<< "\nEnd\n" << endl;
469 
470  return 0;
471 }
472 
473 
474 // ************************************************************************* //
const labelList & patchID() const
Per boundary face label the patch index.
static void listCombineScatter(const List< commsStruct > &comms, List< T > &Value, const int tag, const label comm)
Scatter data. Reverse of combineGather.
#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
A list of face labels.
Definition: faceSet.H:48
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Traits class for primitives.
Definition: pTraits.H:50
virtual const dictionary & dict() const =0
Return dictionary if this entry is a dictionary.
Generic GeometricField class.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
Pre-declare related SubField type.
Definition: Field.H:61
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
label nCells() const
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
dynamicFvMesh & mesh
Pre-declare SubField and related Field type.
Definition: Field.H:57
A class for handling words, derived from string.
Definition: word.H:59
const keyType & keyword() const
Return keyword.
Definition: entry.H:120
const word & constant() const
Return constant name.
Definition: TimePaths.H:124
List< label > labelList
A List of labels.
Definition: labelList.H:56
List< Key > toc() const
Return the table of contents.
Definition: HashTable.C:198
static autoPtr< topoSetSource > New(const word &topoSetSourceType, const polyMesh &mesh, const dictionary &dict)
Return a reference to the selected topoSetSource.
Definition: topoSetSource.C:74
label nFaces() const
label size() const
Return the number of elements in the UList.
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
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
A collection of cell labels.
Definition: cellSet.H:48
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
messageStream Info
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual void applyToSet(const setAction action, topoSet &) const =0
label nInternalFaces() const
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
A keyword and a list of tokens is an &#39;entry&#39;.
Definition: entry.H:65
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451