1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
31 #include "polyBoundaryMeshEntries.H"
32 #include "IOobjectList.H"
35 #include "volFields.H"
36 #include "pointMesh.H"
37 #include "volPointInterpolation.H"
39 #include "vtkFoamReader.h"
40 #include "vtkDataArraySelection.h"
41 #include "vtkUnstructuredGrid.h"
42 #include "vtkPointData.h"
43 #include "vtkCellData.h"
44 #include "vtkFloatArray.h"
45 #include "vtkCharArray.h"
47 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
49 defineTypeNameAndDebug(Foam::vtkFoam, 0);
51 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53 #include "vtkFoamConvertFields.H"
55 void Foam::vtkFoam::SetName
57 vtkUnstructuredGrid* vtkMesh,
61 vtkCharArray* nmArray = vtkCharArray::New();
62 nmArray->SetName("Name");
63 size_t len = strlen(name);
64 nmArray->SetNumberOfTuples(static_cast<vtkIdType>(len)+1);
65 char* copy = nmArray->GetPointer(0);
66 memcpy(copy, name, len);
68 vtkMesh->GetFieldData()->AddArray(nmArray);
73 Foam::string Foam::vtkFoam::padTimeString(const string& ts)
75 return ts + string(" ", max(label(12 - ts.size()), 0));
79 // Pad the patch name string in order to account for dynamic changes
80 // in patch names during topological changes
81 Foam::string Foam::vtkFoam::padPatchString(const string& ps)
83 label n = max(label(50 - ps.size()), 0);
84 return ps + string(" ", n);
88 void Foam::vtkFoam::setSelectedTime
95 instantList Times = runTime.times();
96 int timeIndex = min(max(reader->GetTimeStep() + 1, 0), Times.size()-1);
98 // If this is the first call timeIndex will be 0 ("constant")
99 // so reset to the first time step if one exists and deselect every
100 // element of the selection array
103 timeIndex = min(1, Times.size()-1);
104 reader->GetTimeSelection()->DisableAllArrays();
107 label selectedTimeIndex = -1;
108 label nSelectedTimes = reader->GetTimeSelection()->GetNumberOfArrays();
110 for (label i=nSelectedTimes-1; i>=0; i--)
112 if(reader->GetTimeSelection()->GetArraySetting(i))
114 word timeName = string::validate<word>
116 reader->GetTimeSelection()->GetArrayName(i)
121 if (Times[j].name() == timeName)
123 selectedTimeIndex = j;
131 if (selectedTimeIndex != -1)
133 timeIndex = min(selectedTimeIndex, Times.size()-1);
138 Info<< "Selecting time " << Times[timeIndex].name() << endl;
141 runTime.setTime(Times[timeIndex], timeIndex);
143 Times = runTime.times();
145 reader->SetTimeStepRange(0, max(Times.size()-2, 0));
147 // reset the time steps ...
148 reader->GetTimeSelection()->RemoveAllArrays();
150 int* TimeStepLimits = reader->GetTimeStepLimits();
151 label maxStartTimes = min(Times.size(), TimeStepLimits[0]);
152 label maxNTimes = min(Times.size() - maxStartTimes, TimeStepLimits[1]);
154 for (label i=0; i<maxStartTimes; i++)
156 reader->GetTimeSelection()
157 ->AddArray(padTimeString(Times[i].name()).c_str());
160 if (Times.size() > TimeStepLimits[0] + TimeStepLimits[1])
162 reader->GetTimeSelection()->AddArray(padTimeString("...").c_str());
165 for (label i=Times.size() - maxNTimes; i<Times.size(); i++)
167 reader->GetTimeSelection()
168 ->AddArray(padTimeString(Times[i].name()).c_str());
171 // Disable all the time selections (which are all selected by default) ...
172 reader->GetTimeSelection()->DisableAllArrays();
174 // But maintain the selections made previously
175 if (selectedTimeIndex != -1 && selectedTimeIndex < Times.size())
177 reader->GetTimeSelection()->EnableArray
178 (padTimeString(Times[selectedTimeIndex].name()).c_str());
183 void Foam::vtkFoam::updateSelectedRegions()
187 Info<< "Foam::vtkFoam::updateSelectedRegions()" << endl;
190 label nRegions = reader_->GetRegionSelection()->GetNumberOfArrays();
192 selectedRegions_.setSize(nRegions);
194 // Read the selected patches and add to the region list
195 for (int i=0; i<nRegions; i++)
197 selectedRegions_[i] =
198 reader_->GetRegionSelection()->GetArraySetting(i);
203 void Foam::vtkFoam::convertMesh()
207 Info<< "Foam::vtkFoam::convertMesh()" << endl;
210 const fvMesh& mesh = *meshPtr_;
212 // Read the internal mesh as region 0 if selected
213 if (reader_->GetRegionSelection()->GetArraySetting(0))
215 selectedRegions_[0] = true;
219 vtkUnstructuredGrid::SafeDownCast(reader_->GetOutput(0))
224 selectedRegions_[0] = false;
226 vtkUnstructuredGrid *vtkMesh =
227 vtkUnstructuredGrid::SafeDownCast(reader_->GetOutput(0));
229 vtkMesh->Initialize();
230 SetName(vtkMesh, "(Internal Mesh)");
234 // Read the selected patches and add to the region list
236 polyBoundaryMeshEntries patchEntries
241 dbPtr_().findInstance(polyMesh::meshSubDir, "boundary"),
242 polyMesh::meshSubDir,
251 label regioniLast = 0;
253 // Read in the number Outputs (patch regions) currently being used
254 label currNOutputs = reader_->GetNumberOfOutputs();
256 // Cycle through all the patches in the boundary file for the relevant
258 forAll(patchEntries, entryi)
260 // Number of faces in the current patch (Used to detect dummy patches
262 label nFaces(readLabel(patchEntries[entryi].dict().lookup("nFaces")));
264 // Check to see if the patch is currently a part of the displayed list
267 reader_->GetRegionSelection()->ArrayExists
269 padPatchString(patchEntries[entryi].keyword()).c_str()
275 // Remove patch if it is only a dummy patch in the current
276 // time step with zero faces
277 reader_->GetRegionSelection()->RemoveArrayByName
279 padPatchString(patchEntries[entryi].keyword()).c_str()
284 // A patch already existent in the list and which
285 // continues to exist found
291 // A new patch so far not yet included into the list has been found
296 // Add a new entry to the list of regions
297 reader_->GetRegionSelection()->AddArray
299 padPatchString(patchEntries[entryi].keyword()).c_str()
302 // AddArray automatically enables a new array... disable
304 reader_->GetRegionSelection()->DisableArray
306 padPatchString(patchEntries[entryi].keyword()).c_str()
311 // Avoid Initialization of the same Output twice
312 if (regioni != regioniLast)
314 // Only setup an Output if it has not been setup before
315 if(regioni >= currNOutputs)
317 vtkUnstructuredGrid* ugrid = vtkUnstructuredGrid::New();
318 reader_->SetNthOutput(regioni,ugrid);
321 // Initialize -> Delete memory used, and reset to zero state
322 reader_->GetOutput(regioni)->Initialize();
323 regioniLast = regioni;
327 // Initialize (reset to zero and free) any outputs which are not used
329 if (regioni < currNOutputs)
331 for(label i = (regioni+1); i < currNOutputs;i++)
333 reader_->GetOutput(i)->Initialize();
337 selectedRegions_.setSize(regioni + 1);
341 const polyBoundaryMesh& patches = mesh.boundaryMesh();
343 forAll (patches, patchi)
345 if (patches[patchi].size())
349 if (reader_->GetRegionSelection()->GetArraySetting(regioni))
351 selectedRegions_[regioni] = true;
355 vtkUnstructuredGrid::SafeDownCast
357 reader_->GetOutput(regioni)
363 selectedRegions_[regioni] = false;
365 vtkUnstructuredGrid *vtkMesh =
366 vtkUnstructuredGrid::SafeDownCast
368 reader_->GetOutput(regioni)
371 vtkMesh->Initialize();
375 ('(' + padPatchString(patches[patchi].name()) + ')').c_str()
383 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
385 Foam::vtkFoam::vtkFoam(const char* const FileName, vtkFoamReader* reader)
392 fileName fullCasePath(fileName(FileName).path());
394 if (!isDir(fullCasePath))
399 char* argvStrings[3];
400 argvStrings[0] = new char[9];
401 strcpy(argvStrings[0], "/vtkFoam");
402 argvStrings[1] = new char[6];
403 strcpy(argvStrings[1], "-case");
404 argvStrings[2] = new char[fullCasePath.size()+1];
405 strcpy(argvStrings[2], fullCasePath.c_str());
408 char** argv = &argvStrings[0];
409 argsPtr_.reset(new argList(argc, argv));
411 for(int i = 0; i < argc; i++)
413 delete[] argvStrings[i];
420 Time::controlDictName,
421 argsPtr_().rootPath(),
422 argsPtr_().caseName()
425 dbPtr_().functionObjects().off();
426 setSelectedTime(dbPtr_(), reader_);
430 Info<< "vtkFoam::ExecuteInformation: Initialising outputs" << endl;
433 reader_->GetRegionSelection()->AddArray("Internal Mesh");
435 vtkUnstructuredGrid* ugrid = vtkUnstructuredGrid::New();
436 reader_->SetNthOutput(0, ugrid);
438 reader_->GetOutput(0)->Initialize();
440 polyBoundaryMeshEntries patchEntries
445 dbPtr_().findInstance(polyMesh::meshSubDir, "boundary"),
446 polyMesh::meshSubDir,
455 forAll(patchEntries, entryi)
457 label nFaces(readLabel(patchEntries[entryi].dict().lookup("nFaces")));
463 reader_->GetRegionSelection()->AddArray
465 padPatchString(patchEntries[entryi].keyword()).c_str()
468 vtkUnstructuredGrid* ugrid = vtkUnstructuredGrid::New();
469 reader_->SetNthOutput(regioni, ugrid);
471 reader_->GetOutput(regioni)->Initialize();
475 selectedRegions_.setSize(regioni + 1);
476 selectedRegions_ = true;
482 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
484 Foam::vtkFoam::~vtkFoam()
486 // Do NOT delete meshPtr_ since still referenced somehow.
490 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
492 #include "vtkFoamAddFields.H"
494 void Foam::vtkFoam::UpdateInformation()
498 Info<< "TimeStep = " << reader_->GetTimeStep() << endl;
501 setSelectedTime(dbPtr_(), reader_);
503 // Search for list of objects for this time
504 IOobjectList objects(dbPtr_(), dbPtr_().timeName());
506 addFields<volScalarField>(reader_->GetVolFieldSelection(), objects);
507 addFields<volVectorField>(reader_->GetVolFieldSelection(), objects);
508 addFields<volSphericalTensorField>(reader_->GetVolFieldSelection(), objects);
509 addFields<volSymmTensorField>(reader_->GetVolFieldSelection(), objects);
510 addFields<volTensorField>(reader_->GetVolFieldSelection(), objects);
512 addFields<pointScalarField>(reader_->GetPointFieldSelection(), objects);
513 addFields<pointVectorField>(reader_->GetPointFieldSelection(), objects);
514 addFields<pointSphericalTensorField>(reader_->GetPointFieldSelection(), objects);
515 addFields<pointSymmTensorField>(reader_->GetPointFieldSelection(), objects);
516 addFields<pointTensorField>(reader_->GetPointFieldSelection(), objects);
520 void Foam::vtkFoam::Update()
524 !reader_->GetCacheMesh()
525 || reader_->GetTimeSelection()->GetArraySetting(0)
531 // Clear the current set of selected fields
533 for (label i=0; i<reader_->GetNumberOfOutputs(); i++)
535 vtkUnstructuredGrid *vtkMesh =
536 vtkUnstructuredGrid::SafeDownCast(reader_->GetOutput(i));
538 vtkCellData* cellData = vtkMesh->GetCellData();
539 int numberOfCellArrays = cellData->GetNumberOfArrays();
541 wordList cellFieldNames(numberOfCellArrays);
542 for (int j=0; j<numberOfCellArrays; j++)
544 cellFieldNames[j] = cellData->GetArrayName(j);
547 for (int j=0; j<numberOfCellArrays; j++)
549 cellData->RemoveArray(cellFieldNames[j].c_str());
552 vtkPointData* pointData = vtkMesh->GetPointData();
553 int numberOfPointArrays = pointData->GetNumberOfArrays();
555 wordList pointFieldNames(numberOfPointArrays);
556 for (int j=0; j<numberOfPointArrays; j++)
558 pointFieldNames[j] = pointData->GetArrayName(j);
561 for (int j=0; j<numberOfPointArrays; j++)
563 pointData->RemoveArray(pointFieldNames[j].c_str());
567 // Check to see if the mesh has been created
573 Info<< "Reading Mesh" << endl;
580 fvMesh::defaultRegion,
589 boolList oldSelectedRegions = selectedRegions_;
590 updateSelectedRegions();
593 meshPtr_->readUpdate() != fvMesh::UNCHANGED
594 || oldSelectedRegions != selectedRegions_
603 Info<< "converting fields" << endl;
606 const fvMesh& mesh = *meshPtr_;
608 // Construct interpolation on the raw mesh
609 Foam::pointMesh pMesh(mesh);
611 Foam::volPointInterpolation pInterp(mesh, pMesh);
613 // Search for list of objects for this time
614 Foam::IOobjectList objects(mesh, dbPtr_().timeName());
616 convertVolFields<Foam::scalar>
618 mesh, pInterp, objects, reader_->GetVolFieldSelection()
620 convertVolFields<Foam::vector>
622 mesh, pInterp, objects, reader_->GetVolFieldSelection()
624 convertVolFields<Foam::sphericalTensor>
626 mesh, pInterp, objects, reader_->GetVolFieldSelection()
628 convertVolFields<Foam::symmTensor>
630 mesh, pInterp, objects, reader_->GetVolFieldSelection()
632 convertVolFields<Foam::tensor>
634 mesh, pInterp, objects, reader_->GetVolFieldSelection()
637 convertPointFields<Foam::scalar>
639 mesh, objects, reader_->GetPointFieldSelection()
641 convertPointFields<Foam::vector>
643 mesh, objects, reader_->GetPointFieldSelection()
645 convertPointFields<Foam::sphericalTensor>
647 mesh, objects, reader_->GetPointFieldSelection()
649 convertPointFields<Foam::symmTensor>
651 mesh, objects, reader_->GetPointFieldSelection()
653 convertPointFields<Foam::tensor>
655 mesh, objects, reader_->GetPointFieldSelection()
660 Info<< "done" << endl;
665 // ************************************************************************* //