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