34 namespace functionObjects
41 regionSizeDistribution,
47 template<
class T,
unsigned Size>
68 void Foam::functionObjects::regionSizeDistribution::writeGraph
71 const word& valueName,
75 const wordList valNames(1, valueName);
80 OFstream str(outputPath/formatterPtr_().getFileName(coords, valNames));
82 Info<<
" Writing distribution of " << valueName <<
" to " << str.
name()
87 formatterPtr_().write(coords, valNames, valPtrs, str);
91 void Foam::functionObjects::regionSizeDistribution::writeAlphaFields
99 const scalar maxDropletVol = 1.0/6.0*
pow(maxDiam_, 3);
114 alphaName_ +
"_liquidCore",
115 obr_.time().timeName(),
127 alphaName_ +
"_background",
128 obr_.time().timeName(),
140 label regionI = regions[celli];
141 if (patchRegions.
found(regionI))
143 backgroundAlpha[celli] = 0;
147 liquidCore[celli] = 0;
149 scalar regionVol = regionVolume[regionI];
150 if (regionVol < maxDropletVol)
152 backgroundAlpha[celli] = 0;
156 liquidCore.correctBoundaryConditions();
157 backgroundAlpha.correctBoundaryConditions();
159 Info<<
" Volume of liquid-core = " 162 Info<<
" Volume of background = " 166 Info<<
" Writing liquid-core field to " << liquidCore.name() <<
endl;
168 Info<<
" Writing background field to " << backgroundAlpha.name() <<
endl;
169 backgroundAlpha.
write();
174 Foam::functionObjects::regionSizeDistribution::findPatchRegions
186 label nPatchFaces = 0;
205 regions[faceCells[i]],
221 Foam::functionObjects::regionSizeDistribution::divide
234 result[i] = num[i]/denom[i];
245 void Foam::functionObjects::regionSizeDistribution::writeGraphs
247 const word& fieldName,
260 binSum[indices[i]] += sortedField[i];
269 binSqrSum[indices[i]] +=
Foam::sqr(sortedField[i]);
277 writeGraph(coords, fieldName +
"_sum", binSum);
279 writeGraph(coords, fieldName +
"_avg", binAvg);
281 writeGraph(coords, fieldName +
"_dev", binDev);
286 void Foam::functionObjects::regionSizeDistribution::writeGraphs
288 const word& fieldName,
300 Map<scalar> regionField(regionSum(regions, cellField));
326 Foam::functionObjects::regionSizeDistribution::regionSizeDistribution
334 alphaName_(dict.
lookup(
"field")),
335 patchNames_(dict.
lookup(
"patches"))
337 if (!isA<fvMesh>(obr_))
357 dict.
lookup(
"field") >> alphaName_;
358 dict.
lookup(
"patches") >> patchNames_;
359 dict.
lookup(
"threshold") >> threshold_;
360 dict.
lookup(
"maxDiameter") >> maxDiam_;
363 dict.
lookup(
"nBins") >> nBins_;
364 dict.
lookup(
"fields") >> fields_;
369 if (dict.
found(
"coordinateSystem"))
373 Info<<
"Transforming all vectorFields with coordinate system " 374 << coordSysPtr_().name() <<
endl;
391 const fvMesh& mesh = refCast<const fvMesh>(obr_);
396 Info<<
" Looking up field " << alphaName_ <<
endl;
400 Info<<
" Reading field " << alphaName_ <<
endl;
426 Info<<
" Volume of alpha = " 430 const scalar meshVol =
gSum(mesh.
V());
431 const scalar maxDropletVol = 1.0/6.0*
pow(maxDiam_, 3);
432 const scalar
delta = (maxDiam_-minDiam_)/nBins_;
434 Info<<
" Mesh volume = " << meshVol <<
endl;
435 Info<<
" Maximum droplet diameter = " << maxDiam_ <<
endl;
436 Info<<
" Maximum droplet volume = " << maxDropletVol <<
endl;
446 scalar ownVal = alpha[mesh.
faceOwner()[facei]];
451 (ownVal < threshold_ && neiVal > threshold_)
452 || (ownVal > threshold_ && neiVal < threshold_)
455 blockedFace[facei] =
true;
475 scalar ownVal = ownFld[i];
476 scalar neiVal = nbrFld[i];
480 (ownVal < threshold_ && neiVal > threshold_)
481 || (ownVal > threshold_ && neiVal < threshold_)
484 blockedFace[start+i] =
true;
496 <<
" disconnected regions" <<
endl;
514 Info<<
" Dumping region as volScalarField to " << region.name()
519 region[celli] = regions[celli];
521 region.correctBoundaryConditions();
527 Map<label> patchRegions(findPatchRegions(mesh, regions));
533 Map<scalar> allRegionVolume(regionSum(regions, mesh.
V()));
534 Map<scalar> allRegionAlphaVolume(regionSum(regions, alphaVol));
551 scalar meshSumVol = 0.0;
552 scalar alphaSumVol = 0.0;
561 vIter != allRegionVolume.end()
562 && aIter != allRegionAlphaVolume.
end();
563 ++vIter, ++aIter, ++numIter
572 meshSumVol += vIter();
573 alphaSumVol += aIter();
588 Info<<
" Patch connected regions (liquid core):" <<
endl;
595 label regionI = iter.key();
605 Info<<
" Background regions:" <<
endl;
616 vIter != allRegionVolume.end()
617 && aIter != allRegionAlphaVolume.
end();
624 && vIter() >= maxDropletVol
643 writeAlphaFields(regions, patchRegions, allRegionVolume, alpha);
652 label regionI = vIter.key();
655 patchRegions.
found(regionI)
656 || vIter() >= maxDropletVol
659 allRegionVolume.erase(vIter);
660 allRegionAlphaVolume.
erase(regionI);
661 allRegionNumCells.
erase(regionI);
665 if (allRegionVolume.size())
673 xBin[i] =
point(x, 0, 0);
677 const coordSet coords(
"diameter",
"x", xBin,
mag(xBin));
696 forAll(sortedDiameters, i)
706 labelList indices(sortedDiameters.size());
707 forAll(sortedDiameters, i)
709 indices[i] = (sortedDiameters[i]-minDiam_)/delta;
714 forAll(sortedDiameters, i)
716 binCount[indices[i]] += 1.0;
722 writeGraph(coords,
"count", binCount);
762 const word& fldName = scalarNames[selected[i]];
763 Info<<
" Scalar field " << fldName <<
endl;
768 >(fldName).primitiveField();
791 const word& fldName = vectorNames[selected[i]];
792 Info<<
" Vector field " << fldName <<
endl;
797 >(fldName).primitiveField();
799 if (coordSysPtr_.valid())
801 Info<<
"Transforming vector field " << fldName
802 <<
" with coordinate system " 803 << coordSysPtr_().name()
806 fld = coordSysPtr_().localVector(fld);
void divide(FieldField< Field, Type > &f, const FieldField< Field, Type > &f1, const FieldField< Field, scalar > &f2)
const fvPatch & patch() const
Return patch.
Base class for other coordinate system specifications.
This class separates the mesh into distinct unconnected regions, each of which is then given a label ...
Field< label > labelField
Specialisation of Field<T> for label.
#define forAll(list, i)
Loop across all elements in list.
virtual bool write()
Calculate the regionSizeDistribution and write.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
A class for handling file names.
virtual bool coupled() const
Return true if this patch field is coupled.
An STL-conforming const_iterator.
virtual bool execute()
Do nothing.
errorManipArg< error, int > exit(error &err, const int errNo=1)
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
A 1D vector of objects of type <T> with a fixed size <Size>.
A list of keyword definitions, which are a keyword followed by any number of values (e...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
static const char *const typeName
dimensionedSymmTensor sqr(const dimensionedVector &dv)
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
void size(const label)
Override size to be inconsistent with allocated storage.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
static const char *const componentNames[]
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Plus op for FixedList<scalar>
static bool master(const label communicator=0)
Am I the master process.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
word format(conversionProperties.lookup("format"))
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
const polyPatch & patch() const
Return the polyPatch.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Macros for easy insertion into run-time selection tables.
bool erase(const iterator &)
Erase a hashedEntry specified by given iterator.
bool findStrings(const wordReListMatcher &matcher, const std::string &str)
Return true if string matches one of the regular expressions.
static void mapCombineScatter(const List< commsStruct > &comms, Container &Values, const int tag, const label comm)
Scatter data. Reverse of combineGather.
virtual tmp< Field< Type > > patchNeighbourField() const
Return patchField on the opposite patch of a coupled patch.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
label nRegions() const
Return total number of regions.
bool read(const char *, int32_t &)
static const direction nComponents
Number of components in this vector space.
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))
Holds list of sampling positions.
Type gSum(const FieldField< Field, Type > &f)
label start() const
Return start label of this patch in the polyMesh face list.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
A class for handling words, derived from string.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar cbrt(const dimensionedScalar &ds)
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.
Volume integrate volField creating a volField.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
static void mapCombineGather(const List< commsStruct > &comms, Container &Values, const CombineOp &cop, const int tag, const label comm)
iterator begin()
Iterator set to the beginning of the HashTable.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
virtual bool read(const dictionary &)
Read the regionSizeDistribution data.
bool found(const Key &) const
Return true if hashedEntry is found in table.
static autoPtr< writer > New(const word &writeFormat)
Return a reference to the selected writer.
virtual ~regionSizeDistribution()
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
word name(const complex &)
Return a string representation of a complex.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual const labelList & faceNeighbour() const
Return face neighbour.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
vector point
Point is a vector.
const Key & key() const
Return the Key corresponding to the iterator.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
const labelUList & faceCells() const
Return face-cell addressing.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
addToRunTimeSelectionTable(functionObject, blendingFactor, dictionary)
List< label > sortedToc() const
Return the table of contents as a sorted list.
Mesh data needed to do the Finite Volume discretisation.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
virtual bool write() const
Write using setting from DB.
autoPtr< volScalarField > alphaPtr
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
virtual Ostream & write(const token &)=0
Write next token to stream.
Mesh consisting of general polyhedral cells.
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
A class for managing temporary objects.
virtual const labelList & faceOwner() const
Return face owner.
A patch is a list of labels that address the faces in the global face list.
T & ref() const
Return non-const reference or generate a fatal error.
functionObject base class for writing single files
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
label nInternalFaces() const
const fileName & name() const
Return the name of the stream.
const word & name() const
Return name.
const Time & time() const
Return the top-level database.
label size() const
Return the number of elements in the UPtrList.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.