37 template<
class GeoField>
39 Foam::sampledSets::getOrLoadField(
const word& fieldName)
const 69 void Foam::sampledSets::writeCoordSet
99 void Foam::sampledSets::performAction
101 const VolumeField<Type>&
fld,
105 const word& fieldName =
fld.name();
106 const scalar timeValue =
fld.time().timeOutputValue();
109 autoPtr<interpolation<Type>> interpPtr;
111 if (!samplePointScheme_.empty() && samplePointScheme_ !=
"cell")
117 OFstream* osptr =
nullptr;
119 if (writeAsProbes_ && (request & ACTION_WRITE))
121 osptr = createProbeFile(fieldName);
125 (*osptr) <<
setw(width) << timeValue;
130 Type avgEnsemble =
Zero;
131 label sizeEnsemble = 0;
132 MinMax<Type> limitsEnsemble;
136 const sampledSet&
s = (*this)[seti];
137 const globalIndex& globIdx = globalIndices_[seti];
138 const labelList& globOrder = gatheredSorting_[seti];
140 const word& setName =
s.name();
148 const label celli =
s.cells()[samplei];
149 const label facei =
s.faces()[samplei];
151 if (celli == -1 && facei == -1)
158 values[samplei] = interpPtr().interpolate(
p, celli, facei);
166 const label celli =
s.cells()[samplei];
180 globIdx.gatherInplace(
values);
183 Type avgValue =
Zero;
190 sizeValue =
values.size();
191 limits = MinMax<Type>(
values);
194 avgEnsemble += avgValue;
195 sizeEnsemble += sizeValue;
196 limitsEnsemble += limits;
200 avgValue /= sizeValue;
211 const word resultArg(
'(' + setName +
',' + fieldName +
')');
213 this->setResult(
"average" + resultArg, avgValue);
214 this->setResult(
"min" + resultArg, limits.min());
215 this->setResult(
"max" + resultArg, limits.max());
216 this->setResult(
"size" + resultArg, sizeValue);
220 Info<<
name() <<
' ' << setName <<
" : " << fieldName <<
nl 221 <<
" avg: " << avgValue <<
nl 222 <<
" min: " << limits.min() <<
nl 223 <<
" max: " << limits.max() <<
nl <<
nl;
226 if ((request & ACTION_WRITE) != 0)
232 for (
const Type& val :
values)
234 (*osptr) <<
' ' <<
setw(width) << val;
240 writeCoordSet<Type>(writers_[seti],
values, fieldName);
253 avgEnsemble /= sizeEnsemble;
268 const word resultArg(
'(' + fieldName +
')');
270 this->setResult(
"average" + resultArg, avgEnsemble);
271 this->setResult(
"min" + resultArg, limitsEnsemble.min());
272 this->setResult(
"max" + resultArg, limitsEnsemble.max());
273 this->setResult(
"size" + resultArg, sizeEnsemble);
278 template<
class GeoField>
279 void Foam::sampledSets::performAction
281 const IOobjectList& objects,
288 fieldNames = objects.sortedNames<GeoField>(fieldSelection_);
292 fieldNames = mesh_.thisDb().sortedNames<GeoField>(fieldSelection_);
295 for (
const word& fieldName : fieldNames)
297 tmp<GeoField> tfield = getOrLoadField<GeoField>(fieldName);
301 performAction<typename GeoField::value_type>(tfield(), request);
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
A class for handling file names.
vtk::lineWriter writer(edgeCentres, edgeList::null(), fileName(aMesh.time().globalPath()/"finiteArea-edgesCentres"))
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
static autoPtr< interpolation< Type > > New(const word &interpolationType, const GeometricField< Type, fvPatchField, volMesh > &psi)
Return a reference to the specified interpolation scheme.
constexpr char nl
The newline '\n' character (0x0a)
Ostream & endl(Ostream &os)
Add newline and flush stream.
const T & cref() const
Return const reference to the object or to the contents of a (non-null) managed pointer.
static unsigned int defaultPrecision() noexcept
Return the default precision.
entry * add(entry *entryPtr, bool mergeEntry=false)
Add a new entry.
Ignore writing from objectRegistry::writeObject()
const Time & time() const
Return the top-level database.
static label worldComm
Communicator for all ranks. May differ from commGlobal() if local worlds are in use.
const Type * cfindObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
void write(const word &fieldName, const UList< Type > &field)
Write primitive field of CellData (Poly or Line) or PointData values.
static void broadcast(Type &value, const label comm=UPstream::worldComm)
Broadcast content (contiguous or non-contiguous) to all communicator ranks. Does nothing in non-paral...
IOdictionary propsDict(dictIO)
virtual const objectRegistry & thisDb() const
Return the object registry - resolve conflict polyMesh/lduMesh.
#define forAll(list, i)
Loop across all elements in list.
List< T > values(const HashTable< T, Key, Hash > &tbl, const bool doSort=false)
List of values from HashTable, optionally sorted.
word outputName("finiteArea-edges.obj")
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
Generic templated field type.
A class for handling words, derived from Foam::string.
Base class for writing coordSet(s) and tracks with fields.
static word timeName(const scalar t, const int precision=precision_)
Return a time name for the given scalar time value formatted with the given precision.
static void broadcasts(const label comm, Type &arg1, Args &&... args)
Broadcast multiple items to all communicator ranks. Does nothing in non-parallel. ...
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
List< word > wordList
List of word.
vector point
Point is a vector.
static bool master(const label communicator=worldComm)
True if process corresponds to the master rank in the communicator.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Omanip< int > setw(const int i)
List< label > labelList
A List of labels.
A class for managing temporary objects.
T & emplace(Args &&... args)
Reset with emplace construction. Return reference to the new content.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Defines the attributes of an object for which implicit objectRegistry management is supported...
const fvMesh & mesh_
Reference to the fvMesh.
Do not request registration (bool: false)
static constexpr const zero Zero
Global zero (0)