vtkPVFoam.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-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 
26 #include "vtkPVFoam.H"
27 #include "vtkPVFoamReader.h"
28 
29 // OpenFOAM includes
30 #include "fvMesh.H"
31 #include "Time.H"
32 #include "patchZones.H"
33 #include "collatedFileOperation.H"
34 
35 // VTK includes
36 #include "vtkDataArraySelection.h"
37 #include "vtkMultiBlockDataSet.h"
38 #include "vtkRenderer.h"
39 #include "vtkTextActor.h"
40 #include "vtkTextProperty.h"
41 
42 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 defineTypeNameAndDebug(vtkPVFoam, 0);
47 }
48 
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
54 
55 void Foam::vtkPVFoam::resetCounters()
56 {
57  // Reset array range information (ids and sizes)
58  arrayRangeVolume_.reset();
59  arrayRangePatches_.reset();
60  arrayRangeLagrangian_.reset();
61  arrayRangeCellZones_.reset();
62  arrayRangeFaceZones_.reset();
63  arrayRangePointZones_.reset();
64  arrayRangeCellSets_.reset();
65  arrayRangeFaceSets_.reset();
66  arrayRangePointSets_.reset();
67 }
68 
69 
70 void Foam::vtkPVFoam::reduceMemory()
71 {
72  forAll(regionPolyDecomp_, i)
73  {
74  regionPolyDecomp_[i].clear();
75  }
76 
77  forAll(zonePolyDecomp_, i)
78  {
79  zonePolyDecomp_[i].clear();
80  }
81 
82  forAll(csetPolyDecomp_, i)
83  {
84  csetPolyDecomp_[i].clear();
85  }
86 
87  if (!reader_->GetCacheMesh())
88  {
89  delete meshPtr_;
90  meshPtr_ = nullptr;
91  }
92 }
93 
94 
95 int Foam::vtkPVFoam::setTime(int nRequest, const double requestTimes[])
96 {
97  Time& runTime = dbPtr_();
98 
99  // Get times list
100  instantList Times = runTime.times();
101 
102  int nearestIndex = timeIndex_;
103  for (int requestI = 0; requestI < nRequest; ++requestI)
104  {
105  int index = Time::findClosestTimeIndex(Times, requestTimes[requestI]);
106  if (index >= 0 && index != timeIndex_)
107  {
108  nearestIndex = index;
109  break;
110  }
111  }
112 
113  if (nearestIndex < 0)
114  {
115  nearestIndex = 0;
116  }
117 
118  if (debug)
119  {
120  Info<< "<beg> Foam::vtkPVFoam::setTime(";
121  for (int requestI = 0; requestI < nRequest; ++requestI)
122  {
123  if (requestI)
124  {
125  Info<< ", ";
126  }
127 
128  Info<< requestTimes[requestI];
129  }
130  Info<< ") - previousIndex = " << timeIndex_
131  << ", nearestIndex = " << nearestIndex << endl;
132  }
133 
134 
135  // see what has changed
136  if (timeIndex_ != nearestIndex)
137  {
138  timeIndex_ = nearestIndex;
139  runTime.setTime(Times[nearestIndex], nearestIndex);
140 
141  // the fields change each time
142  fieldsChanged_ = true;
143 
144  if (meshPtr_)
145  {
146  if (meshPtr_->readUpdate() != polyMesh::UNCHANGED)
147  {
148  meshChanged_ = true;
149  }
150  }
151  else
152  {
153  meshChanged_ = true;
154  }
155  }
156 
157  if (debug)
158  {
159  Info<< "<end> Foam::vtkPVFoam::setTime() - selectedTime="
160  << Times[nearestIndex].name() << " index=" << timeIndex_
161  << "/" << Times.size()
162  << " meshChanged=" << Switch(meshChanged_)
163  << " fieldsChanged=" << Switch(fieldsChanged_) << endl;
164  }
165 
166  return nearestIndex;
167 }
168 
169 
170 void Foam::vtkPVFoam::updateMeshPartsStatus()
171 {
172  if (debug)
173  {
174  Info<< "<beg> Foam::vtkPVFoam::updateMeshPartsStatus" << endl;
175  }
176 
177  vtkDataArraySelection* selection = reader_->GetPartSelection();
178  label nElem = selection->GetNumberOfArrays();
179 
180  if (partStatus_.size() != nElem)
181  {
182  partStatus_.setSize(nElem);
183  partStatus_ = false;
184  meshChanged_ = true;
185  }
186 
187  // this needs fixing if we wish to re-use the datasets
188  partDataset_.setSize(nElem);
189  partDataset_ = -1;
190 
191  // Read the selected mesh parts (zones, patches ...) and add to list
192  forAll(partStatus_, partId)
193  {
194  const int setting = selection->GetArraySetting(partId);
195 
196  if (partStatus_[partId] != setting)
197  {
198  partStatus_[partId] = setting;
199  meshChanged_ = true;
200  }
201 
202  if (debug)
203  {
204  Info<< " part[" << partId << "] = "
205  << partStatus_[partId]
206  << " : " << selection->GetArrayName(partId) << endl;
207  }
208  }
209  if (debug)
210  {
211  Info<< "<end> Foam::vtkPVFoam::updateMeshPartsStatus" << endl;
212  }
213 }
214 
215 
216 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
217 
218 Foam::vtkPVFoam::vtkPVFoam
219 (
220  const char* const FileName,
221  vtkPVFoamReader* reader
222 )
223 :
224  reader_(reader),
225  dbPtr_(nullptr),
226  meshPtr_(nullptr),
227  meshRegion_(polyMesh::defaultRegion),
228  meshDir_(polyMesh::meshSubDir),
229  timeIndex_(-1),
230  meshChanged_(true),
231  fieldsChanged_(true),
232  arrayRangeVolume_("unzoned"),
233  arrayRangePatches_("patches"),
234  arrayRangeLagrangian_("lagrangian"),
235  arrayRangeCellZones_("cellZone"),
236  arrayRangeFaceZones_("faceZone"),
237  arrayRangePointZones_("pointZone"),
238  arrayRangeCellSets_("cellSet"),
239  arrayRangeFaceSets_("faceSet"),
240  arrayRangePointSets_("pointSet")
241 {
242  if (debug)
243  {
244  Info<< "Foam::vtkPVFoam::vtkPVFoam - " << FileName << endl;
245  printMemory();
246  }
247 
248  // Make sure not to use the threaded version - it does not like
249  // being loaded as a shared library - static cleanup order is problematic.
250  // For now just disable the threaded writer.
251  fileOperations::collatedFileOperation::maxThreadFileBufferSize = 0;
252 
253 
254  // avoid argList and get rootPath/caseName directly from the file
255  fileName fullCasePath(fileName(FileName).path());
256 
257  if (!isDir(fullCasePath))
258  {
259  return;
260  }
261  if (fullCasePath == ".")
262  {
263  fullCasePath = cwd();
264  }
265 
266  // Set the case as an environment variable - some BCs might use this
267  if (fullCasePath.name().find("processor", 0) == 0)
268  {
269  const fileName globalCase = fullCasePath.path();
270 
271  setEnv("FOAM_CASE", globalCase, true);
272  setEnv("FOAM_CASENAME", globalCase.name(), true);
273  }
274  else
275  {
276  setEnv("FOAM_CASE", fullCasePath, true);
277  setEnv("FOAM_CASENAME", fullCasePath.name(), true);
278  }
279 
280  // look for 'case{region}.OpenFOAM'
281  // could be stringent and insist the prefix match the directory name...
282  // Note: cannot use fileName::name() due to the embedded '{}'
283  string caseName(fileName(FileName).lessExt());
284  string::size_type beg = caseName.find_last_of("/{");
285  string::size_type end = caseName.find('}', beg);
286 
287  if
288  (
289  beg != string::npos && caseName[beg] == '{'
290  && end != string::npos && end == caseName.size()-1
291  )
292  {
293  meshRegion_ = caseName.substr(beg+1, end-beg-1);
294 
295  // some safety
296  if (meshRegion_.empty())
297  {
298  meshRegion_ = polyMesh::defaultRegion;
299  }
300 
301  if (meshRegion_ != polyMesh::defaultRegion)
302  {
303  meshDir_ = meshRegion_/polyMesh::meshSubDir;
304  }
305  }
306 
307  if (debug)
308  {
309  Info<< "fullCasePath=" << fullCasePath << nl
310  << "FOAM_CASE=" << getEnv("FOAM_CASE") << nl
311  << "FOAM_CASENAME=" << getEnv("FOAM_CASENAME") << nl
312  << "region=" << meshRegion_ << endl;
313  }
314 
315  // Create time object
316  dbPtr_.reset
317  (
318  new Time
319  (
320  Time::controlDictName,
321  fileName(fullCasePath.path()),
322  fileName(fullCasePath.name())
323  )
324  );
325 
326  dbPtr_().functionObjects().off();
327 
328  updateInfo();
329 }
330 
331 
332 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
333 
335 {
336  if (debug)
337  {
338  Info<< "<end> Foam::vtkPVFoam::~vtkPVFoam" << endl;
339  }
340 
341  delete meshPtr_;
342 }
343 
344 
345 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
346 
348 {
349  if (debug)
350  {
351  Info<< "<beg> Foam::vtkPVFoam::updateInfo"
352  << " [meshPtr=" << (meshPtr_ ? "set" : "nullptr") << "] timeIndex="
353  << timeIndex_ << endl;
354  }
355 
356  resetCounters();
357 
358  vtkDataArraySelection* partSelection = reader_->GetPartSelection();
359 
360  // there are two ways to ensure we have the correct list of parts:
361  // 1. remove everything and then set particular entries 'on'
362  // 2. build a 'char **' list and call SetArraysWithDefault()
363  //
364  // Nr. 2 has the potential advantage of not touching the modification
365  // time of the vtkDataArraySelection, but the qt/paraview proxy
366  // layer doesn't care about that anyhow.
367 
368  // enable 'internalMesh' on the first call
369  // or preserve the enabled selections
370  stringList enabledEntries;
371  if (!partSelection->GetNumberOfArrays() && !meshPtr_)
372  {
373  enabledEntries.setSize(1);
374  enabledEntries[0] = "internalMesh";
375  }
376  else
377  {
378  enabledEntries = getSelectedArrayEntries(partSelection);
379  }
380 
381  // Clear current mesh parts list
382  partSelection->RemoveAllArrays();
383 
384  // Update mesh parts list - add Lagrangian at the bottom
385  updateInfoInternalMesh(partSelection);
386  updateInfoPatches(partSelection, enabledEntries);
387  updateInfoSets(partSelection);
388  updateInfoZones(partSelection);
389  updateInfoLagrangian(partSelection);
390 
391  // restore the enabled selections
392  setSelectedArrayEntries(partSelection, enabledEntries);
393 
394  if (meshChanged_)
395  {
396  fieldsChanged_ = true;
397  }
398 
399  // Update volume, point and lagrangian fields
400  updateInfoFields<fvPatchField, volMesh>
401  (
402  reader_->GetVolFieldSelection()
403  );
404  updateInfoFields<pointPatchField, pointMesh>
405  (
406  reader_->GetPointFieldSelection()
407  );
408  updateInfoLagrangianFields();
409 
410  if (debug)
411  {
412  // just for debug info
413  getSelectedArrayEntries(partSelection);
414  Info<< "<end> Foam::vtkPVFoam::updateInfo" << endl;
415  }
416 
417 }
418 
419 
420 void Foam::vtkPVFoam::updateFoamMesh()
421 {
422  if (debug)
423  {
424  Info<< "<beg> Foam::vtkPVFoam::updateFoamMesh" << endl;
425  printMemory();
426  }
427 
428  if (!reader_->GetCacheMesh())
429  {
430  delete meshPtr_;
431  meshPtr_ = nullptr;
432  }
433 
434  // Check to see if the OpenFOAM mesh has been created
435  if (!meshPtr_)
436  {
437  if (debug)
438  {
439  Info<< "Creating OpenFOAM mesh for region " << meshRegion_
440  << " at time=" << dbPtr_().timeName()
441  << endl;
442 
443  }
444 
445  meshPtr_ = new fvMesh
446  (
447  IOobject
448  (
449  meshRegion_,
450  dbPtr_().timeName(),
451  dbPtr_(),
453  )
454  );
455 
456  meshChanged_ = true;
457  }
458  else
459  {
460  if (debug)
461  {
462  Info<< "Using existing OpenFOAM mesh" << endl;
463  }
464  }
465 
466  if (debug)
467  {
468  Info<< "<end> Foam::vtkPVFoam::updateFoamMesh" << endl;
469  printMemory();
470  }
471 }
472 
473 
475 (
476  vtkMultiBlockDataSet* output,
477  vtkMultiBlockDataSet* lagrangianOutput
478 )
479 {
480  if (debug)
481  {
482  cout<< "<beg> Foam::vtkPVFoam::Update - output with "
483  << output->GetNumberOfBlocks() << " and "
484  << lagrangianOutput->GetNumberOfBlocks() << " blocks\n";
485  output->Print(cout);
486  lagrangianOutput->Print(cout);
487  printMemory();
488  }
489  reader_->UpdateProgress(0.1);
490 
491  // Set up mesh parts selection(s)
492  updateMeshPartsStatus();
493 
494  reader_->UpdateProgress(0.15);
495 
496  // Update the OpenFOAM mesh
497  updateFoamMesh();
498  reader_->UpdateProgress(0.4);
499 
500  // Convert meshes - start port0 at block=0
501  int blockNo = 0;
502 
503  convertMeshVolume(output, blockNo);
504  convertMeshPatches(output, blockNo);
505  reader_->UpdateProgress(0.6);
506 
507  if (reader_->GetIncludeZones())
508  {
509  convertMeshCellZones(output, blockNo);
510  convertMeshFaceZones(output, blockNo);
511  convertMeshPointZones(output, blockNo);
512  reader_->UpdateProgress(0.65);
513  }
514 
515  if (reader_->GetIncludeSets())
516  {
517  convertMeshCellSets(output, blockNo);
518  convertMeshFaceSets(output, blockNo);
519  convertMeshPointSets(output, blockNo);
520  reader_->UpdateProgress(0.7);
521  }
522 
523  convertMeshLagrangian(lagrangianOutput, blockNo);
524 
525  reader_->UpdateProgress(0.8);
526 
527  // Update fields
528  convertVolFields(output);
529  convertPointFields(output);
530  convertLagrangianFields(lagrangianOutput);
531  if (debug)
532  {
533  Info<< "done reader part" << endl;
534  }
535  reader_->UpdateProgress(0.95);
536 
537  meshChanged_ = fieldsChanged_ = false;
538 }
539 
540 
542 {
543  // reclaim some memory
544  reduceMemory();
545  reader_->UpdateProgress(1.0);
546 }
547 
548 
549 double* Foam::vtkPVFoam::findTimes(int& nTimeSteps)
550 {
551  int nTimes = 0;
552  double* tsteps = nullptr;
553 
554  if (dbPtr_.valid())
555  {
556  Time& runTime = dbPtr_();
557  instantList timeLst = runTime.times();
558 
559  // find the first time for which this mesh appears to exist
560  label timeI = 0;
561  for (; timeI < timeLst.size(); ++timeI)
562  {
563  const word& timeName = timeLst[timeI].name();
564 
565  if
566  (
567  IOobject
568  (
569  "points",
570  timeName,
571  meshDir_,
572  runTime
573  ).typeHeaderOk<pointIOField>(true)
574  )
575  {
576  break;
577  }
578  }
579 
580  nTimes = timeLst.size() - timeI;
581 
582  // skip "constant" time whenever possible
583  if (timeI == 0 && nTimes > 1)
584  {
585  if (timeLst[timeI].name() == runTime.constant())
586  {
587  ++timeI;
588  --nTimes;
589  }
590  }
591 
592 
593  // skip "0/" time if requested and possible
594  if (nTimes > 1 && reader_->GetSkipZeroTime())
595  {
596  if (mag(timeLst[timeI].value()) < SMALL)
597  {
598  ++timeI;
599  --nTimes;
600  }
601  }
602 
603  if (nTimes)
604  {
605  tsteps = new double[nTimes];
606  for (label stepI = 0; stepI < nTimes; ++stepI, ++timeI)
607  {
608  tsteps[stepI] = timeLst[timeI].value();
609  }
610  }
611  }
612  else
613  {
614  if (debug)
615  {
616  cout<< "no valid dbPtr:\n";
617  }
618  }
619 
620  // vector length returned via the parameter
621  nTimeSteps = nTimes;
622 
623  return tsteps;
624 }
625 
626 
628 (
629  vtkRenderer* renderer,
630  const bool show
631 )
632 {
633  if (!meshPtr_)
634  {
635  return;
636  }
637 
638  // always remove old actors first
639 
640  forAll(patchTextActorsPtrs_, patchi)
641  {
642  renderer->RemoveViewProp(patchTextActorsPtrs_[patchi]);
643  patchTextActorsPtrs_[patchi]->Delete();
644  }
645  patchTextActorsPtrs_.clear();
646 
647  if (show)
648  {
649  // get the display patches, strip off any suffix
650  wordHashSet selectedPatches = getSelected
651  (
652  reader_->GetPartSelection(),
653  arrayRangePatches_
654  );
655 
656  if (selectedPatches.empty())
657  {
658  return;
659  }
660 
661  const polyBoundaryMesh& pbMesh = meshPtr_->boundaryMesh();
662 
663  // Find the total number of zones
664  // Each zone will take the patch name
665  // Number of zones per patch ... zero zones should be skipped
666  labelList nZones(pbMesh.size(), 0);
667 
668  // Per global zone number the average face centre position
669  List<DynamicList<point>> zoneCentre(pbMesh.size());
670 
671 
672  // Loop through all patches to determine zones, and centre of each zone
673  forAll(pbMesh, patchi)
674  {
675  const polyPatch& pp = pbMesh[patchi];
676 
677  // Only include the patch if it is selected
678  if (!selectedPatches.found(pp.name()))
679  {
680  continue;
681  }
682 
683  const labelListList& edgeFaces = pp.edgeFaces();
684  const vectorField& n = pp.faceNormals();
685 
686  boolList featEdge(pp.nEdges(), false);
687 
688  forAll(edgeFaces, edgeI)
689  {
690  const labelList& eFaces = edgeFaces[edgeI];
691 
692  if (eFaces.size() == 1)
693  {
694  // Note: could also do ones with > 2 faces but this gives
695  // too many zones for baffles
696  featEdge[edgeI] = true;
697  }
698  else if (mag(n[eFaces[0]] & n[eFaces[1]]) < 0.5)
699  {
700  featEdge[edgeI] = true;
701  }
702  }
703 
704  // Do topological analysis of patch, find disconnected regions
705  patchZones pZones(pp, featEdge);
706 
707  nZones[patchi] = pZones.nZones();
708 
709  labelList zoneNFaces(pZones.nZones(), 0);
710 
711  // Create storage for additional zone centres
712  forAll(zoneNFaces, zoneI)
713  {
714  zoneCentre[patchi].append(Zero);
715  }
716 
717  // Do averaging per individual zone
718  forAll(pp, facei)
719  {
720  label zoneI = pZones[facei];
721  zoneCentre[patchi][zoneI] += pp[facei].centre(pp.points());
722  zoneNFaces[zoneI]++;
723  }
724 
725  forAll(zoneCentre[patchi], zoneI)
726  {
727  zoneCentre[patchi][zoneI] /= zoneNFaces[zoneI];
728  }
729  }
730 
731  // Count number of zones we're actually going to display.
732  // This is truncated to a max per patch
733 
734  const label MAXPATCHZONES = 20;
735 
736  label displayZoneI = 0;
737 
738  forAll(pbMesh, patchi)
739  {
740  displayZoneI += min(MAXPATCHZONES, nZones[patchi]);
741  }
742 
743  if (debug)
744  {
745  Info<< "displayed zone centres = " << displayZoneI << nl
746  << "zones per patch = " << nZones << endl;
747  }
748 
749  // Set the size of the patch labels to max number of zones
750  patchTextActorsPtrs_.setSize(displayZoneI);
751 
752  if (debug)
753  {
754  Info<< "constructing patch labels" << endl;
755  }
756 
757  // Actor index
758  displayZoneI = 0;
759 
760  forAll(pbMesh, patchi)
761  {
762  const polyPatch& pp = pbMesh[patchi];
763 
764  label globalZoneI = 0;
765 
766  // Only selected patches will have a non-zero number of zones
767  label nDisplayZones = min(MAXPATCHZONES, nZones[patchi]);
768  label increment = 1;
769  if (nZones[patchi] >= MAXPATCHZONES)
770  {
771  increment = nZones[patchi]/MAXPATCHZONES;
772  }
773 
774  for (label i = 0; i < nDisplayZones; i++)
775  {
776  if (debug)
777  {
778  Info<< "patch name = " << pp.name() << nl
779  << "anchor = " << zoneCentre[patchi][globalZoneI] << nl
780  << "globalZoneI = " << globalZoneI << endl;
781  }
782 
783  vtkTextActor* txt = vtkTextActor::New();
784 
785  txt->SetInput(pp.name().c_str());
786 
787  // Set text properties
788  vtkTextProperty* tprop = txt->GetTextProperty();
789  tprop->SetFontFamilyToArial();
790  tprop->BoldOff();
791  tprop->ShadowOff();
792  tprop->SetLineSpacing(1.0);
793  tprop->SetFontSize(12);
794  tprop->SetColor(1.0, 0.0, 0.0);
795  tprop->SetJustificationToCentered();
796 
797  // Set text to use 3-D world co-ordinates
798  txt->GetPositionCoordinate()->SetCoordinateSystemToWorld();
799 
800  txt->GetPositionCoordinate()->SetValue
801  (
802  zoneCentre[patchi][globalZoneI].x(),
803  zoneCentre[patchi][globalZoneI].y(),
804  zoneCentre[patchi][globalZoneI].z()
805  );
806 
807  // Add text to each renderer
808  renderer->AddViewProp(txt);
809 
810  // Maintain a list of text labels added so that they can be
811  // removed later
812  patchTextActorsPtrs_[displayZoneI] = txt;
813 
814  globalZoneI += increment;
815  displayZoneI++;
816  }
817  }
818 
819  // Resize the patch names list to the actual number of patch names added
820  patchTextActorsPtrs_.setSize(displayZoneI);
821  }
822 }
823 
824 
825 void Foam::vtkPVFoam::PrintSelf(ostream& os, vtkIndent indent) const
826 {
827  os << indent << "Number of nodes: "
828  << (meshPtr_ ? meshPtr_->nPoints() : 0) << "\n";
829 
830  os << indent << "Number of cells: "
831  << (meshPtr_ ? meshPtr_->nCells() : 0) << "\n";
832 
833  os << indent << "Number of available time steps: "
834  << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << "\n";
835 
836  os << indent << "mesh region: " << meshRegion_ << "\n";
837 }
838 
839 
840 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
string getEnv(const word &)
Return environment variable of given name.
Definition: POSIX.C:104
List< instant > instantList
List of instants.
Definition: instantList.H:42
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
static void printMemory()
Simple memory used debugging information.
double * findTimes(int &nTimeSteps)
Allocate and return a list of selected times.
void updateInfo()
Update.
label nCells() const
void renderPatchNames(vtkRenderer *, const bool show)
Add/remove patch names to/from the view.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void PrintSelf(ostream &, vtkIndent) const
Debug information.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
void Update(vtkMultiBlockDataSet *output, vtkMultiBlockDataSet *lagrangianOutput)
List< bool > boolList
Bool container classes.
Definition: boolList.H:50
scalar y
void CleanUp()
Clean any storage.
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:124
bool isDir(const fileName &, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?
Definition: POSIX.C:536
word timeName
Definition: getTimeIndex.H:3
List< label > labelList
A List of labels.
Definition: labelList.H:56
static const zero Zero
Definition: zero.H:91
HashSet wordHashSet
A HashSet with word keys.
Definition: HashSet.H:207
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:73
static const char nl
Definition: Ostream.H:262
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
List< string > stringList
A List of strings.
Definition: stringList.H:50
void setSize(const label)
Reset size of List.
Definition: List.C:281
label patchi
int setTime(int count, const double requestTimes[])
Set the runTime to the first plausible request time,.
fileName cwd()
Return current working directory path name.
Definition: POSIX.C:248
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
label nPoints() const
Field< vector > vectorField
Specialisation of Field<T> for vector.
~vtkPVFoam()
Destructor.
bool setEnv(const word &name, const std::string &value, const bool overwrite)
Set an environment variable.
Definition: POSIX.C:122
Namespace for OpenFOAM.
static label findClosestTimeIndex(const instantList &, const scalar, const word &constantName="constant")
Search instantList for the time index closest to the given time.
Definition: Time.C:738
IOporosityModelList pZones(mesh)