vtkPVblockMesh.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) 2011-2018 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 "vtkPVblockMesh.H"
27 #include "vtkPVblockMeshReader.h"
28 
29 // OpenFOAM includes
30 #include "blockMesh.H"
31 #include "Time.H"
32 #include "patchZones.H"
33 #include "OStringStream.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(vtkPVblockMesh, 0);
47 }
48 
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
52 void Foam::vtkPVblockMesh::resetCounters()
53 {
54  // Reset mesh part ids and sizes
55  arrayRangeBlocks_.reset();
56  arrayRangeEdges_.reset();
57  arrayRangeCorners_.reset();
58 }
59 
60 
61 void Foam::vtkPVblockMesh::updateInfoBlocks
62 (
63  vtkDataArraySelection* arraySelection
64 )
65 {
66  if (debug)
67  {
68  Info<< "<beg> Foam::vtkPVblockMesh::updateInfoBlocks"
69  << " [meshPtr=" << (meshPtr_ ? "set" : "nullptr") << "]" << endl;
70  }
71 
72  arrayRangeBlocks_.reset( arraySelection->GetNumberOfArrays() );
73 
74  const blockMesh& blkMesh = *meshPtr_;
75 
76  const int nBlocks = blkMesh.size();
77  for (int blockI = 0; blockI < nBlocks; ++blockI)
78  {
79  const blockDescriptor& blockDef = blkMesh[blockI];
80 
81  // Display either blockI as a number or with its name
82  // (looked up from blockMeshDict)
83  OStringStream os;
84  blockDescriptor::write(os, blockI, blkMesh.meshDict());
85  word partName(os.str());
86 
87  // append the (optional) zone name
88  if (!blockDef.zoneName().empty())
89  {
90  partName += " - " + blockDef.zoneName();
91  }
92 
93  // Add blockId and zoneName to GUI list
94  arraySelection->AddArray(partName.c_str());
95  }
96 
97  arrayRangeBlocks_ += nBlocks;
98 
99  if (debug)
100  {
101  // just for debug info
102  getSelectedArrayEntries(arraySelection);
103 
104  Info<< "<end> Foam::vtkPVblockMesh::updateInfoBlocks" << endl;
105  }
106 }
107 
108 
109 void Foam::vtkPVblockMesh::updateInfoEdges
110 (
111  vtkDataArraySelection* arraySelection
112 )
113 {
114  if (debug)
115  {
116  Info<< "<beg> Foam::vtkPVblockMesh::updateInfoEdges"
117  << " [meshPtr=" << (meshPtr_ ? "set" : "nullptr") << "]" << endl;
118  }
119 
120  arrayRangeEdges_.reset( arraySelection->GetNumberOfArrays() );
121 
122  const blockMesh& blkMesh = *meshPtr_;
123  const blockEdgeList& edges = blkMesh.edges();
124 
125  const int nEdges = edges.size();
126  forAll(edges, edgeI)
127  {
128  OStringStream ostr;
129  blockVertex::write(ostr, edges[edgeI].start(), blkMesh.meshDict());
130  ostr<< ":";
131  blockVertex::write(ostr, edges[edgeI].end(), blkMesh.meshDict());
132  ostr << " - " << edges[edgeI].type();
133 
134  // Add "beg:end - type" to GUI list
135  arraySelection->AddArray(ostr.str().c_str());
136  }
137 
138  arrayRangeEdges_ += nEdges;
139 
140  if (debug)
141  {
142  // just for debug info
143  getSelectedArrayEntries(arraySelection);
144 
145  Info<< "<end> Foam::vtkPVblockMesh::updateInfoEdges" << endl;
146  }
147 }
148 
149 
150 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
151 
152 Foam::vtkPVblockMesh::vtkPVblockMesh
153 (
154  const char* const FileName,
155  vtkPVblockMeshReader* reader
156 )
157 :
158  reader_(reader),
159  dbPtr_(nullptr),
160  meshPtr_(nullptr),
161  meshRegion_(polyMesh::defaultRegion),
162  meshDir_(polyMesh::meshSubDir),
163  arrayRangeBlocks_("block"),
164  arrayRangeEdges_("edges"),
165  arrayRangeCorners_("corners")
166 {
167  if (debug)
168  {
169  Info<< "Foam::vtkPVblockMesh::vtkPVblockMesh - "
170  << FileName << endl;
171  }
172 
173  // avoid argList and get rootPath/caseName directly from the file
174  fileName fullCasePath(fileName(FileName).path());
175 
176  if (!isDir(fullCasePath))
177  {
178  return;
179  }
180  if (fullCasePath == ".")
181  {
182  fullCasePath = cwd();
183  }
184 
185  // Set the case as an environment variable - some BCs might use this
186  if (fullCasePath.name().find("processor", 0) == 0)
187  {
188  const fileName globalCase = fullCasePath.path();
189 
190  setEnv("FOAM_CASE", globalCase, true);
191  setEnv("FOAM_CASENAME", globalCase.name(), true);
192  }
193  else
194  {
195  setEnv("FOAM_CASE", fullCasePath, true);
196  setEnv("FOAM_CASENAME", fullCasePath.name(), true);
197  }
198 
199  // look for 'case{region}.OpenFOAM'
200  // could be stringent and insist the prefix match the directory name...
201  // Note: cannot use fileName::name() due to the embedded '{}'
202  string caseName(fileName(FileName).lessExt());
203  string::size_type beg = caseName.find_last_of("/{");
204  string::size_type end = caseName.find('}', beg);
205 
206  if
207  (
208  beg != string::npos && caseName[beg] == '{'
209  && end != string::npos && end == caseName.size()-1
210  )
211  {
212  meshRegion_ = caseName.substr(beg+1, end-beg-1);
213 
214  // some safety
215  if (meshRegion_.empty())
216  {
217  meshRegion_ = polyMesh::defaultRegion;
218  }
219 
220  if (meshRegion_ != polyMesh::defaultRegion)
221  {
222  meshDir_ = meshRegion_/polyMesh::meshSubDir;
223  }
224  }
225 
226  if (debug)
227  {
228  Info<< "fullCasePath=" << fullCasePath << nl
229  << "FOAM_CASE=" << getEnv("FOAM_CASE") << nl
230  << "FOAM_CASENAME=" << getEnv("FOAM_CASENAME") << endl;
231  }
232 
233  // Create time object
234  dbPtr_.reset
235  (
236  new Time
237  (
238  Time::controlDictName,
239  fileName(fullCasePath.path()),
240  fileName(fullCasePath.name())
241  )
242  );
243 
244  dbPtr_().functionObjects().off();
245 
246  updateInfo();
247 }
248 
249 
250 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
251 
253 {
254  if (debug)
255  {
256  Info<< "<end> Foam::vtkPVblockMesh::~vtkPVblockMesh" << endl;
257  }
258 
259  // Hmm. pointNumberTextActors are not getting removed
260  //
261  forAll(pointNumberTextActorsPtrs_, pointi)
262  {
263  pointNumberTextActorsPtrs_[pointi]->Delete();
264  }
265  pointNumberTextActorsPtrs_.clear();
266 
267  delete meshPtr_;
268 }
269 
270 
271 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
272 
274 {
275  if (debug)
276  {
277  Info<< "<beg> Foam::vtkPVblockMesh::updateInfo"
278  << " [meshPtr=" << (meshPtr_ ? "set" : "nullptr") << "] " << endl;
279  }
280 
281  resetCounters();
282 
283  vtkDataArraySelection* blockSelection = reader_->GetBlockSelection();
284  vtkDataArraySelection* edgeSelection = reader_->GetCurvedEdgesSelection();
285 
286  // enable 'internalMesh' on the first call
287  // or preserve the enabled selections
288  stringList enabledParts;
289  stringList enabledEdges;
290  bool firstTime = false;
291  if (!blockSelection->GetNumberOfArrays() && !meshPtr_)
292  {
293  firstTime = true;
294  }
295  else
296  {
297  enabledParts = getSelectedArrayEntries(blockSelection);
298  enabledEdges = getSelectedArrayEntries(edgeSelection);
299  }
300 
301  // Clear current mesh parts list
302  blockSelection->RemoveAllArrays();
303  edgeSelection->RemoveAllArrays();
304 
305  // need a blockMesh
306  updateFoamMesh();
307 
308  // Update mesh parts list
309  updateInfoBlocks( blockSelection );
310 
311  // Update curved edges list
312  updateInfoEdges( edgeSelection );
313 
314  // restore the enabled selections
315  if (!firstTime)
316  {
317  setSelectedArrayEntries(blockSelection, enabledParts);
318  setSelectedArrayEntries(edgeSelection, enabledEdges);
319  }
320 
321  if (debug)
322  {
323  Info<< "<end> Foam::vtkPVblockMesh::updateInfo" << endl;
324  }
325 }
326 
327 
328 void Foam::vtkPVblockMesh::updateFoamMesh()
329 {
330  if (debug)
331  {
332  Info<< "<beg> Foam::vtkPVblockMesh::updateFoamMesh" << endl;
333  }
334 
335  // Check to see if the OpenFOAM mesh has been created
336  if (!meshPtr_)
337  {
338  if (debug)
339  {
340  Info<< "Creating blockMesh at time=" << dbPtr_().timeName()
341  << endl;
342  }
343 
344  // Set path for the blockMeshDict
345  const word dictName("blockMeshDict");
346  fileName dictPath(dbPtr_().system()/dictName);
347 
348  // Check if dictionary is present in the constant directory
349  if
350  (
351  exists
352  (
353  dbPtr_().path()/dbPtr_().constant()
355  )
356  )
357  {
358  dictPath = dbPtr_().constant()/polyMesh::meshSubDir/dictName;
359  }
360 
361  // Store dictionary since is used as database inside blockMesh class
362  // for names of vertices and blocks
363  IOdictionary* meshDictPtr = new IOdictionary
364  (
365  IOobject
366  (
367  dictPath,
368  dbPtr_(),
371  true
372  )
373  );
374  meshDictPtr->store();
375 
376  meshPtr_ = new blockMesh(*meshDictPtr, meshRegion_);
377  }
378 
379 
380  if (debug)
381  {
382  Info<< "<end> Foam::vtkPVblockMesh::updateFoamMesh" << endl;
383  }
384 }
385 
386 
388 (
389  vtkMultiBlockDataSet* output
390 )
391 {
392  reader_->UpdateProgress(0.1);
393 
394  // Set up mesh parts selection(s)
395  updateBoolListStatus(blockStatus_, reader_->GetBlockSelection());
396 
397  // Set up curved edges selection(s)
398  updateBoolListStatus(edgeStatus_, reader_->GetCurvedEdgesSelection());
399 
400  reader_->UpdateProgress(0.2);
401 
402  // Update the OpenFOAM mesh
403  updateFoamMesh();
404  reader_->UpdateProgress(0.5);
405 
406  // Convert mesh element
407  int blockNo = 0;
408 
409  convertMeshCorners(output, blockNo);
410  convertMeshBlocks(output, blockNo);
411  convertMeshEdges(output, blockNo);
412 
413  reader_->UpdateProgress(0.8);
414 
415 }
416 
417 
419 {
420  reader_->UpdateProgress(1.0);
421 }
422 
423 
425 (
426  vtkRenderer* renderer,
427  const bool show
428 )
429 {
430  // always remove old actors first
431 
432  forAll(pointNumberTextActorsPtrs_, pointi)
433  {
434  renderer->RemoveViewProp(pointNumberTextActorsPtrs_[pointi]);
435  pointNumberTextActorsPtrs_[pointi]->Delete();
436  }
437  pointNumberTextActorsPtrs_.clear();
438 
439  if (show && meshPtr_)
440  {
441  const blockMesh& blkMesh = *meshPtr_;
442  const pointField& cornerPts = blkMesh.vertices();
443  const scalar scaleFactor = blkMesh.scaleFactor();
444 
445  pointNumberTextActorsPtrs_.setSize(cornerPts.size());
446  forAll(cornerPts, pointi)
447  {
448  vtkTextActor* txt = vtkTextActor::New();
449 
450  // Display either pointi as a number or with its name
451  // (looked up from blockMeshDict)
452  {
453  OStringStream os;
454  blockVertex::write(os, pointi, blkMesh.meshDict());
455  txt->SetInput(os.str().c_str());
456  }
457 
458  // Set text properties
459  vtkTextProperty* tprop = txt->GetTextProperty();
460  tprop->SetFontFamilyToArial();
461  tprop->BoldOn();
462  tprop->ShadowOff();
463  tprop->SetLineSpacing(1.0);
464  tprop->SetFontSize(14);
465  tprop->SetColor(1.0, 0.0, 1.0);
466  tprop->SetJustificationToCentered();
467 
468  // Set text to use 3-D world co-ordinates
469  txt->GetPositionCoordinate()->SetCoordinateSystemToWorld();
470 
471  txt->GetPositionCoordinate()->SetValue
472  (
473  cornerPts[pointi].x()*scaleFactor,
474  cornerPts[pointi].y()*scaleFactor,
475  cornerPts[pointi].z()*scaleFactor
476  );
477 
478  // Add text to each renderer
479  renderer->AddViewProp(txt);
480 
481  // Maintain a list of text labels added so that they can be
482  // removed later
483  pointNumberTextActorsPtrs_[pointi] = txt;
484  }
485  }
486 }
487 
488 
489 
490 void Foam::vtkPVblockMesh::PrintSelf(ostream& os, vtkIndent indent) const
491 {
492 #if 0
493  os << indent << "Number of nodes: "
494  << (meshPtr_ ? meshPtr_->nPoints() : 0) << "\n";
495 
496  os << indent << "Number of cells: "
497  << (meshPtr_ ? meshPtr_->nCells() : 0) << "\n";
498 
499  os << indent << "Number of available time steps: "
500  << (dbPtr_.valid() ? dbPtr_().times().size() : 0) << endl;
501 #endif
502 }
503 
504 // ************************************************************************* //
string getEnv(const word &)
Return environment variable of given name.
Definition: POSIX.C:96
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
void PrintSelf(ostream &, vtkIndent) const
Debug information.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:312
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
static void write(Ostream &, const label blocki, const dictionary &)
Write block index with dictionary lookup.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
void renderPointNumbers(vtkRenderer *, const bool show)
Add/remove point numbers to/from the view.
scalar y
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
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:528
fileName dictPath
const word dictName("particleTrackDict")
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:73
void Update(vtkMultiBlockDataSet *output)
static const char nl
Definition: Ostream.H:265
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
defineTypeNameAndDebug(combustionModel, 0)
~vtkPVblockMesh()
Destructor.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
List< string > stringList
A List of strings.
Definition: stringList.H:50
void setSize(const label)
Reset size of List.
Definition: List.C:281
void CleanUp()
Clean any storage.
PtrList< blockEdge > blockEdgeList
A PtrList of blockEdges.
Definition: blockEdgeList.H:45
void updateInfo()
Update.
fileName cwd()
Return current working directory path name.
Definition: POSIX.C:240
messageStream Info
static void write(Ostream &, const label, const dictionary &)
Write vertex index with optional name backsubstitution.
Definition: blockVertex.C:120
bool setEnv(const word &name, const std::string &value, const bool overwrite)
Set an environment variable.
Definition: POSIX.C:114
bool exists(const fileName &, const bool checkGzip=true, const bool followLink=true)
Does the name exist (as DIRECTORY or FILE) in the file system?
Definition: POSIX.C:509
Namespace for OpenFOAM.
int system(const std::string &command)
Execute the specified command.
Definition: POSIX.C:1217