timeVaryingMappedFixedValueFvPatchField.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) 2011-2014 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 "AverageIOField.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
37 template<class Type>
40 (
41  const fvPatch& p,
43 )
44 :
46  fieldTableName_(iF.name()),
47  setAverage_(false),
48  perturb_(0),
49  mapperPtr_(NULL),
50  sampleTimes_(0),
51  startSampleTime_(-1),
52  startSampledValues_(0),
53  startAverage_(pTraits<Type>::zero),
54  endSampleTime_(-1),
55  endSampledValues_(0),
56  endAverage_(pTraits<Type>::zero),
57  offset_()
58 {}
59 
60 
61 template<class Type>
64 (
66  const fvPatch& p,
68  const fvPatchFieldMapper& mapper
69 )
70 :
71  fixedValueFvPatchField<Type>(ptf, p, iF, mapper),
72  fieldTableName_(ptf.fieldTableName_),
73  setAverage_(ptf.setAverage_),
74  perturb_(ptf.perturb_),
75  mapMethod_(ptf.mapMethod_),
76  mapperPtr_(NULL),
77  sampleTimes_(0),
78  startSampleTime_(-1),
79  startSampledValues_(0),
80  startAverage_(pTraits<Type>::zero),
81  endSampleTime_(-1),
82  endSampledValues_(0),
83  endAverage_(pTraits<Type>::zero),
84  offset_
85  (
86  ptf.offset_.valid()
87  ? ptf.offset_().clone().ptr()
88  : NULL
89  )
90 {}
91 
92 
93 template<class Type>
96 (
97  const fvPatch& p,
99  const dictionary& dict
100 )
101 :
103  fieldTableName_(iF.name()),
104  setAverage_(readBool(dict.lookup("setAverage"))),
105  perturb_(dict.lookupOrDefault("perturb", 1e-5)),
106  mapMethod_
107  (
108  dict.lookupOrDefault<word>
109  (
110  "mapMethod",
111  "planarInterpolation"
112  )
113  ),
114  mapperPtr_(NULL),
115  sampleTimes_(0),
116  startSampleTime_(-1),
117  startSampledValues_(0),
118  startAverage_(pTraits<Type>::zero),
119  endSampleTime_(-1),
120  endSampledValues_(0),
121  endAverage_(pTraits<Type>::zero),
122  offset_(DataEntry<Type>::New("offset", dict))
123 {
124  if
125  (
126  mapMethod_ != "planarInterpolation"
127  && mapMethod_ != "nearest"
128  )
129  {
131  (
132  "timeVaryingMappedFixedValueFvPatchField<Type>::\n"
133  "timeVaryingMappedFixedValueFvPatchField\n"
134  "(\n"
135  " const fvPatch&\n"
136  " const DimensionedField<Type, volMesh>&\n"
137  " const dictionary&\n"
138  ")\n",
139  dict
140  ) << "mapMethod should be one of 'planarInterpolation'"
141  << ", 'nearest'" << exit(FatalIOError);
142  }
143 
144 
145  dict.readIfPresent("fieldTableName", fieldTableName_);
146 
147  if (dict.found("value"))
148  {
149  fvPatchField<Type>::operator==(Field<Type>("value", dict, p.size()));
150  }
151  else
152  {
153  // Note: we use evaluate() here to trigger updateCoeffs followed
154  // by re-setting of fvatchfield::updated_ flag. This is
155  // so if first use is in the next time step it retriggers
156  // a new update.
157  this->evaluate(Pstream::blocking);
158  }
159 }
160 
161 
162 template<class Type>
165 (
167 )
168 :
170  fieldTableName_(ptf.fieldTableName_),
171  setAverage_(ptf.setAverage_),
172  perturb_(ptf.perturb_),
173  mapMethod_(ptf.mapMethod_),
174  mapperPtr_(NULL),
175  sampleTimes_(ptf.sampleTimes_),
176  startSampleTime_(ptf.startSampleTime_),
177  startSampledValues_(ptf.startSampledValues_),
178  startAverage_(ptf.startAverage_),
179  endSampleTime_(ptf.endSampleTime_),
180  endSampledValues_(ptf.endSampledValues_),
181  endAverage_(ptf.endAverage_),
182  offset_
183  (
184  ptf.offset_.valid()
185  ? ptf.offset_().clone().ptr()
186  : NULL
187  )
188 {}
189 
190 
191 template<class Type>
194 (
197 )
198 :
200  fieldTableName_(ptf.fieldTableName_),
201  setAverage_(ptf.setAverage_),
202  perturb_(ptf.perturb_),
203  mapMethod_(ptf.mapMethod_),
204  mapperPtr_(NULL),
205  sampleTimes_(ptf.sampleTimes_),
206  startSampleTime_(ptf.startSampleTime_),
207  startSampledValues_(ptf.startSampledValues_),
208  startAverage_(ptf.startAverage_),
209  endSampleTime_(ptf.endSampleTime_),
210  endSampledValues_(ptf.endSampledValues_),
211  endAverage_(ptf.endAverage_),
212  offset_
213  (
214  ptf.offset_.valid()
215  ? ptf.offset_().clone().ptr()
216  : NULL
217  )
218 {}
219 
220 
221 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
222 
223 template<class Type>
225 (
226  const fvPatchFieldMapper& m
227 )
228 {
230  if (startSampledValues_.size())
231  {
232  startSampledValues_.autoMap(m);
233  endSampledValues_.autoMap(m);
234  }
235  // Clear interpolator
236  mapperPtr_.clear();
237  startSampleTime_ = -1;
238  endSampleTime_ = -1;
239 }
240 
241 
242 template<class Type>
244 (
245  const fvPatchField<Type>& ptf,
246  const labelList& addr
247 )
248 {
250 
252  refCast<const timeVaryingMappedFixedValueFvPatchField<Type> >(ptf);
253 
254  startSampledValues_.rmap(tiptf.startSampledValues_, addr);
255  endSampledValues_.rmap(tiptf.endSampledValues_, addr);
256 
257  // Clear interpolator
258  mapperPtr_.clear();
259  startSampleTime_ = -1;
260  endSampleTime_ = -1;
261 }
262 
263 
264 template<class Type>
266 {
267  // Initialise
268  if (mapperPtr_.empty())
269  {
270  pointIOField samplePoints
271  (
272  IOobject
273  (
274  "points",
275  this->db().time().constant(),
276  "boundaryData"/this->patch().name(),
277  this->db(),
280  false
281  )
282  );
283 
284  const fileName samplePointsFile = samplePoints.filePath();
285 
286  if (debug)
287  {
288  Info<< "timeVaryingMappedFixedValueFvPatchField :"
289  << " Read " << samplePoints.size() << " sample points from "
290  << samplePointsFile << endl;
291  }
292 
293 
294  // tbd: run-time selection
295  bool nearestOnly =
296  (
297  !mapMethod_.empty()
298  && mapMethod_ != "planarInterpolation"
299  );
300 
301  // Allocate the interpolator
302  mapperPtr_.reset
303  (
305  (
306  samplePoints,
307  this->patch().patch().faceCentres(),
308  perturb_,
309  nearestOnly
310  )
311  );
312 
313  // Read the times for which data is available
314  const fileName samplePointsDir = samplePointsFile.path();
315  sampleTimes_ = Time::findTimes(samplePointsDir);
316 
317  if (debug)
318  {
319  Info<< "timeVaryingMappedFixedValueFvPatchField : In directory "
320  << samplePointsDir << " found times "
322  << endl;
323  }
324  }
325 
326 
327  // Find current time in sampleTimes
328  label lo = -1;
329  label hi = -1;
330 
331  bool foundTime = mapperPtr_().findTime
332  (
333  sampleTimes_,
334  startSampleTime_,
335  this->db().time().value(),
336  lo,
337  hi
338  );
339 
340  if (!foundTime)
341  {
343  (
344  "timeVaryingMappedFixedValueFvPatchField<Type>::checkTable()"
345  ) << "Cannot find starting sampling values for current time "
346  << this->db().time().value() << nl
347  << "Have sampling values for times "
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 
389  // Reread values and interpolate
391  (
392  IOobject
393  (
394  fieldTableName_,
395  this->db().time().constant(),
396  "boundaryData"
397  /this->patch().name()
398  /sampleTimes_[startSampleTime_].name(),
399  this->db(),
402  false
403  )
404  );
405 
406  if (vals.size() != mapperPtr_().sourceSize())
407  {
409  (
410  "timeVaryingMappedFixedValueFvPatchField<Type>::"
411  "checkTable()"
412  ) << "Number of values (" << vals.size()
413  << ") differs from the number of points ("
414  << mapperPtr_().sourceSize()
415  << ") in file " << vals.objectPath() << exit(FatalError);
416  }
417 
418  startAverage_ = vals.average();
419  startSampledValues_ = mapperPtr_().interpolate(vals);
420  }
421  }
422 
423  if (hi != endSampleTime_)
424  {
425  endSampleTime_ = hi;
426 
427  if (endSampleTime_ == -1)
428  {
429  // endTime no longer valid. Might as well clear endValues.
430  if (debug)
431  {
432  Pout<< "checkTable : Clearing endValues" << endl;
433  }
434  endSampledValues_.clear();
435  }
436  else
437  {
438  if (debug)
439  {
440  Pout<< "checkTable : Reading endValues from "
441  << "boundaryData"
442  /this->patch().name()
443  /sampleTimes_[endSampleTime_].name()
444  << endl;
445  }
446 
447  // Reread values and interpolate
449  (
450  IOobject
451  (
452  fieldTableName_,
453  this->db().time().constant(),
454  "boundaryData"
455  /this->patch().name()
456  /sampleTimes_[endSampleTime_].name(),
457  this->db(),
460  false
461  )
462  );
463 
464  if (vals.size() != mapperPtr_().sourceSize())
465  {
467  (
468  "timeVaryingMappedFixedValueFvPatchField<Type>::"
469  "checkTable()"
470  ) << "Number of values (" << vals.size()
471  << ") differs from the number of points ("
472  << mapperPtr_().sourceSize()
473  << ") in file " << vals.objectPath() << exit(FatalError);
474  }
475 
476  endAverage_ = vals.average();
477  endSampledValues_ = mapperPtr_().interpolate(vals);
478  }
479  }
480 }
481 
482 
483 template<class Type>
485 {
486  if (this->updated())
487  {
488  return;
489  }
490 
491 
492  checkTable();
493 
494  // Interpolate between the sampled data
495 
496  Type wantedAverage;
497 
498  if (endSampleTime_ == -1)
499  {
500  // only start value
501  if (debug)
502  {
503  Pout<< "updateCoeffs : Sampled, non-interpolated values"
504  << " from start time:"
505  << sampleTimes_[startSampleTime_].name() << nl;
506  }
507 
508  this->operator==(startSampledValues_);
509  wantedAverage = startAverage_;
510  }
511  else
512  {
513  scalar start = sampleTimes_[startSampleTime_].value();
514  scalar end = sampleTimes_[endSampleTime_].value();
515 
516  scalar s = (this->db().time().value() - start)/(end - start);
517 
518  if (debug)
519  {
520  Pout<< "updateCoeffs : Sampled, interpolated values"
521  << " between start time:"
522  << sampleTimes_[startSampleTime_].name()
523  << " and end time:" << sampleTimes_[endSampleTime_].name()
524  << " with weight:" << s << endl;
525  }
526 
527  this->operator==((1 - s)*startSampledValues_ + s*endSampledValues_);
528  wantedAverage = (1 - s)*startAverage_ + s*endAverage_;
529  }
530 
531  // Enforce average. Either by scaling (if scaling factor > 0.5) or by
532  // offsetting.
533  if (setAverage_)
534  {
535  const Field<Type>& fld = *this;
536 
537  Type averagePsi =
538  gSum(this->patch().magSf()*fld)
539  /gSum(this->patch().magSf());
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  const scalar t = this->db().time().timeOutputValue();
575  this->operator==(*this + offset_->value(t));
576 
577  if (debug)
578  {
579  Pout<< "updateCoeffs : set fixedValue to min:" << gMin(*this)
580  << " max:" << gMax(*this)
581  << " avg:" << gAverage(*this) << endl;
582  }
583 
585 }
586 
587 
588 template<class Type>
590 {
592  os.writeKeyword("setAverage") << setAverage_ << token::END_STATEMENT << nl;
593  if (perturb_ != 1e-5)
594  {
595  os.writeKeyword("perturb") << perturb_ << token::END_STATEMENT << nl;
596  }
597 
598  if (fieldTableName_ != this->dimensionedInternalField().name())
599  {
600  os.writeKeyword("fieldTableName") << fieldTableName_
601  << token::END_STATEMENT << nl;
602  }
603 
604  if
605  (
606  (
607  !mapMethod_.empty()
608  && mapMethod_ != "planarInterpolation"
609  )
610  )
611  {
612  os.writeKeyword("mapMethod") << mapMethod_
613  << token::END_STATEMENT << nl;
614  }
615 
616  offset_->writeData(os);
617 
618  this->writeEntry("value", os);
619 }
620 
621 
622 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
623 
624 } // End namespace Foam
625 
626 // ************************************************************************* //
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
prefixOSstream Pout(cout,"Pout")
Definition: IOstreams.H:53
word name() const
Return file name (part beyond last /)
Definition: fileName.C:206
virtual void operator==(const fvPatchField< Type > &)
Definition: fvPatchField.C:559
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
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))
A class for handling file names.
Definition: fileName.H:69
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
#define FatalIOErrorIn(fn, ios)
Definition: error.H:326
Pre-declare SubField and related Field type.
Definition: Field.H:57
fileName filePath() const
Return complete path + object name if the file exists.
Definition: IOobject.C:320
This boundary conditions interpolates the values from a set of supplied points in space and time...
rDeltaT dimensionedInternalField()
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:248
bool readBool(Istream &)
Definition: boolIO.C:60
error FatalError
A primitive field + average with IO.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: DataEntry.H:52
dimensioned< scalar > mag(const dimensioned< Type > &)
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))
virtual const fileName & name() const
Return the name of the stream.
Definition: OSstream.H:88
void checkTable()
Find boundary data inbetween current time and interpolate.
Type gSum(const FieldField< Field, Type > &f)
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
static wordList timeNames(const instantList &)
Helper: extract words of times.
Type gAverage(const FieldField< Field, Type > &f)
timeVaryingMappedFixedValueFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
virtual label size() const
Return size.
Definition: fvPatch.H:161
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:301
fileName objectPath() const
Return complete path + object name.
Definition: IOobject.H:352
Type gMax(const FieldField< Field, Type > &f)
Interpolates between two sets of unstructured points using 2D Delaunay triangulation. Used in e.g. timeVaryingMapped bcs.
A class for handling words, derived from string.
Definition: word.H:59
fileName path() const
Return directory path name (part before last /)
Definition: fileName.C:287
const Type & average() const
static instantList findTimes(const fileName &, const word &constantName="constant")
Search a given directory for valid time directories.
Definition: findTimes.C:38
messageStream Info
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:354
Type gMin(const FieldField< Field, Type > &f)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
IOerror FatalIOError
Traits class for primitives.
Definition: pTraits.H:50
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define FatalErrorIn(fn)
Definition: error.H:317
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
static const char nl
Definition: Ostream.H:257
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: fvPatchField.C:291
Constant dispersed-phase particle diameter model.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void rmap(const fvPatchField< Type > &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
const double e
Elementary charge.
Definition: doubleFloat.H:78
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
const word & name() const
Return name.
Definition: IOobject.H:252
volScalarField & p
Definition: createFields.H:51
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: fvPatchField.C:230
Foam::fvPatchFieldMapper.