patchProbesTemplates.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 \*---------------------------------------------------------------------------*/
25 
26 #include "patchProbes.H"
27 #include "volFields.H"
28 #include "IOmanip.H"
29 
30 
31 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
32 
33 template<class Type>
34 void Foam::patchProbes::sampleAndWrite
35 (
36  const GeometricField<Type, fvPatchField, volMesh>& vField
37 )
38 {
39  Field<Type> values(sample(vField));
40 
41  if (Pstream::master())
42  {
43  unsigned int w = IOstream::defaultPrecision() + 7;
44  OFstream& probeStream = *probeFilePtrs_[vField.name()];
45 
46  probeStream
47  << setw(w)
48  << vField.time().timeToUserTime(vField.time().value());
49 
50  forAll(values, probei)
51  {
52  probeStream << ' ' << setw(w) << values[probei];
53  }
54  probeStream << endl;
55  }
56 }
57 
58 
59 template<class Type>
60 void Foam::patchProbes::sampleAndWrite
61 (
62  const GeometricField<Type, fvsPatchField, surfaceMesh>& sField
63 )
64 {
65  Field<Type> values(sample(sField));
66 
67  if (Pstream::master())
68  {
69  unsigned int w = IOstream::defaultPrecision() + 7;
70  OFstream& probeStream = *probeFilePtrs_[sField.name()];
71 
72  probeStream
73  << setw(w)
74  << sField.time().timeToUserTime(sField.time().value());
75 
76  forAll(values, probei)
77  {
78  probeStream << ' ' << setw(w) << values[probei];
79  }
80  probeStream << endl;
81  }
82 }
83 
84 
85 template<class Type>
86 void Foam::patchProbes::sampleAndWrite
87 (
88  const fieldGroup<Type>& fields
89 )
90 {
91  forAll(fields, fieldi)
92  {
93  if (loadFromFiles_)
94  {
95  sampleAndWrite
96  (
97  GeometricField<Type, fvPatchField, volMesh>
98  (
99  IOobject
100  (
101  fields[fieldi],
102  mesh_.time().timeName(),
103  mesh_,
106  false
107  ),
108  mesh_
109  )
110  );
111  }
112  else
113  {
114  objectRegistry::const_iterator iter = mesh_.find(fields[fieldi]);
115 
116  if
117  (
118  iter != objectRegistry::end()
119  && iter()->type()
121  )
122  {
123  sampleAndWrite
124  (
126  <GeometricField<Type, fvPatchField, volMesh>>
127  (
128  fields[fieldi]
129  )
130  );
131  }
132  }
133  }
134 }
135 
136 
137 template<class Type>
138 void Foam::patchProbes::sampleAndWriteSurfaceFields
139 (
140  const fieldGroup<Type>& fields
141 )
142 {
143  forAll(fields, fieldi)
144  {
145  if (loadFromFiles_)
146  {
147  sampleAndWrite
148  (
149  GeometricField<Type, fvsPatchField, surfaceMesh>
150  (
151  IOobject
152  (
153  fields[fieldi],
154  mesh_.time().timeName(),
155  mesh_,
158  false
159  ),
160  mesh_
161  )
162  );
163  }
164  else
165  {
166  objectRegistry::const_iterator iter = mesh_.find(fields[fieldi]);
167 
168  if
169  (
170  iter != objectRegistry::end()
171  && iter()->type()
173  )
174  {
175  sampleAndWrite
176  (
178  <GeometricField<Type, fvsPatchField, surfaceMesh>>
179  (
180  fields[fieldi]
181  )
182  );
183  }
184  }
185  }
186 }
187 
188 
189 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
190 
191 template<class Type>
193 Foam::patchProbes::sample
194 (
196 ) const
197 {
198  const Type unsetVal(-VGREAT*pTraits<Type>::one);
199 
200  tmp<Field<Type>> tValues
201  (
202  new Field<Type>(this->size(), unsetVal)
203  );
204 
205  Field<Type>& values = tValues.ref();
206 
208 
209  forAll(*this, probei)
210  {
211  label facei = elementList_[probei];
212 
213  if (facei >= 0)
214  {
215  label patchi = patches.whichPatch(facei);
216  label localFacei = patches[patchi].whichFace(facei);
217  values[probei] = vField.boundaryField()[patchi][localFacei];
218  }
219  }
220 
223 
224  return tValues;
225 }
226 
227 
228 template<class Type>
230 Foam::patchProbes::sample(const word& fieldName) const
231 {
232  return sample
233  (
235  (
236  fieldName
237  )
238  );
239 }
240 
241 
242 template<class Type>
244 Foam::patchProbes::sample
245 (
247 ) const
248 {
249  const Type unsetVal(-VGREAT*pTraits<Type>::one);
250 
251  tmp<Field<Type>> tValues
252  (
253  new Field<Type>(this->size(), unsetVal)
254  );
255 
256  Field<Type>& values = tValues.ref();
257 
259 
260  forAll(*this, probei)
261  {
262  label facei = elementList_[probei];
263 
264  if (facei >= 0)
265  {
266  label patchi = patches.whichPatch(facei);
267  label localFacei = patches[patchi].whichFace(facei);
268  values[probei] = sField.boundaryField()[patchi][localFacei];
269  }
270  }
271 
274 
275  return tValues;
276 }
277 // ************************************************************************* //
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
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
static const char *const typeName
Definition: Field.H:94
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Traits class for primitives.
Definition: pTraits.H:50
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:411
Generic GeometricField class.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
patches[0]
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
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:138
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Pre-declare SubField and related Field type.
Definition: Field.H:57
A class for handling words, derived from string.
Definition: word.H:59
friend class const_iterator
Declare friendship with the const_iterator.
Definition: HashTable.H:192
virtual const word & type() const =0
Runtime type information.
const fvMesh & mesh_
Const reference to fvMesh.
Definition: probes.H:95
Foam::polyBoundaryMesh.
Istream and Ostream manipulators taking arguments.
label size() const
Return the number of elements in the UList.
HashPtrTable< OFstream > probeFilePtrs_
Current open files.
Definition: probes.H:139
label patchi
bool loadFromFiles_
Load fields from files (not from objectRegistry)
Definition: probes.H:98
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
labelList elementList_
Definition: probes.H:133
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
A class for managing temporary objects.
Definition: PtrList.H:54
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243