timeVaryingMappedFixedValuePointPatchField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2012-2017 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
27 #include "Time.H"
28 #include "AverageField.H"
29 #include "IFstream.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class Type>
36 (
37  const pointPatch& p,
39 )
40 :
42  fieldTableName_(iF.name()),
43  setAverage_(false),
44  perturb_(0),
45  mapperPtr_(nullptr),
46  sampleTimes_(0),
47  startSampleTime_(-1),
48  startSampledValues_(0),
49  startAverage_(Zero),
50  endSampleTime_(-1),
51  endSampledValues_(0),
52  endAverage_(Zero),
53  offset_()
54 {}
55 
56 
57 template<class Type>
60 (
61  const pointPatch& p,
63  const dictionary& dict
64 )
65 :
67  fieldTableName_(iF.name()),
68  setAverage_(dict.lookupOrDefault("setAverage", false)),
69  perturb_(dict.lookupOrDefault("perturb", 1e-5)),
70  mapMethod_
71  (
72  dict.lookupOrDefault<word>
73  (
74  "mapMethod",
75  "planarInterpolation"
76  )
77  ),
78  mapperPtr_(nullptr),
79  sampleTimes_(0),
80  startSampleTime_(-1),
81  startSampledValues_(0),
82  startAverage_(Zero),
83  endSampleTime_(-1),
84  endSampledValues_(0),
85  endAverage_(Zero),
86  offset_()
87 {
88  if (dict.found("offset"))
89  {
90  offset_ = Function1<Type>::New("offset", dict);
91  }
92 
93  if
94  (
95  mapMethod_ != "planarInterpolation"
96  && mapMethod_ != "nearest"
97  )
98  {
100  (
101  dict
102  ) << "mapMethod should be one of 'planarInterpolation'"
103  << ", 'nearest'" << exit(FatalIOError);
104  }
105 
106  dict.readIfPresent("fieldTableName", fieldTableName_);
107 
108  if (dict.found("value"))
109  {
111  (
112  Field<Type>("value", dict, p.size())
113  );
114  }
115  else
116  {
117  // Note: use evaluate to do updateCoeffs followed by a reset
118  // of the pointPatchField::updated_ flag. This is
119  // so if first use is in the next time step it retriggers
120  // a new update.
121  pointPatchField<Type>::evaluate(Pstream::commsTypes::blocking);
122  }
123 }
124 
125 
126 template<class Type>
129 (
131  const pointPatch& p,
133  const pointPatchFieldMapper& mapper
134 )
135 :
136  fixedValuePointPatchField<Type>(ptf, p, iF, mapper),
137  fieldTableName_(ptf.fieldTableName_),
138  setAverage_(ptf.setAverage_),
139  perturb_(ptf.perturb_),
140  mapMethod_(ptf.mapMethod_),
141  mapperPtr_(nullptr),
142  sampleTimes_(0),
143  startSampleTime_(-1),
144  startSampledValues_(0),
145  startAverage_(Zero),
146  endSampleTime_(-1),
147  endSampledValues_(0),
148  endAverage_(Zero),
149  offset_(ptf.offset_, false)
150 {}
151 
152 
153 template<class Type>
156 (
158 )
159 :
161  fieldTableName_(ptf.fieldTableName_),
162  setAverage_(ptf.setAverage_),
163  perturb_(ptf.perturb_),
164  mapMethod_(ptf.mapMethod_),
165  mapperPtr_(ptf.mapperPtr_),
166  sampleTimes_(ptf.sampleTimes_),
167  startSampleTime_(ptf.startSampleTime_),
168  startSampledValues_(ptf.startSampledValues_),
169  startAverage_(ptf.startAverage_),
170  endSampleTime_(ptf.endSampleTime_),
171  endSampledValues_(ptf.endSampledValues_),
172  endAverage_(ptf.endAverage_),
173  offset_(ptf.offset_, false)
174 {}
175 
176 
177 template<class Type>
180 (
183 )
184 :
186  fieldTableName_(ptf.fieldTableName_),
187  setAverage_(ptf.setAverage_),
188  perturb_(ptf.perturb_),
189  mapMethod_(ptf.mapMethod_),
190  mapperPtr_(ptf.mapperPtr_),
191  sampleTimes_(ptf.sampleTimes_),
192  startSampleTime_(ptf.startSampleTime_),
193  startSampledValues_(ptf.startSampledValues_),
194  startAverage_(ptf.startAverage_),
195  endSampleTime_(ptf.endSampleTime_),
196  endSampledValues_(ptf.endSampledValues_),
197  endAverage_(ptf.endAverage_),
198  offset_(ptf.offset_, false)
199 {}
200 
201 
202 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
203 
204 template<class Type>
206 (
207  const pointPatchFieldMapper& m
208 )
209 {
211  if (startSampledValues_.size())
212  {
213  startSampledValues_.autoMap(m);
214  endSampledValues_.autoMap(m);
215  }
216  // Clear interpolator
217  mapperPtr_.clear();
218  startSampleTime_ = -1;
219  endSampleTime_ = -1;
220 }
221 
222 
223 template<class Type>
225 (
226  const pointPatchField<Type>& ptf,
227  const labelList& addr
228 )
229 {
231 
233  refCast<const timeVaryingMappedFixedValuePointPatchField<Type>>(ptf);
234 
235  startSampledValues_.rmap(tiptf.startSampledValues_, addr);
236  endSampledValues_.rmap(tiptf.endSampledValues_, addr);
237 
238  // Clear interpolator
239  mapperPtr_.clear();
240  startSampleTime_ = -1;
241  endSampleTime_ = -1;
242 }
243 
244 
245 template<class Type>
247 {
248  // Initialise
249  if (startSampleTime_ == -1 && endSampleTime_ == -1)
250  {
251  const polyMesh& pMesh = this->patch().boundaryMesh().mesh()();
252 
253  // Read the initial point position
254  pointField meshPts;
255 
256  if (pMesh.pointsInstance() == pMesh.facesInstance())
257  {
258  meshPts = pointField(pMesh.points(), this->patch().meshPoints());
259  }
260  else
261  {
262  // Load points from facesInstance
263  if (debug)
264  {
265  Info<< "Reloading points0 from " << pMesh.facesInstance()
266  << endl;
267  }
268 
269  pointIOField points0
270  (
271  IOobject
272  (
273  "points",
274  pMesh.facesInstance(),
275  polyMesh::meshSubDir,
276  pMesh,
277  IOobject::MUST_READ,
278  IOobject::NO_WRITE,
279  false
280  )
281  );
282  meshPts = pointField(points0, this->patch().meshPoints());
283  }
284 
285  // Reread values and interpolate
286  fileName samplePointsFile
287  (
288  this->db().time().constant()
289  /"boundaryData"
290  /this->patch().name()
291  /"points"
292  );
293 
294  pointField samplePoints((IFstream(samplePointsFile)()));
295 
296  // tbd: run-time selection
297  bool nearestOnly =
298  (
299  !mapMethod_.empty()
300  && mapMethod_ != "planarInterpolation"
301  );
302 
303  // Allocate the interpolator
304  mapperPtr_.reset
305  (
307  (
308  samplePoints,
309  meshPts,
310  perturb_,
311  nearestOnly
312  )
313  );
314 
315  // Read the times for which data is available
316 
317  const fileName samplePointsDir = samplePointsFile.path();
318  sampleTimes_ = Time::findTimes(samplePointsDir);
319 
320  if (debug)
321  {
322  Info<< "timeVaryingMappedFixedValuePointPatchField : In directory "
323  << samplePointsDir << " found times "
324  << pointToPointPlanarInterpolation::timeNames(sampleTimes_)
325  << endl;
326  }
327  }
328 
329  // Find current time in sampleTimes
330  label lo = -1;
331  label hi = -1;
332 
333  bool foundTime = mapperPtr_().findTime
334  (
335  sampleTimes_,
336  startSampleTime_,
337  this->db().time().value(),
338  lo,
339  hi
340  );
341 
342  if (!foundTime)
343  {
345  << "Cannot find starting sampling values for current time "
346  << this->db().time().value() << nl
347  << "Have sampling values for times "
348  << pointToPointPlanarInterpolation::timeNames(sampleTimes_) << nl
349  << "In directory "
350  << this->db().time().constant()/"boundaryData"/this->patch().name()
351  << "\n on patch " << this->patch().name()
352  << " of field " << fieldTableName_
353  << exit(FatalError);
354  }
355 
356 
357  // Update sampled data fields.
358 
359  if (lo != startSampleTime_)
360  {
361  startSampleTime_ = lo;
362 
363  if (startSampleTime_ == endSampleTime_)
364  {
365  // No need to reread since are end values
366  if (debug)
367  {
368  Pout<< "checkTable : Setting startValues to (already read) "
369  << "boundaryData"
370  /this->patch().name()
371  /sampleTimes_[startSampleTime_].name()
372  << endl;
373  }
374  startSampledValues_ = endSampledValues_;
375  startAverage_ = endAverage_;
376  }
377  else
378  {
379  if (debug)
380  {
381  Pout<< "checkTable : Reading startValues from "
382  << "boundaryData"
383  /this->patch().name()
384  /sampleTimes_[lo].name()
385  << endl;
386  }
387 
388  // Reread values and interpolate
389  fileName valsFile
390  (
391  this->db().time().constant()
392  /"boundaryData"
393  /this->patch().name()
394  /sampleTimes_[startSampleTime_].name()
395  /fieldTableName_
396  );
397 
398  Field<Type> vals;
399 
400  if (setAverage_)
401  {
402  AverageField<Type> avals((IFstream(valsFile)()));
403  vals = avals;
404  startAverage_ = avals.average();
405  }
406  else
407  {
408  IFstream(valsFile)() >> vals;
409  }
410 
411  if (vals.size() != mapperPtr_().sourceSize())
412  {
414  << "Number of values (" << vals.size()
415  << ") differs from the number of points ("
416  << mapperPtr_().sourceSize()
417  << ") in file " << valsFile << exit(FatalError);
418  }
419 
420  startSampledValues_ = mapperPtr_().interpolate(vals);
421  }
422  }
423 
424  if (hi != endSampleTime_)
425  {
426  endSampleTime_ = hi;
427 
428  if (endSampleTime_ == -1)
429  {
430  // endTime no longer valid. Might as well clear endValues.
431  if (debug)
432  {
433  Pout<< "checkTable : Clearing endValues" << endl;
434  }
435  endSampledValues_.clear();
436  }
437  else
438  {
439  if (debug)
440  {
441  Pout<< "checkTable : Reading endValues from "
442  << "boundaryData"
443  /this->patch().name()
444  /sampleTimes_[endSampleTime_].name()
445  << endl;
446  }
447 
448  // Reread values and interpolate
449  fileName valsFile
450  (
451  this->db().time().constant()
452  /"boundaryData"
453  /this->patch().name()
454  /sampleTimes_[endSampleTime_].name()
455  /fieldTableName_
456  );
457 
458  Field<Type> vals;
459 
460  if (setAverage_)
461  {
462  AverageField<Type> avals((IFstream(valsFile)()));
463  vals = avals;
464  endAverage_ = avals.average();
465  }
466  else
467  {
468  IFstream(valsFile)() >> vals;
469  }
470 
471  if (vals.size() != mapperPtr_().sourceSize())
472  {
474  << "Number of values (" << vals.size()
475  << ") differs from the number of points ("
476  << mapperPtr_().sourceSize()
477  << ") in file " << valsFile << exit(FatalError);
478  }
479 
480  endSampledValues_ = mapperPtr_().interpolate(vals);
481  }
482  }
483 }
484 
485 
486 template<class Type>
488 {
489  if (this->updated())
490  {
491  return;
492  }
493 
494  checkTable();
495 
496  // Interpolate between the sampled data
497 
498  Type wantedAverage;
499 
500  if (endSampleTime_ == -1)
501  {
502  // only start value
503  if (debug)
504  {
505  Pout<< "updateCoeffs : Sampled, non-interpolated values"
506  << " from start time:"
507  << sampleTimes_[startSampleTime_].name() << nl;
508  }
509 
510  this->operator==(startSampledValues_);
511  wantedAverage = startAverage_;
512  }
513  else
514  {
515  scalar start = sampleTimes_[startSampleTime_].value();
516  scalar end = sampleTimes_[endSampleTime_].value();
517 
518  scalar s = (this->db().time().value()-start)/(end-start);
519 
520  if (debug)
521  {
522  Pout<< "updateCoeffs : Sampled, interpolated values"
523  << " between start time:"
524  << sampleTimes_[startSampleTime_].name()
525  << " and end time:" << sampleTimes_[endSampleTime_].name()
526  << " with weight:" << s << endl;
527  }
528 
529  this->operator==((1-s)*startSampledValues_ + s*endSampledValues_);
530  wantedAverage = (1-s)*startAverage_ + s*endAverage_;
531  }
532 
533  // Enforce average. Either by scaling (if scaling factor > 0.5) or by
534  // offsetting.
535  if (setAverage_)
536  {
537  const Field<Type>& fld = *this;
538 
539  Type averagePsi = gAverage(fld);
540 
541  if (debug)
542  {
543  Pout<< "updateCoeffs :"
544  << " actual average:" << averagePsi
545  << " wanted average:" << wantedAverage
546  << endl;
547  }
548 
549  if (mag(averagePsi) < VSMALL)
550  {
551  // Field too small to scale. Offset instead.
552  const Type offset = wantedAverage - averagePsi;
553  if (debug)
554  {
555  Pout<< "updateCoeffs :"
556  << " offsetting with:" << offset << endl;
557  }
558  this->operator==(fld+offset);
559  }
560  else
561  {
562  const scalar scale = mag(wantedAverage)/mag(averagePsi);
563 
564  if (debug)
565  {
566  Pout<< "updateCoeffs :"
567  << " scaling with:" << scale << endl;
568  }
569  this->operator==(scale*fld);
570  }
571  }
572 
573  // Apply offset to mapped values
574  if (offset_.valid())
575  {
576  const scalar t = this->db().time().timeOutputValue();
577  this->operator==(*this + offset_->value(t));
578  }
579 
580  if (debug)
581  {
582  Pout<< "updateCoeffs : set fixedValue to min:" << gMin(*this)
583  << " max:" << gMax(*this)
584  << " avg:" << gAverage(*this) << endl;
585  }
586 
588 }
589 
590 
591 template<class Type>
593 (
594  Ostream& os
595 ) const
596 {
598 
599  this->writeEntryIfDifferent(os, "setAverage", Switch(false), setAverage_);
600 
601  this->writeEntryIfDifferent(os, "perturb", scalar(1e-5), perturb_);
602 
603  this->writeEntryIfDifferent
604  (
605  os,
606  "fieldTable",
607  this->internalField().name(),
608  fieldTableName_
609  );
610 
611  this->writeEntryIfDifferent
612  (
613  os,
614  "mapMethod",
615  word("planarInterpolation"),
616  mapMethod_
617  );
618 
619  if (offset_.valid())
620  {
621  offset_->writeData(os);
622  }
623 }
624 
625 
626 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:88
dictionary dict
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const word & name() const
Return name.
Definition: IOobject.H:291
virtual void rmap(const pointPatchField< Type > &, const labelList &)
Reverse map the given PointPatchField onto.
A class for handling file names.
Definition: fileName.H:69
void checkTable()
Find boundary data inbetween current time and interpolate.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:801
A time-varying form of a mapped fixed value boundary condition.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
A primitive field with a separate average value.
Definition: AverageField.H:49
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
A FixedValue boundary condition for pointField.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Foam::pointPatchFieldMapper.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
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.
Definition: Switch.H:60
label size() const
Return size.
timeVaryingMappedFixedValuePointPatchField(const pointPatch &, const DimensionedField< Type, pointMesh > &)
Construct from patch and internal field.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1011
Abstract base class for point-mesh patch fields.
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
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.
Definition: pointFieldFwd.H:42
const dimensionedScalar e
Elementary charge.
Definition: doubleFloat.H:78
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:795
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))
Pre-declare SubField and related Field type.
Definition: Field.H:57
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
word name() const
Return file name (part beyond last /)
Definition: fileName.C:180
const polyMesh & mesh() const
Return the mesh reference.
static const zero Zero
Definition: zero.H:91
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
Type gMax(const FieldField< Field, Type > &f)
Input from file stream.
Definition: IFstream.H:81
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Constant dispersed-phase particle diameter model.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:56
Type gAverage(const FieldField< Field, Type > &f)
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
fileName path() const
Return directory path name (part before last /)
Definition: fileName.C:250
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
volScalarField & p
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const Type & average() const
Definition: AverageField.C:61
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Interpolates between two sets of unstructured points using 2D Delaunay triangulation. Used in e.g. timeVaryingMapped bcs.
IOerror FatalIOError