FreeFOAM The Cross-Platform CFD Toolkit
sampledSets.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 <OpenFOAM/dictionary.H>
28 #include <OpenFOAM/Time.H>
29 #include <finiteVolume/volFields.H>
30 #include <OpenFOAM/ListListOps.H>
31 #include <OpenFOAM/SortableList.H>
33 
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(sampledSets, 0);
40 }
41 
42 bool Foam::sampledSets::verbose_ = false;
43 
44 
45 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
46 
47 bool Foam::sampledSets::checkFieldTypes()
48 {
49  wordList fieldTypes(fieldNames_.size());
50 
51  // check files for a particular time
52  if (loadFromFiles_)
53  {
54  forAll(fieldNames_, fieldi)
55  {
56  IOobject io
57  (
58  fieldNames_[fieldi],
59  mesh_.time().timeName(),
60  mesh_,
63  false
64  );
65 
66  if (io.headerOk())
67  {
68  fieldTypes[fieldi] = io.headerClassName();
69  }
70  else
71  {
72  fieldTypes[fieldi] = "(notFound)";
73  }
74  }
75  }
76  else
77  {
78  // check objectRegistry
79  forAll(fieldNames_, fieldi)
80  {
82  mesh_.find(fieldNames_[fieldi]);
83 
84  if (iter != mesh_.objectRegistry::end())
85  {
86  fieldTypes[fieldi] = iter()->type();
87  }
88  else
89  {
90  fieldTypes[fieldi] = "(notFound)";
91  }
92  }
93  }
94 
95 
96  label nFields = 0;
97 
98  // classify fieldTypes
99  nFields += grep(scalarFields_, fieldTypes);
100  nFields += grep(vectorFields_, fieldTypes);
101  nFields += grep(sphericalTensorFields_, fieldTypes);
102  nFields += grep(symmTensorFields_, fieldTypes);
103  nFields += grep(tensorFields_, fieldTypes);
104 
105  if (Pstream::master())
106  {
107  if (debug)
108  {
109  Pout<< "timeName = " << mesh_.time().timeName() << nl
110  << "scalarFields " << scalarFields_ << nl
111  << "vectorFields " << vectorFields_ << nl
112  << "sphTensorFields " << sphericalTensorFields_ << nl
113  << "symTensorFields " << symmTensorFields_ <<nl
114  << "tensorFields " << tensorFields_ <<nl;
115  }
116 
117  if (nFields > 0)
118  {
119  if (debug)
120  {
121  Pout<< "Creating directory "
122  << outputPath_/mesh_.time().timeName()
123  << nl << endl;
124  }
125 
126  mkDir(outputPath_/mesh_.time().timeName());
127  }
128  }
129 
130  return nFields > 0;
131 }
132 
133 
134 void Foam::sampledSets::combineSampledSets
135 (
136  PtrList<coordSet>& masterSampledSets,
137  labelListList& indexSets
138 )
139 {
140  // Combine sampleSets from processors. Sort by curveDist. Return
141  // ordering in indexSets.
142  // Note: only master results are valid
143 
144  masterSampledSets_.clear();
145  masterSampledSets_.setSize(size());
146  indexSets_.setSize(size());
147 
148  const PtrList<sampledSet>& sampledSets = *this;
149 
150  forAll(sampledSets, seti)
151  {
152  const sampledSet& samplePts = sampledSets[seti];
153 
154  // Collect data from all processors
155  List<List<point> > gatheredPts(Pstream::nProcs());
156  gatheredPts[Pstream::myProcNo()] = samplePts;
157  Pstream::gatherList(gatheredPts);
158 
159  List<labelList> gatheredSegments(Pstream::nProcs());
160  gatheredSegments[Pstream::myProcNo()] = samplePts.segments();
161  Pstream::gatherList(gatheredSegments);
162 
163  List<scalarList> gatheredDist(Pstream::nProcs());
164  gatheredDist[Pstream::myProcNo()] = samplePts.curveDist();
165  Pstream::gatherList(gatheredDist);
166 
167 
168  // Combine processor lists into one big list.
169  List<point> allPts
170  (
171  ListListOps::combine<List<point> >
172  (
173  gatheredPts, accessOp<List<point> >()
174  )
175  );
176  labelList allSegments
177  (
178  ListListOps::combine<labelList>
179  (
180  gatheredSegments, accessOp<labelList>()
181  )
182  );
183  scalarList allCurveDist
184  (
185  ListListOps::combine<scalarList>
186  (
187  gatheredDist, accessOp<scalarList>()
188  )
189  );
190 
191 
192  if (Pstream::master() && allCurveDist.size() == 0)
193  {
194  WarningIn("sampledSets::combineSampledSets(..)")
195  << "Sample set " << samplePts.name()
196  << " has zero points." << endl;
197  }
198 
199  // Sort curveDist and use to fill masterSamplePts
200  SortableList<scalar> sortedDist(allCurveDist);
201  indexSets[seti] = sortedDist.indices();
202 
203  // Get reference point (note: only master has all points)
204  point refPt;
205 
206  if (allPts.size())
207  {
208  refPt = samplePts.getRefPoint(allPts);
209  }
210  else
211  {
212  refPt = vector::zero;
213  }
214 
215 
216  masterSampledSets.set
217  (
218  seti,
219  new coordSet
220  (
221  samplePts.name(),
222  samplePts.axis(),
223  UIndirectList<point>(allPts, indexSets[seti]),
224  refPt
225  )
226  );
227  }
228 }
229 
230 
231 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
232 
233 Foam::sampledSets::sampledSets
234 (
235  const word& name,
236  const objectRegistry& obr,
237  const dictionary& dict,
238  const bool loadFromFiles
239 )
240 :
242  name_(name),
243  mesh_(refCast<const fvMesh>(obr)),
244  loadFromFiles_(loadFromFiles),
245  outputPath_(fileName::null),
246  searchEngine_(mesh_, true),
247  fieldNames_(),
248  interpolationScheme_(word::null),
249  writeFormat_(word::null)
250 {
251  if (Pstream::parRun())
252  {
253  outputPath_ = mesh_.time().path()/".."/name_;
254  }
255  else
256  {
257  outputPath_ = mesh_.time().path()/name_;
258  }
259  if (mesh_.name() != fvMesh::defaultRegion)
260  {
261  outputPath_ = outputPath_/mesh_.name();
262  }
263 
264  read(dict);
265 }
266 
267 
268 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
269 
271 {}
272 
273 
274 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
275 
276 void Foam::sampledSets::verbose(const bool verbosity)
277 {
278  verbose_ = verbosity;
279 }
280 
281 
283 {
284  // Do nothing - only valid on write
285 }
286 
287 
289 {
290  // Do nothing - only valid on write
291 }
292 
293 
295 {
296  if (size() && checkFieldTypes())
297  {
298  sampleAndWrite(scalarFields_);
299  sampleAndWrite(vectorFields_);
300  sampleAndWrite(sphericalTensorFields_);
301  sampleAndWrite(symmTensorFields_);
302  sampleAndWrite(tensorFields_);
303  }
304 }
305 
306 
308 {
309  dict_ = dict;
310 
311  fieldNames_ = wordList(dict_.lookup("fields"));
312 
313  interpolationScheme_ = "cell";
314  dict_.readIfPresent("interpolationScheme", interpolationScheme_);
315 
316  dict_.lookup("setFormat") >> writeFormat_;
317 
318  scalarFields_.clear();
319  vectorFields_.clear();
320  sphericalTensorFields_.clear();
321  symmTensorFields_.clear();
322  tensorFields_.clear();
323 
324  PtrList<sampledSet> newList
325  (
326  dict_.lookup("sets"),
327  sampledSet::iNew(mesh_, searchEngine_)
328  );
329  transfer(newList);
330  combineSampledSets(masterSampledSets_, indexSets_);
331 
332  if (Pstream::master() && debug)
333  {
334  Pout<< "sample fields:" << fieldNames_ << nl
335  << "sample sets:" << nl << "(" << nl;
336 
337  forAll(*this, si)
338  {
339  Pout << " " << operator[](si) << endl;
340  }
341  Pout << ")" << endl;
342  }
343 }
344 
345 
347 {
348  // reset interpolation
349  pointMesh::Delete(mesh_);
351 
352  searchEngine_.correct();
353 
354  PtrList<sampledSet> newList
355  (
356  dict_.lookup("sets"),
357  sampledSet::iNew(mesh_, searchEngine_)
358  );
359  transfer(newList);
360  combineSampledSets(masterSampledSets_, indexSets_);
361 }
362 
363 
365 {
366  correct();
367 }
368 
369 
371 {
372  correct();
373 }
374 
375 
377 {
378  if (state != polyMesh::UNCHANGED)
379  {
380  correct();
381  }
382 }
383 
384 
385 // ************************ vim: set sw=4 sts=4 et: ************************ //