patchProbesTemplates.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-2022 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().userTimeValue();
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().userTimeValue();
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  objectRegistry::const_iterator iter = mesh_.find(fields[fieldi]);
94 
95  if
96  (
97  iter != objectRegistry::end()
98  && iter()->type()
100  )
101  {
102  sampleAndWrite
103  (
105  <GeometricField<Type, fvPatchField, volMesh>>
106  (
107  fields[fieldi]
108  )
109  );
110  }
111  }
112 }
113 
114 
115 template<class Type>
116 void Foam::patchProbes::sampleAndWriteSurfaceFields
117 (
118  const fieldGroup<Type>& fields
119 )
120 {
121  forAll(fields, fieldi)
122  {
123  objectRegistry::const_iterator iter = mesh_.find(fields[fieldi]);
124 
125  if
126  (
127  iter != objectRegistry::end()
128  && iter()->type()
130  )
131  {
132  sampleAndWrite
133  (
135  <GeometricField<Type, fvsPatchField, surfaceMesh>>
136  (
137  fields[fieldi]
138  )
139  );
140  }
141  }
142 }
143 
144 
145 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
146 
147 template<class Type>
149 Foam::patchProbes::sample
150 (
152 ) const
153 {
154  const Type unsetVal(-vGreat*pTraits<Type>::one);
155 
156  tmp<Field<Type>> tValues
157  (
158  new Field<Type>(this->size(), unsetVal)
159  );
160 
161  Field<Type>& values = tValues.ref();
162 
164 
165  forAll(*this, probei)
166  {
167  label facei = elementList_[probei];
168 
169  if (facei >= 0)
170  {
171  label patchi = patches.whichPatch(facei);
172  label localFacei = patches[patchi].whichFace(facei);
173  values[probei] = vField.boundaryField()[patchi][localFacei];
174  }
175  }
176 
179 
180  return tValues;
181 }
182 
183 
184 template<class Type>
186 Foam::patchProbes::sample(const word& fieldName) const
187 {
188  return sample
189  (
191  (
192  fieldName
193  )
194  );
195 }
196 
197 
198 template<class Type>
200 Foam::patchProbes::sample
201 (
203 ) const
204 {
205  const Type unsetVal(-vGreat*pTraits<Type>::one);
206 
207  tmp<Field<Type>> tValues
208  (
209  new Field<Type>(this->size(), unsetVal)
210  );
211 
212  Field<Type>& values = tValues.ref();
213 
215 
216  forAll(*this, probei)
217  {
218  label facei = elementList_[probei];
219 
220  if (facei >= 0)
221  {
222  label patchi = patches.whichPatch(facei);
223  label localFacei = patches[patchi].whichFace(facei);
224  values[probei] = sField.boundaryField()[patchi][localFacei];
225  }
226  }
227 
230 
231  return tValues;
232 }
233 // ************************************************************************* //
const fvPatchList & patches
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:453
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:434
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
static const char *const typeName
Definition: Field.H:105
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:458
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Traits class for primitives.
Definition: pTraits.H:50
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
Generic GeometricField class.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
static void listCombineGather(const List< commsStruct > &comms, List< T > &Value, const CombineOp &cop, const int tag, const label comm)
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
Pre-declare SubField and related Field type.
Definition: Field.H:56
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:197
virtual const word & type() const =0
Runtime type information.
const fvMesh & mesh_
Const reference to fvMesh.
Definition: probes.H:91
Foam::polyBoundaryMesh.
Istream and Ostream manipulators taking arguments.
HashPtrTable< OFstream > probeFilePtrs_
Current open files.
Definition: probes.H:132
label patchi
labelList elementList_
Definition: probes.H:126
Omanip< int > setw(const int i)
Definition: IOmanip.H:199
A class for managing temporary objects.
Definition: PtrList.H:53
label size() const
Return the number of elements in the UList.
Definition: ListI.H:171
label whichPatch(const label faceIndex) const
Return patch index for a given face label.