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