Go to the documentation of this file.
37 const word& redirectType,
38 const word& entryName,
44 dictConstructed_(
true),
45 setAverage_(
dict.getOrDefault(
"setAverage",
false)),
46 fieldTableName_(entryName),
47 perturb_(
dict.getOrDefault<scalar>(
"perturb", 1
e-5)),
48 pointsName_(
dict.getOrDefault<
word>(
"points",
"points")),
60 startSampledValues_(0),
69 mapMethod_ !=
"planarInterpolation"
70 && mapMethod_ !=
"nearest"
74 <<
"mapMethod should be one of 'planarInterpolation'"
78 dict.readIfPresent(
"fieldTable", fieldTableName_);
86 const word& entryName,
88 const word& fieldTableName,
93 dictConstructed_(
false),
94 setAverage_(
dict.getOrDefault(
"setAverage",
false)),
95 fieldTableName_(fieldTableName),
96 perturb_(
dict.getOrDefault<scalar>(
"perturb", 1
e-5)),
97 pointsName_(
dict.getOrDefault<
word>(
"points",
"points")),
103 "planarInterpolation"
108 startSampleTime_(-1),
109 startSampledValues_(0),
112 endSampledValues_(0),
118 mapMethod_ !=
"planarInterpolation"
119 && mapMethod_ !=
"nearest"
123 <<
"mapMethod should be one of 'planarInterpolation'"
147 dictConstructed_(rhs.dictConstructed_),
148 setAverage_(rhs.setAverage_),
149 fieldTableName_(rhs.fieldTableName_),
150 perturb_(rhs.perturb_),
151 pointsName_(rhs.pointsName_),
152 mapMethod_(rhs.mapMethod_),
153 mapperPtr_(rhs.mapperPtr_.clone()),
154 sampleTimes_(rhs.sampleTimes_),
155 startSampleTime_(rhs.startSampleTime_),
156 startSampledValues_(rhs.startSampledValues_),
157 startAverage_(rhs.startAverage_),
158 endSampleTime_(rhs.endSampleTime_),
159 endSampledValues_(rhs.endSampledValues_),
160 endAverage_(rhs.endAverage_),
161 offset_(rhs.offset_.clone())
175 if (startSampledValues_.size())
177 startSampledValues_.autoMap(mapper);
180 if (endSampledValues_.size())
182 endSampledValues_.autoMap(mapper);
187 startSampleTime_ = -1;
202 refCast<const PatchFunction1Types::MappedFile<Type>>(pf1);
204 if (tiptf.startSampledValues_.size())
206 startSampledValues_.setSize(this->size());
207 startSampledValues_.
rmap(tiptf.startSampledValues_, addr);
210 if (tiptf.endSampledValues_.size())
212 endSampledValues_.setSize(this->size());
213 endSampledValues_.rmap(tiptf.endSampledValues_, addr);
218 startSampleTime_ = -1;
229 const polyMesh&
mesh = this->patch_.boundaryMesh().mesh();
259 <<
"Read " << samplePoints.size() <<
" sample points from "
260 << samplePointsFile <<
endl;
267 && mapMethod_ !=
"planarInterpolation"
271 if (this->faceValues())
278 this->localPosition(this->patch_.faceCentres()),
291 this->localPosition(this->patch_.localPoints()),
300 const fileName samplePointsDir = samplePointsFile.path();
301 sampleTimes_ = Time::findTimes(samplePointsDir);
305 << samplePointsDir <<
" found times "
306 << pointToPointPlanarInterpolation::timeNames(sampleTimes_)
315 bool foundTime = mapperPtr_().findTime
327 <<
"Cannot find starting sampling values for index "
329 <<
"Have sampling values for "
330 << pointToPointPlanarInterpolation::timeNames(sampleTimes_) <<
nl
332 << time.
constant()/
"boundaryData"/this->patch_.name()
333 <<
"\n on patch " << this->patch_.name()
334 <<
" of field " << fieldTableName_
341 if (lo != startSampleTime_)
343 startSampleTime_ = lo;
345 if (startSampleTime_ == endSampleTime_)
350 Pout<<
"checkTable : Setting startValues to (already read) "
353 /sampleTimes_[startSampleTime_].
name()
356 startSampledValues_ = endSampledValues_;
357 startAverage_ = endAverage_;
363 Pout<<
"checkTable : Reading startValues from "
366 /sampleTimes_[lo].
name()
372 const fileName valsFile
378 /sampleTimes_[startSampleTime_].name()
392 const rawIOField<Type> vals(io, setAverage_);
395 startAverage_ = vals.average();
398 if (vals.size() != mapperPtr_().sourceSize())
401 <<
"Number of values (" << vals.size()
402 <<
") differs from the number of points ("
403 << mapperPtr_().sourceSize()
407 startSampledValues_ = mapperPtr_().interpolate(vals);
411 if (hi != endSampleTime_)
415 if (endSampleTime_ == -1)
420 Pout<<
"checkTable : Clearing endValues" <<
endl;
422 endSampledValues_.clear();
428 Pout<<
"checkTable : Reading endValues from "
431 /sampleTimes_[endSampleTime_].
name()
442 /sampleTimes_[endSampleTime_].name()
456 const rawIOField<Type> vals(io, setAverage_);
459 endAverage_ = vals.average();
462 if (vals.size() != mapperPtr_().sourceSize())
465 <<
"Number of values (" << vals.size()
466 <<
") differs from the number of points ("
467 << mapperPtr_().sourceSize()
471 endSampledValues_ = mapperPtr_().interpolate(vals);
487 auto&
fld = tfld.ref();
490 if (endSampleTime_ == -1)
495 Pout<<
"MappedFile<Type>::value : Sampled, non-interpolated values"
496 <<
" from start time:"
497 << sampleTimes_[startSampleTime_].
name() <<
nl;
500 fld = startSampledValues_;
501 wantedAverage = startAverage_;
505 scalar start = sampleTimes_[startSampleTime_].value();
506 scalar
end = sampleTimes_[endSampleTime_].value();
508 scalar
s = (
x - start)/(
end - start);
512 Pout<<
"MappedFile<Type>::value : Sampled, interpolated values"
513 <<
" between start time:"
514 << sampleTimes_[startSampleTime_].
name()
515 <<
" and end time:" << sampleTimes_[endSampleTime_].
name()
516 <<
" with weight:" <<
s <<
endl;
519 fld = ((1 -
s)*startSampledValues_ +
s*endSampledValues_);
520 wantedAverage = (1 -
s)*startAverage_ +
s*endAverage_;
528 if (this->faceValues())
540 Pout<<
"MappedFile<Type>::value :"
541 <<
" actual average:" << averagePsi
542 <<
" wanted average:" << wantedAverage
546 if (
mag(averagePsi) < VSMALL)
549 const Type offset = wantedAverage - averagePsi;
552 Pout<<
"MappedFile<Type>::value :"
553 <<
" offsetting with:" << offset <<
endl;
559 const scalar scale =
mag(wantedAverage)/
mag(averagePsi);
563 Pout<<
"MappedFile<Type>::value :"
564 <<
" scaling with:" << scale <<
endl;
573 fld += offset_->value(
x);
578 Pout<<
"MappedFile<Type>::value : set fixedValue to min:" <<
gMin(
fld)
608 os.writeEntry(
"setAverage", setAverage_);
611 os.writeEntryIfDifferent<scalar>(
"perturb", 1
e-5, perturb_);
613 os.writeEntryIfDifferent<
word>(
"points",
"points", pointsName_);
615 os.writeEntryIfDifferent<
word>
618 "planarInterpolation",
624 offset_->writeData(
os);
639 if (dictConstructed_)
int debug
Static debugging option.
virtual void autoMap(const FieldMapper &mapper)
Map (and resize as needed) from self given a mapping object.
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
A class for handling file names.
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))
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Type gAverage(const FieldField< Field, Type > &f)
virtual void writeData(Ostream &os) const
Write in dictionary format.
static std::string name(const std::string &str)
Return basename (part beyond last /), including its extension.
const bool writeData(pdfDictionary.get< bool >("writeData"))
Abstract base class to hold the Field mapping addressing and weights.
Interpolates between two sets of unstructured points using 2D Delaunay triangulation....
Ostream & endl(Ostream &os)
Add newline and flush stream.
Type gSum(const FieldField< Field, Type > &f)
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Mesh consisting of general polyhedral cells.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Like IOField but falls back to raw IFstream if no header found. Optionally reads average value....
A patch is a list of labels that address the faces in the global face list.
virtual const fileName & name() const
Get the name of the stream.
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))
virtual tmp< Field< Type > > integrate(const scalar x1, const scalar x2) const
Integrate between two values.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
OBJstream os(runTime.globalPath()/outputName)
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Patch value mapping from a set of values stored in a file and a set of unstructured points using the ...
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
virtual void rmap(const PatchFunction1< Type > &pf1, const labelList &addr)
Reverse map the given PatchFunction1 onto this PatchFunction1.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
fileName globalPath() const
Return global path for the case.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
MappedFile(const polyPatch &pp, const word &redirectType, const word &entryName, const dictionary &dict, const bool faceValues=true)
Construct from entry name and dictionary.
void writeEntries(Ostream &os) const
Write coefficient entries in dictionary format.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Type gMin(const FieldField< Field, Type > &f)
const word & constant() const
Return constant name.
Type gMax(const FieldField< Field, Type > &f)
virtual tmp< Field< Type > > value(const scalar) const
Return MappedFile value.