FreeFOAM The Cross-Platform CFD Toolkit
sampledSetsTemplates.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) 1991-2010 OpenCFD Ltd.
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 "sampledSets.H"
27 #include <finiteVolume/volFields.H>
28 #include <OpenFOAM/ListListOps.H>
29 
30 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31 
32 template <class Type>
33 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
34 (
35  const word& interpolationScheme,
36  const GeometricField<Type, fvPatchField, volMesh>& field,
37  const PtrList<sampledSet>& samplers
38 )
39 :
40  List<Field<Type> >(samplers.size()),
41  name_(field.name())
42 {
43  autoPtr<interpolation<Type> > interpolator
44  (
45  interpolation<Type>::New(interpolationScheme, field)
46  );
47 
48  forAll(samplers, seti)
49  {
50  Field<Type>& values = this->operator[](seti);
51  const sampledSet& samples = samplers[seti];
52 
53  values.setSize(samples.size());
54  forAll(samples, samplei)
55  {
56  const point& samplePt = samples[samplei];
57  label celli = samples.cells()[samplei];
58  label facei = samples.faces()[samplei];
59 
60  if (celli == -1 && facei == -1)
61  {
62  // Special condition for illegal sampling points
63  values[samplei] = pTraits<Type>::max;
64  }
65  else
66  {
67  values[samplei] = interpolator().interpolate
68  (
69  samplePt,
70  celli,
71  facei
72  );
73  }
74  }
75  }
76 }
77 
78 
79 template <class Type>
80 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
81 (
82  const GeometricField<Type, fvPatchField, volMesh>& field,
83  const PtrList<sampledSet>& samplers
84 )
85 :
86  List<Field<Type> >(samplers.size()),
87  name_(field.name())
88 {
89  forAll(samplers, seti)
90  {
91  Field<Type>& values = this->operator[](seti);
92  const sampledSet& samples = samplers[seti];
93 
94  values.setSize(samples.size());
95  forAll(samples, samplei)
96  {
97  label celli = samples.cells()[samplei];
98 
99  if (celli ==-1)
100  {
101  values[samplei] = pTraits<Type>::max;
102  }
103  else
104  {
105  values[samplei] = field[celli];
106  }
107  }
108  }
109 }
110 
111 
112 template <class Type>
113 Foam::sampledSets::volFieldSampler<Type>::volFieldSampler
114 (
115  const List<Field<Type> >& values,
116  const word& name
117 )
118 :
119  List<Field<Type> >(values),
120  name_(name)
121 {}
122 
123 
124 template<class Type>
125 Foam::label Foam::sampledSets::grep
126 (
127  fieldGroup<Type>& fieldList,
128  const wordList& fieldTypes
129 ) const
130 {
131  fieldList.setSize(fieldNames_.size());
132  label nFields = 0;
133 
134  forAll(fieldNames_, fieldi)
135  {
136  if
137  (
138  fieldTypes[fieldi]
140  )
141  {
142  fieldList[nFields] = fieldNames_[fieldi];
143  nFields++;
144  }
145  }
146 
147  fieldList.setSize(nFields);
148 
149  return nFields;
150 }
151 
152 
153 template<class Type>
154 void Foam::sampledSets::writeSampleFile
155 (
156  const coordSet& masterSampleSet,
157  const PtrList<volFieldSampler<Type> >& masterFields,
158  const label seti,
159  const fileName& timeDir,
160  const writer<Type>& formatter
161 )
162 {
163  wordList valueSetNames(masterFields.size());
164  List<const Field<Type>*> valueSets(masterFields.size());
165 
166  forAll(masterFields, fieldi)
167  {
168  valueSetNames[fieldi] = masterFields[fieldi].name();
169  valueSets[fieldi] = &masterFields[fieldi][seti];
170  }
171 
172  fileName fName
173  (
174  timeDir/formatter.getFileName(masterSampleSet, valueSetNames)
175  );
176 
177  OFstream ofs(fName);
178  if (ofs.opened())
179  {
180  formatter.write
181  (
182  masterSampleSet,
183  valueSetNames,
184  valueSets,
185  ofs
186  );
187  }
188  else
189  {
190  WarningIn
191  (
192  "void Foam::sampledSets::writeSampleFile"
193  "("
194  "const coordSet&, "
195  "const PtrList<volFieldSampler<Type> >&, "
196  "const label, "
197  "const fileName&, "
198  "const writer<Type>&"
199  ")"
200  ) << "File " << ofs.name() << " could not be opened. "
201  << "No data will be written" << endl;
202  }
203 }
204 
205 
206 template<class T>
207 void Foam::sampledSets::combineSampledValues
208 (
209  const PtrList<volFieldSampler<T> >& sampledFields,
210  const labelListList& indexSets,
211  PtrList<volFieldSampler<T> >& masterFields
212 )
213 {
214  forAll(sampledFields, fieldi)
215  {
216  List<Field<T> > masterValues(indexSets.size());
217 
218  forAll(indexSets, seti)
219  {
220  // Collect data from all processors
221  List<Field<T> > gatheredData(Pstream::nProcs());
222  gatheredData[Pstream::myProcNo()] = sampledFields[fieldi][seti];
223  Pstream::gatherList(gatheredData);
224 
225  if (Pstream::master())
226  {
227  Field<T> allData
228  (
229  ListListOps::combine<Field<T> >
230  (
231  gatheredData,
232  Foam::accessOp<Field<T> >()
233  )
234  );
235 
236  masterValues[seti] = UIndirectList<T>
237  (
238  allData,
239  indexSets[seti]
240  )();
241  }
242  }
243 
244  masterFields.set
245  (
246  fieldi,
247  new volFieldSampler<T>
248  (
249  masterValues,
250  sampledFields[fieldi].name()
251  )
252  );
253  }
254 }
255 
256 
257 template<class Type>
258 void Foam::sampledSets::sampleAndWrite
259 (
260  fieldGroup<Type>& fields
261 )
262 {
263  if (fields.size())
264  {
265  bool interpolate = interpolationScheme_ != "cell";
266 
267  // Create or use existing writer
268  if (fields.formatter.empty())
269  {
270  fields.formatter = writer<Type>::New(writeFormat_);
271  }
272 
273  // Storage for interpolated values
274  PtrList<volFieldSampler<Type> > sampledFields(fields.size());
275 
276  forAll(fields, fieldi)
277  {
278  if (Pstream::master() && verbose_)
279  {
280  Pout<< "sampledSets::sampleAndWrite: "
281  << fields[fieldi] << endl;
282  }
283 
284  if (loadFromFiles_)
285  {
286  GeometricField<Type, fvPatchField, volMesh> vf
287  (
288  IOobject
289  (
290  fields[fieldi],
291  mesh_.time().timeName(),
292  mesh_,
295  false
296  ),
297  mesh_
298  );
299 
300  if (interpolate)
301  {
302  sampledFields.set
303  (
304  fieldi,
305  new volFieldSampler<Type>
306  (
307  interpolationScheme_,
308  vf,
309  *this
310  )
311  );
312  }
313  else
314  {
315  sampledFields.set
316  (
317  fieldi,
318  new volFieldSampler<Type>(vf, *this)
319  );
320  }
321  }
322  else
323  {
324  if (interpolate)
325  {
326  sampledFields.set
327  (
328  fieldi,
329  new volFieldSampler<Type>
330  (
331  interpolationScheme_,
332  mesh_.lookupObject
333  <GeometricField<Type, fvPatchField, volMesh> >
334  (fields[fieldi]),
335  *this
336  )
337  );
338  }
339  else
340  {
341  sampledFields.set
342  (
343  fieldi,
344  new volFieldSampler<Type>
345  (
346  mesh_.lookupObject
347  <GeometricField<Type, fvPatchField, volMesh> >
348  (fields[fieldi]),
349  *this
350  )
351  );
352  }
353  }
354  }
355 
356  // Combine sampled fields from processors.
357  // Note: only master results are valid
358 
359  PtrList<volFieldSampler<Type> > masterFields(sampledFields.size());
360  combineSampledValues(sampledFields, indexSets_, masterFields);
361 
362  if (Pstream::master())
363  {
364  forAll(masterSampledSets_, seti)
365  {
366  writeSampleFile
367  (
368  masterSampledSets_[seti],
369  masterFields,
370  seti,
371  outputPath_/mesh_.time().timeName(),
372  fields.formatter()
373  );
374  }
375  }
376  }
377 }
378 
379 
380 // ************************ vim: set sw=4 sts=4 et: ************************ //