externalCoupledMixedFvPatchField.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) 2013 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 "fvPatchFieldMapper.H"
28 #include "volFields.H"
29 #include "IFstream.H"
30 #include "globalIndex.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 template<class Type>
35 Foam::word Foam::externalCoupledMixedFvPatchField<Type>::lockName = "OpenFOAM";
36 
37 template<class Type>
39 Foam::externalCoupledMixedFvPatchField<Type>::patchKey = "# Patch: ";
40 
41 
42 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
43 
44 template<class Type>
46 {
47  word regionName(this->dimensionedInternalField().mesh().name());
48  if (regionName == polyMesh::defaultRegion)
49  {
50  regionName = ".";
51  }
52 
53  fileName result(commsDir_/regionName);
54  result.clean();
55 
56  return result;
57 }
58 
59 
60 template<class Type>
62 (
63  const labelList& patchIDs
64 )
65 {
66  const volFieldType& cvf =
67  static_cast<const volFieldType&>(this->dimensionedInternalField());
68 
69  volFieldType& vf = const_cast<volFieldType&>(cvf);
70 
71  typename volFieldType::GeometricBoundaryField& bf = vf.boundaryField();
72 
73  // number of patches can be different in parallel...
74  label nPatch = bf.size();
75  reduce(nPatch, maxOp<label>());
76 
77  offsets_.setSize(nPatch);
78  forAll(offsets_, i)
79  {
80  offsets_[i].setSize(Pstream::nProcs());
81  offsets_[i] = 0;
82  }
83 
84  // set the master patch
85  forAll(patchIDs, i)
86  {
87  label patchI = patchIDs[i];
88 
89  patchType& pf = refCast<patchType>(bf[patchI]);
90 
91  offsets_[patchI][Pstream::myProcNo()] = pf.size();
92 
93  if (i == 0)
94  {
95  pf.master() = true;
96  }
97  else
98  {
99  pf.master() = false;
100  }
101  }
102 
103  // set the patch offsets
104  int tag = Pstream::msgType() + 1;
105  forAll(offsets_, i)
106  {
107  Pstream::gatherList(offsets_[i], tag);
108  Pstream::scatterList(offsets_[i], tag);
109  }
110 
111  label patchOffset = 0;
112  forAll(offsets_, patchI)
113  {
114  label sumOffset = 0;
115  List<label>& procOffsets = offsets_[patchI];
116 
117  forAll(procOffsets, procI)
118  {
119  label o = procOffsets[procI];
120  if (o > 0)
121  {
122  procOffsets[procI] = patchOffset + sumOffset;
123  sumOffset += o;
124  }
125  }
126  patchOffset += sumOffset;
127  }
128 }
129 
130 
131 template<class Type>
133 (
134  OFstream& osPoints,
135  OFstream& osFaces
136 ) const
137 {
138  int tag = Pstream::msgType() + 1;
139 
140  const label procI = Pstream::myProcNo();
141  const polyPatch& p = this->patch().patch();
142  const polyMesh& mesh = p.boundaryMesh().mesh();
143 
144  labelList pointToGlobal;
145  labelList uniquePointIDs;
146  (void)mesh.globalData().mergePoints
147  (
148  p.meshPoints(),
149  p.meshPointMap(),
150  pointToGlobal,
151  uniquePointIDs
152  );
153 
154  List<pointField> allPoints(Pstream::nProcs());
155  allPoints[procI] = pointField(mesh.points(), uniquePointIDs);
156  Pstream::gatherList(allPoints, tag);
157 
158  List<faceList> allFaces(Pstream::nProcs());
159  faceList& patchFaces = allFaces[procI];
160  patchFaces = p.localFaces();
161  forAll(patchFaces, faceI)
162  {
163  inplaceRenumber(pointToGlobal, patchFaces[faceI]);
164  }
165 
166  Pstream::gatherList(allFaces, tag);
167 
168  if (Pstream::master())
169  {
170  pointField pts
171  (
172  ListListOps::combine<pointField>(allPoints, accessOp<pointField>())
173  );
174 
175  // write points
176  osPoints << patchKey.c_str() << this->patch().name() << pts << endl;
177 
178  faceList fcs
179  (
180  ListListOps::combine<faceList>(allFaces, accessOp<faceList>())
181  );
182 
183  // write faces
184  osFaces<< patchKey.c_str() << this->patch().name() << fcs << endl;
185  }
186 }
187 
188 
189 template<class Type>
191 {
192  return fileName(baseDir()/(lockName + ".lock"));
193 }
194 
195 
196 template<class Type>
198 {
199  if (!master_ || !Pstream::master())
200  {
201  return;
202  }
203 
204  const fileName fName(lockFile());
205  IFstream is(fName);
206 
207  // only create lock file if it doesn't already exist
208  if (!is.good())
209  {
210  if (log_)
211  {
212  Info<< type() << ": creating lock file" << endl;
213  }
214 
215  OFstream os(fName);
216  os << "lock file";
217  os.flush();
218  }
219 }
220 
221 
222 template<class Type>
224 {
225  if (!master_ || !Pstream::master())
226  {
227  return;
228  }
229 
230  if (log_)
231  {
232  Info<< type() << ": removing lock file" << endl;
233  }
234 
235  rm(lockFile());
236 }
237 
238 
239 template<class Type>
241 {
242  // only wait on master patch
243 
244  const volFieldType& cvf =
245  static_cast<const volFieldType&>(this->dimensionedInternalField());
246 
247  const typename volFieldType::GeometricBoundaryField& bf =
248  cvf.boundaryField();
249 
250  forAll(coupledPatchIDs_, i)
251  {
252  label patchI = coupledPatchIDs_[i];
253 
254  const patchType& pf = refCast<const patchType>(bf[patchI]);
255 
256  if (pf.master())
257  {
258  pf.wait();
259  break;
260  }
261  }
262 }
263 
264 
265 template<class Type>
267 {
268  const fileName fName(lockFile());
269  label found = 0;
270  label totalTime = 0;
271 
272  if (log_)
273  {
274  Info<< type() << ": beginning wait for lock file " << fName << endl;
275  }
276 
277  while (found == 0)
278  {
279  if (Pstream::master())
280  {
281  if (totalTime > timeOut_)
282  {
284  (
285  "void "
286  "Foam::externalCoupledMixedFvPatchField<Type>::wait() "
287  "const"
288  )
289  << "Wait time exceeded time out time of " << timeOut_
290  << " s" << abort(FatalError);
291  }
292 
293  IFstream is(fName);
294 
295  if (is.good())
296  {
297  found++;
298 
299  if (log_)
300  {
301  Info<< type() << ": found lock file " << fName << endl;
302  }
303  }
304  else
305  {
306  sleep(waitInterval_);
307  totalTime += waitInterval_;
308 
309  if (log_)
310  {
311  Info<< type() << ": wait time = " << totalTime << endl;
312  }
313  }
314  }
315 
316  // prevent other procs from racing ahead
317  reduce(found, sumOp<label>());
318  }
319 }
320 
321 
322 template<class Type>
324 (
325  IFstream& is
326 ) const
327 {
328  if (!is.good())
329  {
331  (
332  "void Foam::externalCoupledMixedFvPatchField<Type>::"
333  "initialiseRead"
334  "("
335  "IFstream&"
336  ") const"
337  )
338  << "Unable to open data transfer file " << is.name()
339  << " for patch " << this->patch().name()
340  << exit(FatalError);
341  }
342 
343  label offset = offsets_[this->patch().index()][Pstream::myProcNo()];
344 
345  // scan forward to start reading at relevant line/position
346  string line;
347  for (label i = 0; i < offset; i++)
348  {
349  if (is.good())
350  {
351  is.getLine(line);
352  }
353  else
354  {
356  (
357  "void Foam::externalCoupledMixedFvPatchField<Type>::"
358  "initialiseRead"
359  "("
360  "IFstream&"
361  ") const"
362  )
363  << "Unable to scan forward to appropriate read position for "
364  << "data transfer file " << is.name()
365  << " for patch " << this->patch().name()
366  << exit(FatalError);
367  }
368  }
369 }
370 
371 
372 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
373 
374 template<class Type>
376 (
377  const fileName& transferFile
378 )
379 {
380  // read data passed back from external source
381  IFstream is(transferFile + ".in");
382 
383  // pre-process the input transfer file
384  initialiseRead(is);
385 
386  // read data from file
387  forAll(this->patch(), faceI)
388  {
389  if (is.good())
390  {
391  is >> this->refValue()[faceI]
392  >> this->refGrad()[faceI]
393  >> this->valueFraction()[faceI];
394  }
395  else
396  {
398  (
399  "void Foam::externalCoupledMixedFvPatchField<Type>::readData"
400  "("
401  "const fileName&"
402  ")"
403  )
404  << "Insufficient data for patch "
405  << this->patch().name()
406  << " in file " << is.name()
407  << exit(FatalError);
408  }
409  }
410 
411  initialised_ = true;
412 
413  // update the value from the mixed condition
415 }
416 
417 
418 template<class Type>
420 (
421  const fileName& transferFile
422 ) const
423 {
424  if (!master_)
425  {
426  return;
427  }
428 
429  OFstream os(transferFile);
430 
431  writeHeader(os);
432 
433  const volFieldType& cvf =
434  static_cast<const volFieldType&>(this->dimensionedInternalField());
435 
436  const typename volFieldType::GeometricBoundaryField& bf =
437  cvf.boundaryField();
438 
439  forAll(coupledPatchIDs_, i)
440  {
441  label patchI = coupledPatchIDs_[i];
442 
443  const patchType& pf = refCast<const patchType>(bf[patchI]);
444 
445  pf.transferData(os);
446  }
447 }
448 
449 
450 template<class Type>
452 (
453  OFstream& os
454 ) const
455 {
456  os << "# Values: magSf value snGrad" << endl;
457 }
458 
459 
460 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
461 
462 template<class Type>
465 (
466  const fvPatch& p,
468 )
469 :
471  commsDir_("unknown-commsDir"),
472  fName_("unknown-fName"),
473  waitInterval_(0),
474  timeOut_(0),
475  calcFrequency_(0),
476  initByExternal_(false),
477  log_(false),
478  master_(false),
479  offsets_(),
480  initialised_(false),
481  coupledPatchIDs_()
482 {
483  this->refValue() = pTraits<Type>::zero;
484  this->refGrad() = pTraits<Type>::zero;
485  this->valueFraction() = 0.0;
486 }
487 
488 
489 template<class Type>
492 (
494  const fvPatch& p,
496  const fvPatchFieldMapper& mapper
497 )
498 :
499  mixedFvPatchField<Type>(ptf, p, iF, mapper),
500  commsDir_(ptf.commsDir_),
501  fName_(ptf.fName_),
502  waitInterval_(ptf.waitInterval_),
503  timeOut_(ptf.timeOut_),
504  calcFrequency_(ptf.calcFrequency_),
505  initByExternal_(ptf.initByExternal_),
506  log_(ptf.log_),
507  master_(ptf.master_),
508  offsets_(ptf.offsets_),
509  initialised_(ptf.initialised_),
510  coupledPatchIDs_(ptf.coupledPatchIDs_)
511 {}
512 
513 
514 template<class Type>
517 (
518  const fvPatch& p,
520  const dictionary& dict
521 )
522 :
524  commsDir_(dict.lookup("commsDir")),
525  fName_(dict.lookup("fileName")),
526  waitInterval_(dict.lookupOrDefault("waitInterval", 1)),
527  timeOut_(dict.lookupOrDefault("timeOut", 100*waitInterval_)),
528  calcFrequency_(dict.lookupOrDefault("calcFrequency", 1)),
529  initByExternal_(readBool(dict.lookup("initByExternal"))),
530  log_(dict.lookupOrDefault("log", false)),
531  master_(true),
532  offsets_(),
533  initialised_(false),
534  coupledPatchIDs_()
535 {
536  if (dict.found("value"))
537  {
539  (
540  Field<Type>("value", dict, p.size())
541  );
542  }
543  else
544  {
545  fvPatchField<Type>::operator=(this->patchInternalField());
546  }
547 
548  commsDir_.expand();
549 
550  if (Pstream::master())
551  {
552  mkDir(baseDir());
553  }
554 
555  if (!initByExternal_)
556  {
557  createLockFile();
558  }
559 
560  // initialise as a fixed value
561  this->refValue() = *this;
562  this->refGrad() = pTraits<Type>::zero;
563  this->valueFraction() = 1.0;
564 }
565 
566 
567 template<class Type>
570 (
572 )
573 :
575  commsDir_(ecmpf.commsDir_),
576  fName_(ecmpf.fName_),
577  waitInterval_(ecmpf.waitInterval_),
578  timeOut_(ecmpf.timeOut_),
579  calcFrequency_(ecmpf.calcFrequency_),
580  initByExternal_(ecmpf.initByExternal_),
581  log_(ecmpf.log_),
582  master_(ecmpf.master_),
583  offsets_(ecmpf.offsets_),
584  initialised_(ecmpf.initialised_),
585  coupledPatchIDs_(ecmpf.coupledPatchIDs_)
586 {}
587 
588 
589 template<class Type>
592 (
595 )
596 :
597  mixedFvPatchField<Type>(ecmpf, iF),
598  commsDir_(ecmpf.commsDir_),
599  fName_(ecmpf.fName_),
600  waitInterval_(ecmpf.waitInterval_),
601  timeOut_(ecmpf.timeOut_),
602  calcFrequency_(ecmpf.calcFrequency_),
603  initByExternal_(ecmpf.initByExternal_),
604  log_(ecmpf.log_),
605  master_(ecmpf.master_),
606  offsets_(ecmpf.offsets_),
607  initialised_(ecmpf.initialised_),
608  coupledPatchIDs_(ecmpf.coupledPatchIDs_)
609 {}
610 
611 
612 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
613 
614 template<class Type>
617 {
618  removeLockFile();
619 }
620 
621 
622 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
623 
624 template<class Type>
626 (
627  const fileName& transferFile
628 )
629 {
630  if (initialised_)
631  {
632  return;
633  }
634 
635  const volFieldType& cvf =
636  static_cast<const volFieldType&>(this->dimensionedInternalField());
637 
638  volFieldType& vf = const_cast<volFieldType&>(cvf);
639 
641 
642  // identify all coupled patches
643  DynamicList<label> coupledPatchIDs(bf.size());
644 
645  forAll(bf, patchI)
646  {
647  if (isA<patchType>(bf[patchI]))
648  {
649  coupledPatchIDs.append(patchI);
650  }
651  }
652 
653  coupledPatchIDs_.transfer(coupledPatchIDs);
654 
655 
656  // initialise by external solver, or just set the master patch
657  if (initByExternal_)
658  {
659  // remove lock file, signalling external source to execute
660 // removeLockFile();
661 
662  forAll(coupledPatchIDs_, i)
663  {
664  label patchI = coupledPatchIDs_[i];
665 
666  patchType& pf = refCast<patchType>(bf[patchI]);
667 
668  pf.setMaster(coupledPatchIDs_);
669  }
670 
671 
672  // wait for initial data to be made available
673  startWait();
674 
675  // read the initial data
676  if (master_)
677  {
678  forAll(coupledPatchIDs_, i)
679  {
680  label patchI = coupledPatchIDs_[i];
681 
682  patchType& pf = refCast<patchType>(bf[patchI]);
683 
684  pf.readData(transferFile);
685  }
686  }
687  }
688  else
689  {
690  setMaster(coupledPatchIDs_);
691  }
692 
693  initialised_ = true;
694 }
695 
696 
697 template<class Type>
699 (
700  const Pstream::commsTypes comms
701 )
702 {
703  if (!initialised_ || this->db().time().timeIndex() % calcFrequency_ == 0)
704  {
705  const fileName transferFile(baseDir()/fName_);
706 
707  // initialise the coupling
708  initialise(transferFile);
709 
710  // write data for external source
711  writeData(transferFile + ".out");
712 
713  // remove lock file, signalling external source to execute
714  removeLockFile();
715 
716  // wait for response
717  startWait();
718 
719  if (master_ && Pstream::master())
720  {
721  // remove old data file from OpenFOAM
722  rm(transferFile + ".out");
723  }
724 
725  // read data passed back from external source
726  readData(transferFile);
727 
728  // create lock file for external source
729  createLockFile();
730  }
731 }
732 
733 
734 template<class Type>
736 (
737  OFstream& os
738 ) const
739 {
740  if (log_)
741  {
742  Info<< type() << ": writing data to " << os.name() << endl;
743  }
744 
745  if (Pstream::parRun())
746  {
747  int tag = Pstream::msgType() + 1;
748 
749  List<Field<scalar> > magSfs(Pstream::nProcs());
750  magSfs[Pstream::myProcNo()].setSize(this->patch().size());
751  magSfs[Pstream::myProcNo()] = this->patch().magSf();
752  Pstream::gatherList(magSfs, tag);
753 
754  List<Field<Type> > values(Pstream::nProcs());
755  values[Pstream::myProcNo()].setSize(this->patch().size());
756  values[Pstream::myProcNo()] = this->refValue();
757  Pstream::gatherList(values, tag);
758 
759  List<Field<Type> > snGrads(Pstream::nProcs());
760  snGrads[Pstream::myProcNo()].setSize(this->patch().size());
761  snGrads[Pstream::myProcNo()] = this->snGrad();
762  Pstream::gatherList(snGrads, tag);
763 
764  if (Pstream::master())
765  {
766  forAll(values, procI)
767  {
768  const Field<scalar>& magSf = magSfs[procI];
769  const Field<Type>& value = values[procI];
770  const Field<Type>& snGrad = snGrads[procI];
771 
772  forAll(magSf, faceI)
773  {
774  os << magSf[faceI] << token::SPACE
775  << value[faceI] << token::SPACE
776  << snGrad[faceI] << nl;
777  }
778  }
779 
780  os.flush();
781  }
782  }
783  else
784  {
785  const Field<scalar>& magSf(this->patch().magSf());
786  const Field<Type>& value(this->refValue());
787  const Field<Type> snGrad(this->snGrad());
788 
789  forAll(magSf, faceI)
790  {
791  os << magSf[faceI] << token::SPACE
792  << value[faceI] << token::SPACE
793  << snGrad[faceI] << nl;
794  }
795 
796  os.flush();
797  }
798 }
799 
800 
801 template<class Type>
803 {
804  const volFieldType& cvf =
805  static_cast<const volFieldType&>(this->dimensionedInternalField());
806 
807  const typename volFieldType::GeometricBoundaryField& bf =
808  cvf.boundaryField();
809 
810  OFstream osPoints(baseDir()/"patchPoints");
811  OFstream osFaces(baseDir()/"patchFaces");
812 
813  if (log_)
814  {
815  Info<< "writing collated patch points to: " << osPoints.name() << nl
816  << "writing collated patch faces to: " << osFaces.name() << endl;
817  }
818 
819  forAll(bf, patchI)
820  {
821  if (isA<patchType>(bf[patchI]))
822  {
823  const patchType& pf = refCast<const patchType>(bf[patchI]);
824 
825  pf.writeGeometry(osPoints, osFaces);
826  }
827  }
828 }
829 
830 
831 template<class Type>
833 {
835 
836  os.writeKeyword("commsDir") << commsDir_ << token::END_STATEMENT << nl;
837  os.writeKeyword("fileName") << fName_ << token::END_STATEMENT << nl;
838  os.writeKeyword("waitInterval") << waitInterval_ << token::END_STATEMENT
839  << nl;
840  os.writeKeyword("timeOut") << timeOut_ << token::END_STATEMENT << nl;
841  os.writeKeyword("calcFrequency") << calcFrequency_ << token::END_STATEMENT
842  << nl;
843  os.writeKeyword("initByExternal") << initByExternal_ << token::END_STATEMENT
844  << nl;
845  os.writeKeyword("log") << log_ << token::END_STATEMENT << nl;
846 
847  this->writeEntry("value", os);
848 }
849 
850 
851 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Output to file stream.
Definition: OFstream.H:81
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
bool readBool(Istream &)
Definition: boolIO.C:60
This boundary condition provides a base class for &#39;mixed&#39; type boundary conditions, i.e. conditions that mix fixed value and patch-normal gradient conditions.
const fileName & name() const
Return the name of the stream.
Definition: IFstream.H:116
label size() const
Return the number of elements in the PtrList.
Definition: PtrListI.H:32
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label timeIndex
Definition: getTimeIndex.H:4
A class for handling words, derived from string.
Definition: word.H:59
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
A class for handling character strings derived from std::string.
Definition: string.H:74
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Foam::fvPatchFieldMapper.
messageStream Info
List< face > faceList
Definition: faceListFwd.H:43
externalCoupledMixedFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &)
Construct from patch and internal field.
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Input from file stream.
Definition: IFstream.H:81
bool rm(const fileName &)
Remove a file, returning true if successful otherwise false.
Definition: POSIX.C:955
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
runTime write()
void writeGeometry() const
Write the geometry to the comms dir.
virtual void readData(const fileName &transferFile)
Read data from external source.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
dictionary dict
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
virtual void flush()
Flush stream.
Definition: OSstream.C:246
volScalarField & p
Definition: createFields.H:51
This boundary condition provides an interface to an external application. Values are transferred as p...
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:118
commsTypes
Types of communications.
Definition: UPstream.H:64
#define forAll(list, i)
Definition: UList.H:421
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Pre-declare SubField and related Field type.
Definition: Field.H:57
errorManip< error > abort(error &err)
Definition: errorManip.H:131
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
const bool writeData(readBool(pdfDictionary.lookup("writeData")))
tmp< pointField > allPoints(const Triangulation &t)
Extract all points in vertex-index order.
rDeltaT dimensionedInternalField()
Generic GeometricField class.
virtual void evaluate(const Pstream::commsTypes commsType=Pstream::blocking)
Evaluate the patch field.
unsigned int sleep(const unsigned int)
Sleep for the specified number of seconds.
Definition: POSIX.C:1055
Traits class for primitives.
Definition: pTraits.H:50
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
List< label > labelList
A List of labels.
Definition: labelList.H:56
A class for handling file names.
Definition: fileName.H:69
Foam::word regionName
virtual void transferData(OFstream &os) const
Transfer data for external source.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual void writeData(const fileName &transferFile) const
Write data for external source - calls transferData.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
bool mkDir(const fileName &, mode_t=0777)
Make a directory and return an error if it could not be created.
Definition: POSIX.C:420
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:300
void inplaceRenumber(const labelUList &oldToNew, ListType &)
Inplace renumber the values of a list.
virtual void writeHeader(OFstream &os) const
Write header to transfer file.