regionSizeDistribution.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) 2013-2015 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 "regionSizeDistribution.H"
27 #include "volFields.H"
28 #include "regionSplit.H"
29 #include "fvcVolumeIntegrate.H"
30 #include "mathematicalConstants.H"
31 #include "stringListOps.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(regionSizeDistribution, 0);
38 
39  //- Plus op for FixedList<scalar>
40  template<class T, unsigned Size>
42  {
43  public:
44  void operator()
45  (
47  const FixedList<T, Size>& y
48  ) const
49  {
50  forAll(x, i)
51  {
52  x[i] += y[i];
53  }
54  }
55  };
56 }
57 
58 
59 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
60 
61 void Foam::regionSizeDistribution::writeGraph
62 (
63  const coordSet& coords,
64  const word& valueName,
65  const scalarField& values
66 ) const
67 {
68  const wordList valNames(1, valueName);
69 
70  fileName outputPath = baseTimeDir();
71  mkDir(outputPath);
72 
73  OFstream str(outputPath/formatterPtr_().getFileName(coords, valNames));
74 
75  Info<< "Writing distribution of " << valueName << " to " << str.name()
76  << endl;
77 
78  List<const scalarField*> valPtrs(1);
79  valPtrs[0] = &values;
80  formatterPtr_().write(coords, valNames, valPtrs, str);
81 }
82 
83 
84 void Foam::regionSizeDistribution::writeAlphaFields
85 (
86  const regionSplit& regions,
87  const Map<label>& patchRegions,
88  const Map<scalar>& regionVolume,
89  const volScalarField& alpha
90 ) const
91 {
92  const scalar maxDropletVol = 1.0/6.0*pow(maxDiam_, 3);
93 
94  // Split alpha field
95  // ~~~~~~~~~~~~~~~~~
96  // Split into
97  // - liquidCore : region connected to inlet patches
98  // - per region a volume : for all other regions
99  // - backgroundAlpha : remaining alpha
100 
101 
102  // Construct field
103  volScalarField liquidCore
104  (
105  IOobject
106  (
107  alphaName_ + "_liquidCore",
108  obr_.time().timeName(),
109  obr_,
111  ),
112  alpha,
114  );
115 
116  volScalarField backgroundAlpha
117  (
118  IOobject
119  (
120  alphaName_ + "_background",
121  obr_.time().timeName(),
122  obr_,
124  ),
125  alpha,
127  );
128 
129 
130  // Knock out any cell not in patchRegions
131  forAll(liquidCore, cellI)
132  {
133  label regionI = regions[cellI];
134  if (patchRegions.found(regionI))
135  {
136  backgroundAlpha[cellI] = 0;
137  }
138  else
139  {
140  liquidCore[cellI] = 0;
141 
142  scalar regionVol = regionVolume[regionI];
143  if (regionVol < maxDropletVol)
144  {
145  backgroundAlpha[cellI] = 0;
146  }
147  }
148  }
149  liquidCore.correctBoundaryConditions();
150  backgroundAlpha.correctBoundaryConditions();
151 
152  Info<< " Volume of liquid-core = "
153  << fvc::domainIntegrate(liquidCore).value()
154  << endl;
155  Info<< " Volume of background = "
156  << fvc::domainIntegrate(backgroundAlpha).value()
157  << endl;
158 
159  Info<< "Writing liquid-core field to " << liquidCore.name() << endl;
160  liquidCore.write();
161  Info<< "Writing background field to " << backgroundAlpha.name() << endl;
162  backgroundAlpha.write();
163 }
164 
165 
166 Foam::Map<Foam::label> Foam::regionSizeDistribution::findPatchRegions
167 (
168  const polyMesh& mesh,
169  const regionSplit& regions
170 ) const
171 {
172  // Mark all regions starting at patches
173  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
174 
175  // Count number of patch faces (just for initial sizing)
176  const labelHashSet patchIDs(mesh.boundaryMesh().patchSet(patchNames_));
177 
178  label nPatchFaces = 0;
179  forAllConstIter(labelHashSet, patchIDs, iter)
180  {
181  nPatchFaces += mesh.boundaryMesh()[iter.key()].size();
182  }
183 
184 
185  Map<label> patchRegions(nPatchFaces);
186  forAllConstIter(labelHashSet, patchIDs, iter)
187  {
188  const polyPatch& pp = mesh.boundaryMesh()[iter.key()];
189 
190  // Collect all regions on the patch
191  const labelList& faceCells = pp.faceCells();
192 
193  forAll(faceCells, i)
194  {
195  patchRegions.insert
196  (
197  regions[faceCells[i]],
198  Pstream::myProcNo() // dummy value
199  );
200  }
201  }
202 
203 
204  // Make sure all the processors have the same set of regions
205  Pstream::mapCombineGather(patchRegions, minEqOp<label>());
206  Pstream::mapCombineScatter(patchRegions);
207 
208  return patchRegions;
209 }
210 
211 
212 Foam::tmp<Foam::scalarField> Foam::regionSizeDistribution::divide
213 (
214  const scalarField& num,
215  const scalarField& denom
216 )
217 {
218  tmp<scalarField> tresult(new scalarField(num.size()));
219  scalarField& result = tresult();
220 
221  forAll(denom, i)
222  {
223  if (denom[i] != 0)
224  {
225  result[i] = num[i]/denom[i];
226  }
227  else
228  {
229  result[i] = 0.0;
230  }
231  }
232  return tresult;
233 }
234 
235 
236 void Foam::regionSizeDistribution::writeGraphs
237 (
238  const word& fieldName, // name of field
239  const labelList& indices, // index of bin for each region
240  const scalarField& sortedField, // per region field data
241  const scalarField& binCount, // per bin number of regions
242  const coordSet& coords // graph data for bins
243 ) const
244 {
245  if (Pstream::master())
246  {
247  // Calculate per-bin average
248  scalarField binSum(nBins_, 0.0);
249  forAll(sortedField, i)
250  {
251  binSum[indices[i]] += sortedField[i];
252  }
253 
254  scalarField binAvg(divide(binSum, binCount));
255 
256  // Per bin deviation
257  scalarField binSqrSum(nBins_, 0.0);
258  forAll(sortedField, i)
259  {
260  binSqrSum[indices[i]] += Foam::sqr(sortedField[i]);
261  }
262  scalarField binDev
263  (
264  sqrt(divide(binSqrSum, binCount) - Foam::sqr(binAvg))
265  );
266 
267  // Write average
268  writeGraph(coords, fieldName + "_sum", binSum);
269  // Write average
270  writeGraph(coords, fieldName + "_avg", binAvg);
271  // Write deviation
272  writeGraph(coords, fieldName + "_dev", binDev);
273  }
274 }
275 
276 
277 void Foam::regionSizeDistribution::writeGraphs
278 (
279  const word& fieldName, // name of field
280  const scalarField& cellField, // per cell field data
281  const regionSplit& regions, // per cell the region(=droplet)
282  const labelList& sortedRegions, // valid regions in sorted order
283  const scalarField& sortedNormalisation,
284 
285  const labelList& indices, // per region index of bin
286  const scalarField& binCount, // per bin number of regions
287  const coordSet& coords // graph data for bins
288 ) const
289 {
290  // Sum on a per-region basis. Parallel reduced.
291  Map<scalar> regionField(regionSum(regions, cellField));
292 
293  // Extract in region order
294  scalarField sortedField
295  (
296  sortedNormalisation
297  * extractData
298  (
299  sortedRegions,
300  regionField
301  )
302  );
303 
304  writeGraphs
305  (
306  fieldName, // name of field
307  indices, // index of bin for each region
308  sortedField, // per region field data
309  binCount, // per bin number of regions
310  coords // graph data for bins
311  );
312 }
313 
314 
315 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
316 
317 Foam::regionSizeDistribution::regionSizeDistribution
318 (
319  const word& name,
320  const objectRegistry& obr,
321  const dictionary& dict,
322  const bool loadFromFiles
323 )
324 :
325  functionObjectFile(obr, name, typeName),
326  name_(name),
327  obr_(obr),
328  active_(true),
329  alphaName_(dict.lookup("field")),
330  patchNames_(dict.lookup("patches"))
331 {
332  // Check if the available mesh is an fvMesh, otherwise deactivate
333  if (isA<fvMesh>(obr_))
334  {
335  read(dict);
336  }
337  else
338  {
339  active_ = false;
340  WarningIn
341  (
342  "regionSizeDistribution::regionSizeDistribution"
343  "("
344  "const word&, "
345  "const objectRegistry&, "
346  "const dictionary&, "
347  "const bool"
348  ")"
349  ) << "No fvMesh available, deactivating " << name_ << nl
350  << endl;
351  }
352 }
353 
354 
355 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
356 
358 {}
359 
360 
361 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
362 
364 {
365  if (active_)
366  {
367  dict.lookup("field") >> alphaName_;
368  dict.lookup("patches") >> patchNames_;
369  dict.lookup("threshold") >> threshold_;
370  dict.lookup("maxDiameter") >> maxDiam_;
371  minDiam_ = 0.0;
372  dict.readIfPresent("minDiameter", minDiam_);
373  dict.lookup("nBins") >> nBins_;
374  dict.lookup("fields") >> fields_;
375 
376  word format(dict.lookup("setFormat"));
377  formatterPtr_ = writer<scalar>::New(format);
378 
379  if (dict.found("coordinateSystem"))
380  {
381  coordSysPtr_.reset(new coordinateSystem(obr_, dict));
382 
383  Info<< "Transforming all vectorFields with coordinate system "
384  << coordSysPtr_().name() << endl;
385  }
386  }
387 }
388 
389 
391 {
392  // Do nothing - only valid on write
393 }
394 
395 
397 {
398  // Do nothing - only valid on write
399 }
400 
401 
403 {
404  // Do nothing - only valid on write
405 }
406 
407 
409 {
410  if (active_)
411  {
412  Info<< type() << " " << name_ << " output:" << nl;
413 
414  const fvMesh& mesh = refCast<const fvMesh>(obr_);
415 
417  if (obr_.foundObject<volScalarField>(alphaName_))
418  {
419  Info<< " Looking up field " << alphaName_ << endl;
420  }
421  else
422  {
423  Info<< " Reading field " << alphaName_ << endl;
424  alphaPtr.reset
425  (
426  new volScalarField
427  (
428  IOobject
429  (
430  alphaName_,
431  mesh.time().timeName(),
432  mesh,
435  ),
436  mesh
437  )
438  );
439  }
440 
441 
442  const volScalarField& alpha =
443  (
444  alphaPtr.valid()
445  ? alphaPtr()
446  : obr_.lookupObject<volScalarField>(alphaName_)
447  );
448 
449  Info<< " Volume of alpha = "
450  << fvc::domainIntegrate(alpha).value()
451  << endl;
452 
453  const scalar meshVol = gSum(mesh.V());
454  const scalar maxDropletVol = 1.0/6.0*pow(maxDiam_, 3);
455  const scalar delta = (maxDiam_-minDiam_)/nBins_;
456 
457  Info<< " Mesh volume = " << meshVol << endl;
458  Info<< " Maximum droplet diameter = " << maxDiam_ << endl;
459  Info<< " Maximum droplet volume = " << maxDropletVol << endl;
460 
461 
462  // Determine blocked faces
463  boolList blockedFace(mesh.nFaces(), false);
464  label nBlocked = 0;
465 
466  {
467  for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
468  {
469  scalar ownVal = alpha[mesh.faceOwner()[faceI]];
470  scalar neiVal = alpha[mesh.faceNeighbour()[faceI]];
471 
472  if
473  (
474  (ownVal < threshold_ && neiVal > threshold_)
475  || (ownVal > threshold_ && neiVal < threshold_)
476  )
477  {
478  blockedFace[faceI] = true;
479  nBlocked++;
480  }
481  }
482 
483  // Block coupled faces
484  forAll(alpha.boundaryField(), patchI)
485  {
486  const fvPatchScalarField& fvp = alpha.boundaryField()[patchI];
487  if (fvp.coupled())
488  {
489  tmp<scalarField> townFld(fvp.patchInternalField());
490  const scalarField& ownFld = townFld();
491  tmp<scalarField> tnbrFld(fvp.patchNeighbourField());
492  const scalarField& nbrFld = tnbrFld();
493 
494  label start = fvp.patch().patch().start();
495 
496  forAll(ownFld, i)
497  {
498  scalar ownVal = ownFld[i];
499  scalar neiVal = nbrFld[i];
500 
501  if
502  (
503  (ownVal < threshold_ && neiVal > threshold_)
504  || (ownVal > threshold_ && neiVal < threshold_)
505  )
506  {
507  blockedFace[start+i] = true;
508  nBlocked++;
509  }
510  }
511  }
512  }
513  }
514 
515 
516  regionSplit regions(mesh, blockedFace);
517 
518  Info<< " Determined " << regions.nRegions()
519  << " disconnected regions" << endl;
520 
521 
522  if (debug)
523  {
524  volScalarField region
525  (
526  IOobject
527  (
528  "region",
529  mesh.time().timeName(),
530  mesh,
533  ),
534  mesh,
535  dimensionedScalar("zero", dimless, 0)
536  );
537  Info<< " Dumping region as volScalarField to " << region.name()
538  << endl;
539 
540  forAll(regions, cellI)
541  {
542  region[cellI] = regions[cellI];
543  }
544  region.correctBoundaryConditions();
545  region.write();
546  }
547 
548 
549  // Determine regions connected to supplied patches
550  Map<label> patchRegions(findPatchRegions(mesh, regions));
551 
552 
553 
554  // Sum all regions
555  const scalarField alphaVol(alpha.internalField()*mesh.V());
556  Map<scalar> allRegionVolume(regionSum(regions, mesh.V()));
557  Map<scalar> allRegionAlphaVolume(regionSum(regions, alphaVol));
558  Map<label> allRegionNumCells
559  (
560  regionSum
561  (
562  regions,
563  labelField(mesh.nCells(), 1.0)
564  )
565  );
566 
567  if (debug)
568  {
569  Info<< " " << token::TAB << "Region"
570  << token::TAB << "Volume(mesh)"
571  << token::TAB << "Volume(" << alpha.name() << "):"
572  << token::TAB << "nCells"
573  << endl;
574  scalar meshSumVol = 0.0;
575  scalar alphaSumVol = 0.0;
576  label nCells = 0;
577 
578  Map<scalar>::const_iterator vIter = allRegionVolume.begin();
579  Map<scalar>::const_iterator aIter = allRegionAlphaVolume.begin();
580  Map<label>::const_iterator numIter = allRegionNumCells.begin();
581  for
582  (
583  ;
584  vIter != allRegionVolume.end()
585  && aIter != allRegionAlphaVolume.end();
586  ++vIter, ++aIter, ++numIter
587  )
588  {
589  Info<< " " << token::TAB << vIter.key()
590  << token::TAB << vIter()
591  << token::TAB << aIter()
592  << token::TAB << numIter()
593  << endl;
594 
595  meshSumVol += vIter();
596  alphaSumVol += aIter();
597  nCells += numIter();
598  }
599  Info<< " " << token::TAB << "Total:"
600  << token::TAB << meshSumVol
601  << token::TAB << alphaSumVol
602  << token::TAB << nCells
603  << endl;
604  Info<< endl;
605  }
606 
607 
608 
609 
610  {
611  Info<< " Patch connected regions (liquid core):" << endl;
612  Info<< token::TAB << " Region"
613  << token::TAB << "Volume(mesh)"
614  << token::TAB << "Volume(" << alpha.name() << "):"
615  << endl;
616  forAllConstIter(Map<label>, patchRegions, iter)
617  {
618  label regionI = iter.key();
619  Info<< " " << token::TAB << iter.key()
620  << token::TAB << allRegionVolume[regionI]
621  << token::TAB << allRegionAlphaVolume[regionI] << endl;
622 
623  }
624  Info<< endl;
625  }
626 
627  {
628  Info<< " Background regions:" << endl;
629  Info<< " " << token::TAB << "Region"
630  << token::TAB << "Volume(mesh)"
631  << token::TAB << "Volume(" << alpha.name() << "):"
632  << endl;
633  Map<scalar>::const_iterator vIter = allRegionVolume.begin();
634  Map<scalar>::const_iterator aIter = allRegionAlphaVolume.begin();
635 
636  for
637  (
638  ;
639  vIter != allRegionVolume.end()
640  && aIter != allRegionAlphaVolume.end();
641  ++vIter, ++aIter
642  )
643  {
644  if
645  (
646  !patchRegions.found(vIter.key())
647  && vIter() >= maxDropletVol
648  )
649  {
650  Info<< " " << token::TAB << vIter.key()
651  << token::TAB << vIter()
652  << token::TAB << aIter() << endl;
653  }
654  }
655  Info<< endl;
656  }
657 
658 
659 
660  // Split alpha field
661  // ~~~~~~~~~~~~~~~~~
662  // Split into
663  // - liquidCore : region connected to inlet patches
664  // - per region a volume : for all other regions
665  // - backgroundAlpha : remaining alpha
666  writeAlphaFields(regions, patchRegions, allRegionVolume, alpha);
667 
668 
669  // Extract droplet-only allRegionVolume, i.e. delete liquid core
670  // (patchRegions) and background regions from maps.
671  // Note that we have to use mesh volume (allRegionVolume) and not
672  // allRegionAlphaVolume since background might not have alpha in it.
673  forAllIter(Map<scalar>, allRegionVolume, vIter)
674  {
675  label regionI = vIter.key();
676  if
677  (
678  patchRegions.found(regionI)
679  || vIter() >= maxDropletVol
680  )
681  {
682  allRegionVolume.erase(vIter);
683  allRegionAlphaVolume.erase(regionI);
684  allRegionNumCells.erase(regionI);
685  }
686  }
687 
688  if (allRegionVolume.size())
689  {
690  // Construct mids of bins for plotting
691  pointField xBin(nBins_);
692 
693  scalar x = 0.5*delta;
694  forAll(xBin, i)
695  {
696  xBin[i] = point(x, 0, 0);
697  x += delta;
698  }
699 
700  const coordSet coords("diameter", "x", xBin, mag(xBin));
701 
702 
703  // Get in region order the alpha*volume and diameter
704  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
705 
706  const labelList sortedRegions = allRegionAlphaVolume.sortedToc();
707 
708  scalarField sortedVols
709  (
710  extractData
711  (
712  sortedRegions,
713  allRegionAlphaVolume
714  )
715  );
716 
717  // Calculate the diameters
718  scalarField sortedDiameters(sortedVols.size());
719  forAll(sortedDiameters, i)
720  {
721  sortedDiameters[i] = Foam::cbrt
722  (
723  sortedVols[i]
725  );
726  }
727 
728  // Determine the bin index for all the diameters
729  labelList indices(sortedDiameters.size());
730  forAll(sortedDiameters, i)
731  {
732  indices[i] = (sortedDiameters[i]-minDiam_)/delta;
733  }
734 
735  // Calculate the counts per diameter bin
736  scalarField binCount(nBins_, 0.0);
737  forAll(sortedDiameters, i)
738  {
739  binCount[indices[i]] += 1.0;
740  }
741 
742  // Write counts
743  if (Pstream::master())
744  {
745  writeGraph(coords, "count", binCount);
746  }
747 
748  // Write to screen
749  {
750  Info<< " Bins:" << endl;
751  Info<< " " << token::TAB << "Bin"
752  << token::TAB << "Min diameter"
753  << token::TAB << "Count:"
754  << endl;
755 
756  scalar diam = 0.0;
757  forAll(binCount, binI)
758  {
759  Info<< " " << token::TAB << binI
760  << token::TAB << diam
761  << token::TAB << binCount[binI] << endl;
762  diam += delta;
763  }
764  Info<< endl;
765  }
766 
767 
768  // Write average and deviation of droplet volume.
769  writeGraphs
770  (
771  "volume", // name of field
772  indices, // per region the bin index
773  sortedVols, // per region field data
774  binCount, // per bin number of regions
775  coords // graph data for bins
776  );
777 
778  // Collect some more field
779  {
780  wordList scalarNames(obr_.names(volScalarField::typeName));
781  labelList selected = findStrings(fields_, scalarNames);
782 
783  forAll(selected, i)
784  {
785  const word& fldName = scalarNames[selected[i]];
786  Info<< " Scalar field " << fldName << endl;
787 
788  const scalarField& fld = obr_.lookupObject
789  <
791  >(fldName).internalField();
792 
793  writeGraphs
794  (
795  fldName, // name of field
796  alphaVol*fld, // per cell field data
797 
798  regions, // per cell the region(=droplet)
799  sortedRegions, // valid regions in sorted order
800  1.0/sortedVols, // per region normalisation
801 
802  indices, // index of bin for each region
803  binCount, // per bin number of regions
804  coords // graph data for bins
805  );
806  }
807  }
808  {
809  wordList vectorNames(obr_.names(volVectorField::typeName));
810  labelList selected = findStrings(fields_, vectorNames);
811 
812  forAll(selected, i)
813  {
814  const word& fldName = vectorNames[selected[i]];
815  Info<< " Vector field " << fldName << endl;
816 
817  vectorField fld = obr_.lookupObject
818  <
820  >(fldName).internalField();
821 
822  if (coordSysPtr_.valid())
823  {
824  Info<< "Transforming vector field " << fldName
825  << " with coordinate system "
826  << coordSysPtr_().name()
827  << endl;
828 
829  fld = coordSysPtr_().localVector(fld);
830  }
831 
832 
833  // Components
834 
835  for (direction cmp = 0; cmp < vector::nComponents; cmp++)
836  {
837  writeGraphs
838  (
839  fldName + vector::componentNames[cmp],
840  alphaVol*fld.component(cmp),// per cell field data
841 
842  regions, // per cell the region(=droplet)
843  sortedRegions, // valid regions in sorted order
844  1.0/sortedVols, // per region normalisation
845 
846  indices, // index of bin for each region
847  binCount, // per bin number of regions
848  coords // graph data for bins
849  );
850  }
851 
852  // Magnitude
853  writeGraphs
854  (
855  fldName + "mag", // name of field
856  alphaVol*mag(fld), // per cell field data
857 
858  regions, // per cell the region(=droplet)
859  sortedRegions, // valid regions in sorted order
860  1.0/sortedVols, // per region normalisation
861 
862  indices, // index of bin for each region
863  binCount, // per bin number of regions
864  coords // graph data for bins
865  );
866  }
867  }
868  }
869  }
870 }
871 
872 
873 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
Output to file stream.
Definition: OFstream.H:81
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
static autoPtr< writer > New(const word &writeFormat)
Return a reference to the selected writer.
Definition: writer.C:35
unsigned char direction
Definition: direction.H:43
static const char *const typeName
Definition: Field.H:94
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
virtual void timeSet()
Called when time was set at the end of the Time::operator++.
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
label nFaces() const
vector point
Point is a vector.
Definition: point.H:41
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
bool erase(const iterator &)
Erase a hashedEntry specified by given iterator.
dimensioned< scalar > mag(const dimensioned< Type > &)
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
An STL-conforming const_iterator.
Definition: HashTable.H:470
Holds list of sampling positions.
Definition: coordSet.H:49
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:578
Base class for other coordinate system specifications.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
static const char * componentNames[]
Definition: Vector.H:79
label nRegions() const
Return total number of regions.
Definition: regionSplit.H:218
List< label > sortedToc() const
Return the table of contents as a sorted list.
const labelUList & faceCells() const
Return face-cell addressing.
Definition: polyPatch.C:316
#define forAllIter(Container, container, iter)
Definition: UList.H:440
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
A class for handling words, derived from string.
Definition: word.H:59
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 word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
dimensionedScalar cbrt(const dimensionedScalar &ds)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
InternalField & internalField()
Return internal field.
label nCells() const
messageStream Info
dynamicFvMesh & mesh
A 1D vector of objects of type <T> with a fixed size <Size>.
Definition: FixedList.H:53
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual bool coupled() const
Return true if this patch field is coupled.
Definition: fvPatchField.H:342
Namespace for OpenFOAM.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:39
Type gSum(const FieldField< Field, Type > &f)
virtual bool write() const
Write using setting from DB.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:143
word format(conversionProperties.lookup("format"))
Number of components in this vector space.
Definition: VectorSpace.H:88
virtual tmp< Field< Type > > patchNeighbourField() const
Return patchField on the opposite patch of a coupled patch.
Definition: fvPatchField.H:411
static void mapCombineScatter(const List< commsStruct > &comms, Container &Values, const int tag, const label comm)
Scatter data. Reverse of combineGather.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define WarningIn(functionName)
Report a warning using Foam::Warning.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:215
bool findStrings(const wordReListMatcher &matcher, const std::string &str)
Return true if string matches one of the regular expressions.
Definition: stringListOps.H:52
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:118
#define forAll(list, i)
Definition: UList.H:421
virtual Ostream & write(const token &)=0
Write next token to stream.
scalar delta
virtual void write()
Calculate the regionSizeDistribution and write.
static void mapCombineGather(const List< commsStruct > &comms, Container &Values, const CombineOp &cop, const int tag, const label comm)
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
const word & name() const
Return name.
Definition: IOobject.H:260
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1073
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:404
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool usePatchGroups=true) const
Return the set of patch IDs corresponding to the given names.
virtual void end()
Execute at the final time-loop, currently does nothing.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
label nInternalFaces() const
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
autoPtr< volScalarField > alphaPtr
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
scalar y
Registry of regIOobjects.
A class for handling file names.
Definition: fileName.H:69
virtual void execute()
Execute, currently does nothing.
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
Operations on lists of strings.
Base class for output file data handling.
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:398
bool read(const char *, int32_t &)
Definition: int32IO.C:87
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
iterator begin()
Iterator set to the beginning of the HashTable.
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:420
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1079
const Key & key() const
Return the Key corresponding to the iterator.
Definition: HashTableI.H:262
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:300
Volume integrate volField creating a volField.
Plus op for FixedList<scalar>
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
bool found(const Key &) const
Return true if hashedEntry is found in table.
Definition: HashTable.C:109
A class for managing temporary objects.
Definition: PtrList.H:118
This class separates the mesh into distinct unconnected regions, each of which is then given a label ...
Definition: regionSplit.H:116
conserve internalField()+
defineTypeNameAndDebug(combustionModel, 0)
void divide(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
virtual void read(const dictionary &)
Read the regionSizeDistribution data.