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);
178 endSampledValues_.autoMap(mapper);
182 startSampleTime_ = -1;
197 refCast<const PatchFunction1Types::MappedFile<Type>>(pf1);
199 startSampledValues_.
rmap(tiptf.startSampledValues_, addr);
200 endSampledValues_.rmap(tiptf.endSampledValues_, addr);
204 startSampleTime_ = -1;
215 const polyMesh&
mesh = this->patch_.boundaryMesh().mesh();
245 <<
"Read " << samplePoints.size() <<
" sample points from "
246 << samplePointsFile <<
endl;
253 && mapMethod_ !=
"planarInterpolation"
257 if (this->faceValues())
264 this->localPosition(this->patch_.faceCentres()),
277 this->localPosition(this->patch_.localPoints()),
286 const fileName samplePointsDir = samplePointsFile.path();
287 sampleTimes_ = Time::findTimes(samplePointsDir);
291 << samplePointsDir <<
" found times "
292 << pointToPointPlanarInterpolation::timeNames(sampleTimes_)
301 bool foundTime = mapperPtr_().findTime
313 <<
"Cannot find starting sampling values for index "
315 <<
"Have sampling values for "
316 << pointToPointPlanarInterpolation::timeNames(sampleTimes_) <<
nl
318 << time.
constant()/
"boundaryData"/this->patch_.name()
319 <<
"\n on patch " << this->patch_.name()
320 <<
" of field " << fieldTableName_
327 if (lo != startSampleTime_)
329 startSampleTime_ = lo;
331 if (startSampleTime_ == endSampleTime_)
336 Pout<<
"checkTable : Setting startValues to (already read) "
339 /sampleTimes_[startSampleTime_].
name()
342 startSampledValues_ = endSampledValues_;
343 startAverage_ = endAverage_;
349 Pout<<
"checkTable : Reading startValues from "
352 /sampleTimes_[lo].
name()
358 const fileName valsFile
364 /sampleTimes_[startSampleTime_].name()
378 const rawIOField<Type> vals(io, setAverage_);
381 startAverage_ = vals.average();
384 if (vals.size() != mapperPtr_().sourceSize())
387 <<
"Number of values (" << vals.size()
388 <<
") differs from the number of points ("
389 << mapperPtr_().sourceSize()
393 startSampledValues_ = mapperPtr_().interpolate(vals);
397 if (hi != endSampleTime_)
401 if (endSampleTime_ == -1)
406 Pout<<
"checkTable : Clearing endValues" <<
endl;
408 endSampledValues_.clear();
414 Pout<<
"checkTable : Reading endValues from "
417 /sampleTimes_[endSampleTime_].
name()
428 /sampleTimes_[endSampleTime_].name()
442 const rawIOField<Type> vals(io, setAverage_);
445 endAverage_ = vals.average();
448 if (vals.size() != mapperPtr_().sourceSize())
451 <<
"Number of values (" << vals.size()
452 <<
") differs from the number of points ("
453 << mapperPtr_().sourceSize()
457 endSampledValues_ = mapperPtr_().interpolate(vals);
473 auto&
fld = tfld.ref();
476 if (endSampleTime_ == -1)
481 Pout<<
"MappedFile<Type>::value : Sampled, non-interpolated values"
482 <<
" from start time:"
483 << sampleTimes_[startSampleTime_].
name() <<
nl;
486 fld = startSampledValues_;
487 wantedAverage = startAverage_;
491 scalar start = sampleTimes_[startSampleTime_].value();
492 scalar
end = sampleTimes_[endSampleTime_].value();
494 scalar
s = (
x - start)/(
end - start);
498 Pout<<
"MappedFile<Type>::value : Sampled, interpolated values"
499 <<
" between start time:"
500 << sampleTimes_[startSampleTime_].
name()
501 <<
" and end time:" << sampleTimes_[endSampleTime_].
name()
502 <<
" with weight:" <<
s <<
endl;
505 fld = ((1 -
s)*startSampledValues_ +
s*endSampledValues_);
506 wantedAverage = (1 -
s)*startAverage_ +
s*endAverage_;
514 if (this->faceValues())
526 Pout<<
"MappedFile<Type>::value :"
527 <<
" actual average:" << averagePsi
528 <<
" wanted average:" << wantedAverage
532 if (
mag(averagePsi) < VSMALL)
535 const Type offset = wantedAverage - averagePsi;
538 Pout<<
"MappedFile<Type>::value :"
539 <<
" offsetting with:" << offset <<
endl;
545 const scalar scale =
mag(wantedAverage)/
mag(averagePsi);
549 Pout<<
"MappedFile<Type>::value :"
550 <<
" scaling with:" << scale <<
endl;
559 fld += offset_->value(
x);
564 Pout<<
"MappedFile<Type>::value : set fixedValue to min:" <<
gMin(
fld)
596 if (dictConstructed_)
620 "planarInterpolation",
626 offset_->writeData(os);
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.
word name(const complex &c)
Return string representation of complex.
virtual const fileName & name() const
Return 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,...
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 file.
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.
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
fileName globalPath() const
Return global path for the case.
MappedFile(const polyPatch &pp, const word &redirectType, const word &entryName, const dictionary &dict, const bool faceValues=true)
Construct from entry name and dictionary.
#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.