40 const fileName fieldFileName
42 dataDir_/timeName/sampleName_/fieldTableName_
45 const fileName typeFieldFileName
47 dataDir_/timeName/sampleName_
48 /pTraits<Type>::typeName + Field<Type>::typeName
56 else if (
exists(typeFieldFileName))
58 return typeFieldFileName;
63 <<
"Cannot find field file " 64 << fieldFileName <<
" " << typeFieldFileName
67 return fileName::null;
76 if (mapperPtr_.empty())
79 const fileName samplePointsFile(dataDir_/pointsName_);
81 pointField samplePoints((IFstream(samplePointsFile)()));
85 Info<<
"timeVaryingMappedFvPatchField :" 86 <<
" Read " << samplePoints.size() <<
" sample points from " 87 << samplePointsFile <<
endl;
92 const bool nearestOnly
95 && mapMethod_ !=
"planarInterpolation" 101 new pointToPointPlanarInterpolation
104 patch_.patch().faceCentres(),
111 sampleTimes_ = Time::findTimes(dataDir_);
115 Info<<
"timeVaryingMappedFvPatchField : In directory " 116 << dataDir_ <<
" found times " 117 << pointToPointPlanarInterpolation::timeNames(sampleTimes_)
127 bool foundTime = mapperPtr_().findTime
139 <<
"Cannot find starting sampling values for current time " 140 << time().value() <<
nl 141 <<
"Have sampling values for times " 142 << pointToPointPlanarInterpolation::timeNames(sampleTimes_) <<
nl 143 <<
"In directory " << dataDir_ <<
" of field " << fieldTableName_
150 if (lo != startSampleTime_)
152 startSampleTime_ = lo;
154 if (startSampleTime_ == endSampleTime_)
159 Pout<<
"checkTable : Setting startValues to (already read) " 160 << dataDir_/sampleTimes_[startSampleTime_].
name()
163 startSampledValues_ = endSampledValues_;
164 startAverage_ = endAverage_;
170 Pout<<
"checkTable : Reading startValues from " 171 << dataDir_/sampleTimes_[lo].
name()
176 const fileName valsFile
178 findFieldFile(sampleTimes_[startSampleTime_].
name())
185 AverageField<Type> avals((IFstream(valsFile)()));
187 startAverage_ = avals.average();
191 IFstream(valsFile)() >> vals;
194 if (vals.size() != mapperPtr_().sourceSize())
197 <<
"Number of values (" << vals.size()
198 <<
") differs from the number of points (" 199 << mapperPtr_().sourceSize()
203 startSampledValues_ = mapperPtr_().interpolate(vals);
207 if (hi != endSampleTime_)
211 if (endSampleTime_ == -1)
216 Pout<<
"checkTable : Clearing endValues" <<
endl;
218 endSampledValues_.clear();
224 Pout<<
"checkTable : Reading endValues from " 225 << dataDir_/sampleTimes_[endSampleTime_].
name()
230 const fileName valsFile
232 findFieldFile(sampleTimes_[endSampleTime_].
name())
239 AverageField<Type> avals((IFstream(valsFile)()));
241 endAverage_ = avals.average();
245 IFstream(valsFile)() >> vals;
248 if (vals.size() != mapperPtr_().sourceSize())
251 <<
"Number of values (" << vals.size()
252 <<
") differs from the number of points (" 253 << mapperPtr_().sourceSize()
257 endSampledValues_ = mapperPtr_().interpolate(vals);
269 const word& fieldName
273 fieldTableName_(word::null),
274 dataDir_(time().constant()/
"boundaryData"/p.
name()),
275 pointsName_(
"points"),
276 sampleName_(word::null),
281 startSampleTime_(-1),
282 startSampledValues_(0),
285 endSampledValues_(0),
296 const word& fieldName
306 time().constant()/
"boundaryData"/p.
name()
318 "planarInterpolation" 323 startSampleTime_(-1),
324 startSampledValues_(0),
327 endSampledValues_(0),
332 pointsName_.expand();
333 sampleName_.expand();
335 if (dict.
found(
"offset"))
342 mapMethod_ !=
"planarInterpolation" 343 && mapMethod_ !=
"nearest" 349 ) <<
"mapMethod should be one of 'planarInterpolation'" 363 fieldTableName_(ptf.fieldTableName_),
364 dataDir_(ptf.dataDir_),
365 pointsName_(ptf.pointsName_),
366 sampleName_(ptf.sampleName_),
367 setAverage_(ptf.setAverage_),
368 perturb_(ptf.perturb_),
369 mapMethod_(ptf.mapMethod_),
371 sampleTimes_(ptf.sampleTimes_),
372 startSampleTime_(ptf.startSampleTime_),
373 startSampledValues_(ptf.startSampledValues_),
374 startAverage_(ptf.startAverage_),
375 endSampleTime_(ptf.endSampleTime_),
376 endSampledValues_(ptf.endSampledValues_),
377 endAverage_(ptf.endAverage_),
378 offset_(ptf.offset_,
false)
390 if (startSampledValues_.size())
392 m(startSampledValues_, startSampledValues_);
393 m(endSampledValues_, endSampledValues_);
397 startSampleTime_ = -1;
409 startSampledValues_.rmap(tiptf.startSampledValues_, addr);
410 endSampledValues_.rmap(tiptf.endSampledValues_, addr);
414 startSampleTime_ = -1;
431 if (endSampleTime_ == -1)
436 Pout<<
"updateCoeffs : Sampled, non-interpolated values" 437 <<
" from start time:" 438 << sampleTimes_[startSampleTime_].
name() <<
nl;
441 fld = startSampledValues_;
442 wantedAverage = startAverage_;
446 const scalar start = sampleTimes_[startSampleTime_].value();
447 const scalar end = sampleTimes_[endSampleTime_].value();
449 const scalar
s = (time().value() - start)/(end - start);
453 Pout<<
"updateCoeffs : Sampled, interpolated values" 454 <<
" between start time:" 455 << sampleTimes_[startSampleTime_].
name()
456 <<
" and end time:" << sampleTimes_[endSampleTime_].
name()
457 <<
" with weight:" << s <<
endl;
460 fld = (1 -
s)*startSampledValues_ + s*endSampledValues_;
461 wantedAverage = (1 -
s)*startAverage_ + s*endAverage_;
468 Type averagePsi =
gSum(patch_.magSf()*
fld)/
gSum(patch_.magSf());
472 Pout<<
"updateCoeffs :" 473 <<
" actual average:" << averagePsi
474 <<
" wanted average:" << wantedAverage
478 if (
mag(averagePsi) < vSmall)
481 const Type offset = wantedAverage - averagePsi;
484 Pout<<
"updateCoeffs :" 485 <<
" offsetting with:" << offset <<
endl;
491 const scalar scale =
mag(wantedAverage)/
mag(averagePsi);
495 Pout<<
"updateCoeffs :" 496 <<
" scaling with:" << scale <<
endl;
505 const scalar t = time().timeOutputValue();
506 fld += offset_->value(t);
511 Pout<<
"updateCoeffs : set fixedValue to min:" <<
gMin(fld)
512 <<
" max:" <<
gMax(fld)
530 time().constant()/
"boundaryData"/patch_.
name(),
539 writeEntry(os,
"fieldTable", fieldTableName_);
545 word(
"planarInterpolation"),
virtual const fileName & name() const
Return the name of the stream.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
bool exists(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist (as directory or file) in the file system?
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.
errorManipArg< error, int > exit(error &err, const int errNo=1)
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.
Type gMin(const FieldField< Field, Type > &f)
tmp< Field< Type > > map()
Return the current mapped patch field.
T & ref() const
Return non-const reference or generate a fatal error.
void writeEntryIfDifferent(Ostream &os, const word &entryName, const EntryType &value1, const EntryType &value2)
Helper function to write the keyword and entry only if the.
Ostream & endl(Ostream &os)
Add newline and flush stream.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
void write(Ostream &) const
Write.
Patch field mapper which interpolates the values from a set of supplied points in space and time...
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){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
vectorField pointField
pointField is a vectorField.
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))
Type gSum(const FieldField< Field, Type > &f)
void rmap(const timeVaryingMappedFvPatchField< Type > &, const labelList &)
Reverse map the given timeVaryingMappedFvPatchField.
Pre-declare SubField and related Field type.
const dimensionedScalar & e
Elementary charge.
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
virtual const fileName & name() const
Return the name of the stream.
Foam::fvPatchFieldMapper.
word name() const
Return file name (part beyond last /)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Type gMax(const FieldField< Field, Type > &f)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
word name(const complex &)
Return a string representation of a complex.
timeVaryingMappedFvPatchField(const fvPatch &, const word &fieldName)
Construct from patch and field name.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Type gAverage(const FieldField< Field, Type > &f)
prefixOSstream Pout(cout, "Pout")
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual const word & name() const
Return name.
A class for managing temporary objects.