fieldToCell.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2021 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 \*---------------------------------------------------------------------------*/
25 
26 #include "fieldToCell.H"
27 #include "polyMesh.H"
28 #include "cellSet.H"
29 #include "Time.H"
30 #include "IFstream.H"
31 #include "fieldDictionary.H"
32 #include "volFields.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(fieldToCell, 0);
40  addToRunTimeSelectionTable(topoSetSource, fieldToCell, word);
41 }
42 
43 
44 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
45 
46 void Foam::fieldToCell::applyToSet
47 (
48  const topoSetSource::setAction action,
49  const scalarField& field,
50  topoSet& set
51 ) const
52 {
53  Info<< " Field min:" << min(field)
54  << " max:" << max(field) << endl;
55 
56  if ((action == topoSetSource::NEW) || (action == topoSetSource::ADD))
57  {
58  Info<< " Adding all cells with value of field " << fieldName_
59  << " within range " << min_ << ".." << max_ << endl;
60 
61  forAll(field, celli)
62  {
63  if (field[celli] >= min_ && field[celli] <= max_)
64  {
65  set.insert(celli);
66  }
67  }
68  }
69  else if (action == topoSetSource::DELETE)
70  {
71  Info<< " Removing all cells with value of field " << fieldName_
72  << " within range " << min_ << ".." << max_ << endl;
73 
74  forAll(field, celli)
75  {
76  if (field[celli] >= min_ && field[celli] <= max_)
77  {
78  set.erase(celli);
79  }
80  }
81  }
82 }
83 
84 
85 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
86 
88 (
89  const polyMesh& mesh,
90  const word& fieldName,
91  const scalar min,
92  const scalar max
93 )
94 :
95  topoSetSource(mesh),
96  fieldName_(fieldName),
97  min_(min),
98  max_(max)
99 {}
100 
101 
103 (
104  const polyMesh& mesh,
105  const dictionary& dict
106 )
107 :
108  topoSetSource(mesh),
109  fieldName_(dict.lookup("field")),
110  min_(dict.lookup<scalar>("min")),
111  max_(dict.lookup<scalar>("max"))
112 {}
113 
114 
115 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
116 
118 {}
119 
120 
121 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
122 
123 void Foam::fieldToCell::applyToSet
124 (
125  const topoSetSource::setAction action,
126  topoSet& set
127 ) const
128 {
130  (
131  fieldName_,
132  mesh().time().timeName(),
133  mesh(),
136  false
137  );
138 
139  if (!fieldObject.IOobject::headerOk())
140  {
142  << "Cannot read field " << fieldName_
143  << " from time " << mesh().time().timeName() << endl;
144  }
145  else if (fieldObject.IOobject::headerOk())
146  {
147  IFstream str(fieldObject.filePath());
148 
149  // Read dictionary
150  const fieldDictionary fieldDict
151  (
152  fieldObject,
154  );
155 
156  const scalarField internalVals
157  (
158  "internalField",
159  fieldDict, mesh().nCells()
160  );
161 
162  applyToSet(action, internalVals, set);
163  }
164  else
165  {
167  (
168  fieldName_,
169  mesh().time().timeName(),
170  mesh(),
173  false
174  );
175 
176  if (fieldObject.IOobject::headerOk())
177  {
178  IFstream str(fieldObject.filePath());
179 
180  // Read dictionary
181  const fieldDictionary fieldDict
182  (
183  fieldObject,
185  );
186 
187  const vectorField internalVals
188  (
189  "internalField",
190  fieldDict, mesh().nCells()
191  );
192 
193  applyToSet(action, mag(internalVals), set);
194  }
195  else
196  {
198  << "Cannot handle fields of type " << fieldObject.headerClassName()
199  << endl;
200  }
201  }
202 }
203 
204 
205 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
virtual ~fieldToCell()
Destructor.
Definition: fieldToCell.C:117
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
static const char *const typeName
Definition: Field.H:105
Templated form of IOobject providing type information for file reading and header type checking...
Definition: IOobject.H:537
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
IOobject fieldObject(fieldNames[var2field[nVar]], runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE)
fvMesh & mesh
fileName filePath() const
Return the path for the file for this Type.
Definition: IOobject.H:572
Macros for easy insertion into run-time selection tables.
Base class of a source for a topoSet.
Definition: topoSetSource.H:63
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
setAction
Enumeration defining the valid actions.
Definition: topoSetSource.H:82
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
word timeName
Definition: getTimeIndex.H:3
Read field as dictionary (without mesh).
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
Input from file stream.
Definition: IFstream.H:81
General set of labels of mesh quantity (points, cells, faces).
Definition: topoSet.H:61
#define WarningInFunction
Report a warning using Foam::Warning.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
fieldToCell(const polyMesh &mesh, const word &fieldName, const scalar min, const scalar max)
Construct from components.
Definition: fieldToCell.C:88
const word & headerClassName() const
Return name of the class name read from header.
Definition: IOobject.H:321
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864