/*========================================================================= Program: Visualization Toolkit Module: vtkOpenFOAMReader.cxx Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen All rights reserved. See Copyright.txt or http://www.kitware.com/Copyright.htm for details. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the above copyright notice for more information. =========================================================================*/ // Thanks to Terry Jordan (terry.jordan@sa.netl.doe.gov) of SAIC // at the National Energy Technology Laboratory who originally developed this class. // // -------- // Takuya Oshima of Niigata University, Japan (oshima@eng.niigata-u.ac.jp) // provided the major bulk of improvements (rewrite) that made the reader // truly functional. // // Token-based FoamFile format lexer/parser, // performance/stability/compatibility enhancements, gzipped file // support, lagrangian field support, variable timestep support, // builtin cell-to-point filter, pointField support, polyhedron // decomposition support, multiregion support, // parallelization support for // decomposed cases in conjunction with vtkPOpenFOAMReader etc. // // -------- // Philippose Rajan (sarith@rocketmail.com) // provided various adjustments // // * GUI Based selection of mesh regions and fields available in the case // * Minor bug fixes / Strict memory allocation checks // * Minor performance enhancements // // -------- // Mark Olesen (OpenCFD Ltd.) www.openfoam.com // provided various bugfixes, improvements, cleanup // // --------------------------------------------------------------------------- // // Bugs or support questions should be addressed to the discourse forum // https://discourse.paraview.org/ and/or KitWare // // --------------------------------------------------------------------------- // OpenFOAM mesh files (serial), located under constant/polyMesh/ // // https://www.openfoam.com/documentation/user-guide/mesh-description.php // // - points (type: vectorField) // * x,y,z values // // - faces (type: faceList or faceCompactList) // * a list of list of nodes. // Either stored as such, or as offsets and content // // - owner (type: labelList) // * the 'owner' cell for each face. // // - neighbour (type: labelList) // * for 'neighbour' cell for each internal face. // // - boundary (type: polyBoundaryMesh) // * list of patches with startFace/nFaces for external boundary regions // // The owner cell always has a lower number than neighbour. // The face points outwards from owner to neighbour. // // To construct the internal (volume) mesh // - require points, faces, owner/neighbour. // Construct cells from owner/neighbour + faces. // // To construct the boundary mesh // - require points, faces. // The owners from the boundary faces are size (owner_list - neighbour_list). // // To construct cell zones, cell sets // - similar requirements as internal mesh // // To construct face zones, face sets // - require points, faces, owners // // To construct point zones, point sets // - require points only // // --------------------------------------------------------------------------- // Patch/mesh selection naming // single region: // - internalMesh // - group/... // - patch/... // - lagrangian/... // // multi-region: // - /regionName/internalMesh // - /regionName/group/... // - /regionName/patch/... // - /regionName/lagrangian/... // // Prefixed with "/regionName/" to provide unambiguous names. For example, // - "lagrangian/..." (lagrangian on default region) // - "/lagrangian/..." (mesh region called 'lagrangian' - silly, but accept) // // --------------------------------------------------------------------------- // Hijack the CRC routine of zlib to omit CRC check for gzipped files // (on OSes other than Windows where the mechanism doesn't work due // to pre-bound DLL symbols) if set to 1, or not (set to 0). Affects // performance by about 3% - 4%. #define VTK_FOAMFILE_OMIT_CRCCHECK 0 // The input/output buffer sizes for zlib in bytes. #define VTK_FOAMFILE_INBUFSIZE (16384) #define VTK_FOAMFILE_OUTBUFSIZE (131072) #define VTK_FOAMFILE_INCLUDE_STACK_SIZE (10) #if defined(_MSC_VER) #define _CRT_SECURE_NO_WARNINGS 1 // No strtoll on msvc: #define strtoll _strtoi64 #endif #if VTK_FOAMFILE_OMIT_CRCCHECK #define ZLIB_INTERNAL #endif // For possible future extension of linehead-aware directives #define VTK_FOAMFILE_RECOGNIZE_LINEHEAD 0 // List time directories according to system/controlDict #define VTK_FOAMFILE_LIST_TIMEDIRS_BY_CONTROLDICT 1 // Ignore things like 'U_0' restart files. // This could also be made part of the GUI properties #define VTK_FOAMFILE_IGNORE_FIELD_RESTART 1 // Support for finiteArea #define VTK_FOAMFILE_FINITE_AREA 0 // Support extra decomposition of polyhedral cells #define VTK_FOAMFILE_DECOMPOSE_POLYHEDRA 1 //------------------------------------------------------------------------------ // Developer option to debug the reader states #define VTK_FOAMFILE_DEBUG 0 // Similar to vtkErrorMacro etc. #if VTK_FOAMFILE_DEBUG #define vtkFoamDebug(x) \ do \ { \ std::cerr << "" x; \ } while (false) #else #define vtkFoamDebug(x) \ do \ { \ } while (false) #endif // VTK_FOAMFILE_DEBUG //------------------------------------------------------------------------------ #include "vtkOpenFOAMReader.h" #include "vtk_zlib.h" #include "vtksys/RegularExpression.hxx" #include "vtksys/SystemTools.hxx" #include "vtkAssume.h" #include "vtkCellArray.h" #include "vtkCellData.h" #include "vtkCharArray.h" #include "vtkCollection.h" #include "vtkDataArraySelection.h" #include "vtkDirectory.h" #include "vtkDoubleArray.h" #include "vtkFloatArray.h" #include "vtkHexahedron.h" #include "vtkInformation.h" #include "vtkInformationVector.h" #include "vtkIntArray.h" #include "vtkMultiBlockDataSet.h" #include "vtkObjectFactory.h" #include "vtkPointData.h" #include "vtkPoints.h" #include "vtkPolyData.h" #include "vtkPolygon.h" #include "vtkPolyhedron.h" #include "vtkPyramid.h" #include "vtkQuad.h" #include "vtkSmartPointer.h" #include "vtkSortDataArray.h" #include "vtkStreamingDemandDrivenPipeline.h" #include "vtkStringArray.h" #include "vtkTetra.h" #include "vtkTriangle.h" #include "vtkTypeInt32Array.h" #include "vtkTypeInt64Array.h" #include "vtkTypeInt8Array.h" #include "vtkTypeTraits.h" #include "vtkTypeUInt8Array.h" #include "vtkUnstructuredGrid.h" #include "vtkVertex.h" #include "vtkWedge.h" #if !(defined(_WIN32) && !defined(__CYGWIN__) || defined(__LIBCATAMOUNT__)) #include // For getpwnam(), getpwuid() #include #include // For getuid() #endif #include #include // For isalnum(), isdigit(), isspace() #include // For abs() #include #include #include #include #include #include #include #include #include #if VTK_FOAMFILE_OMIT_CRCCHECK uLong ZEXPORT crc32(uLong, const Bytef*, uInt) { return 0; } #endif // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // vtkStandardNewMacro(vtkOpenFOAMReader); #if VTK_FOAMFILE_FINITE_AREA // The name for finiteArea mesh static constexpr const char* const NAME_AREAMESH = "areaMesh"; #endif // The name for finiteVolume internal mesh (unzoned) static constexpr const char* const NAME_INTERNALMESH = "internalMesh"; // Index is "constant" time static constexpr int TIMEINDEX_CONSTANT = -1; // Index has not been visited static constexpr int TIMEINDEX_UNVISITED = -2; //------------------------------------------------------------------------------ // Local Functions namespace { // True if data array uses 64-bit representation for its storage bool Is64BitArray(const vtkDataArray* array) { return (array && array->GetElementComponentSize() == 8); } // Given a data array and a flag indicating whether 64 bit labels are used, // lookup and return a single element in the array. The data array must // be either a vtkTypeInt32Array or vtkTypeInt64Array. vtkTypeInt64 GetLabelValue(const vtkDataArray* array, vtkIdType idx, bool use64BitLabels) { if (!use64BitLabels) { vtkTypeInt64 result = static_cast(static_cast(array)->GetValue(idx)); assert(result >= -1); // some arrays store -1 == 'uninitialized'. return result; } else { vtkTypeInt64 result = static_cast(array)->GetValue(idx); assert(result >= -1); // some arrays store -1 == 'uninitialized'. return result; } } // Setter analogous to the above getter. void SetLabelValue(vtkDataArray* array, vtkIdType idx, vtkTypeInt64 value, bool use64BitLabels) { if (!use64BitLabels) { assert(static_cast(value) >= 0); static_cast(array)->SetValue(idx, static_cast(value)); } else { assert(value >= 0); static_cast(array)->SetValue(idx, value); } } // Another helper for appending an id to a list void AppendLabelValue(vtkDataArray* array, vtkTypeInt64 val, bool use64BitLabels) { if (!use64BitLabels) { assert(static_cast(val) >= 0); static_cast(array)->InsertNextValue(static_cast(val)); } else { assert(val >= 0); static_cast(array)->InsertNextValue(val); } } // Append unique string to list void appendUniq(vtkStringArray* list, vtkStringArray* items) { for (int i = 0; i < items->GetNumberOfTuples(); ++i) { vtkStdString& str = items->GetValue(i); if (list->LookupValue(str) == -1) { list->InsertNextValue(str); } } } // Tuple remapping for symmTensor ordering // OpenFOAM [XX XY XZ YY YZ ZZ] // VTK uses [XX YY ZZ XY YZ XZ] template void remapFoamSymmTensor(T data[]) { std::swap(data[1], data[3]); // swap XY <-> YY std::swap(data[2], data[5]); // swap XZ <-> ZZ } // Generic tuple remapping is a no-op template void remapFoamTuple(T[]) { } // Remapping for symmTensor (float) template <> void remapFoamTuple(float data[]) { ::remapFoamSymmTensor(data); } // Remapping for symmTensor (double) template <> void remapFoamTuple(double data[]) { ::remapFoamSymmTensor(data); } } // End anonymous namespace //------------------------------------------------------------------------------ // Forward Declarations struct vtkFoamDict; struct vtkFoamEntry; struct vtkFoamEntryValue; struct vtkFoamFile; struct vtkFoamIOobject; struct vtkFoamToken; //------------------------------------------------------------------------------ // class vtkFoamError // for exception-carrying object or general place to collect errors struct vtkFoamError : public std::string { vtkFoamError& operator<<(const std::string& str) { this->std::string::operator+=(str); return *this; } vtkFoamError& operator<<(const char* str) { this->std::string::operator+=(str); return *this; } template vtkFoamError& operator<<(const T& val) { std::ostringstream os; os << val; this->std::string::operator+=(os.str()); return *this; } }; //------------------------------------------------------------------------------ // Some storage containers // Manage a list of pointers template struct vtkFoamPtrList : public std::vector { private: typedef std::vector Superclass; // Plain 'delete' each entry void DeleteAll() { for (T* ptr : *this) { delete ptr; } } public: // Inherit all constructors using std::vector::vector; // Default construct vtkFoamPtrList() = default; // No copy construct/assignment vtkFoamPtrList(const vtkFoamPtrList&) = delete; void operator=(const vtkFoamPtrList&) = delete; // Destructor - delete each entry ~vtkFoamPtrList() { DeleteAll(); } // Remove top element, deleting its pointer void remove_back() { if (!Superclass::empty()) { delete Superclass::back(); Superclass::pop_back(); } } // Clear list, delete all elements void clear() { DeleteAll(); Superclass::clear(); } }; // Manage a list of vtkDataObject pointers template struct vtkFoamDataArrayVector : public std::vector { private: typedef std::vector Superclass; // Invoke vtkDataObject Delete() on each (non-null) entry void DeleteAll() { for (ObjectT* ptr : *this) { if (ptr) { ptr->Delete(); } } } public: // Destructor - invoke vtkDataObject Delete() on each entry ~vtkFoamDataArrayVector() { DeleteAll(); } // Remove top element, invoking vtkDataObject Delete() on it void remove_back() { if (!Superclass::empty()) { ObjectT* ptr = Superclass::back(); if (ptr) { ptr->Delete(); } Superclass::pop_back(); } } // Clear list, invoking vtkDataObject Delete() on each element void clear() { DeleteAll(); Superclass::clear(); } }; // Forward Declarations typedef vtkFoamDataArrayVector vtkFoamLabelArrayVector; //------------------------------------------------------------------------------ // A std::vector-like data structure where the data // lies on the stack. If the requested size in the // resize method is larger than N, the class allocates // the array on the heap. // // Unlike std::vector, the array is not default initialized // and behaves more like std::array in that manner. // // Since this simple structure is largely used for scratch space, // it allocates on growth, but not on shrinking. // It has both copying and non-copying reserve/resize methods. template struct vtkFoamStackVector { typedef T value_type; /** * Default construct, zero length and default capacity */ vtkFoamStackVector() = default; /** * Construct with specified length */ explicit vtkFoamStackVector(std::size_t len) { this->fast_resize(len); } ~vtkFoamStackVector() { if (ptr != stck) { delete[] ptr; } } bool empty() const noexcept { return !size_; } std::size_t size() const noexcept { return size_; } std::size_t capacity() const noexcept { return capacity_; } T* data() noexcept { return ptr; } const T* data() const noexcept { return ptr; } T* begin() noexcept { return ptr; } T* end() noexcept { return (ptr + size_); } const T* begin() const noexcept { return ptr; } const T* end() const noexcept { return (ptr + size_); } T& operator[](std::size_t pos) { return ptr[pos]; } const T& operator[](std::size_t pos) const { return ptr[pos]; } // Reserve space, retaining old values on growth. Uses doubling strategy. void copy_reserve(std::size_t len) { _reserve(len, false); } // Resize, retaining old values on growth. Uses doubling strategy. void copy_resize(std::size_t len) { _reserve(len, false); size_ = len; } // Faster reserve space, may discard old values on growth. Uses doubling strategy. void fast_reserve(std::size_t len) { _reserve(len, true); } // Faster resize, may discard old values on growth. Uses doubling strategy. void fast_resize(std::size_t len) { _reserve(len, true); size_ = len; } private: T stck[N]; T* ptr = stck; std::size_t capacity_ = N; std::size_t size_ = 0; // Reserve space, using doubling strategy. // Fast (non-copying) or copy/move old values on growth. void _reserve(std::size_t len, bool fast) { if (capacity_ < len) { while (capacity_ < len) { capacity_ *= 2; } if (fast) { if (ptr != stck) { delete[] ptr; } ptr = new T[capacity_]; } else { T* old = ptr; ptr = new T[capacity_]; for (size_t i = 0; i < size_; ++i) { ptr[i] = std::move(old[i]); } if (old != stck) { delete[] old; } } } } }; //------------------------------------------------------------------------------ // struct vtkFoamLabelListList - details in the implementation class struct vtkFoamLabelListList { using CellType = vtkFoamStackVector; virtual ~vtkFoamLabelListList() = default; virtual size_t GetLabelSize() const = 0; // in bytes bool IsLabel64() const { return this->GetLabelSize() == 8; } virtual vtkIdType GetNumberOfElements() const = 0; virtual vtkDataArray* GetOffsetsArray() = 0; virtual vtkDataArray* GetDataArray() = 0; virtual void ResizeExact(vtkIdType numElem, vtkIdType numValues) = 0; virtual void ResizeData(vtkIdType numValues) = 0; // Fill offsets with zero virtual void ResetOffsets() = 0; virtual vtkTypeInt64 GetBeginOffset(vtkIdType i) const = 0; virtual vtkTypeInt64 GetEndOffset(vtkIdType i) const = 0; virtual vtkIdType GetSize(vtkIdType i) const = 0; virtual void SetOffset(vtkIdType i, vtkIdType val) = 0; virtual void IncrementOffset(vtkIdType i) = 0; // Combine assignment of the new offset and accessing the data virtual void* WritePointer(vtkIdType cellId, vtkIdType dataOffset, vtkIdType elemLength) = 0; virtual vtkTypeInt64 GetValue(vtkIdType bodyIndex) const = 0; virtual void SetValue(vtkIdType bodyIndex, vtkTypeInt64 val) = 0; virtual vtkTypeInt64 GetValue(vtkIdType cellId, vtkIdType subIndex) const = 0; virtual void SetValue(vtkIdType cellId, vtkIdType subIndex, vtkTypeInt64 val) = 0; virtual void InsertValue(vtkIdType bodyIndex, vtkTypeInt64 val) = 0; virtual void GetCell(vtkIdType i, CellType& cell) const = 0; }; //------------------------------------------------------------------------------ // struct vtkFoamLabelListListImpl (implementation for vtkFoamLabelListList) // This is roughly comparable to an OpenFOAM CompactListList and largely // mirrors what the new vtkCellArray (2020: VTK_CELL_ARRAY_V2) now does. // It contains packed data and a table of offsets // template struct vtkFoamLabelListListImpl : public vtkFoamLabelListList { private: ArrayT* Offsets; ArrayT* Data; public: using LabelArrayType = ArrayT; using LabelType = typename ArrayT::ValueType; // Default construct vtkFoamLabelListListImpl() : Offsets(LabelArrayType::New()) , Data(LabelArrayType::New()) { } // Construct a shallow copy from base class explicit vtkFoamLabelListListImpl(const vtkFoamLabelListList& rhs) : Offsets(nullptr) , Data(nullptr) { assert("Require same element representation." && this->IsLabel64() == rhs.IsLabel64()); const auto& rhsCast = static_cast&>(rhs); this->Offsets = rhsCast.Offsets; this->Data = rhsCast.Data; this->Offsets->Register(nullptr); // ref count the copy this->Data->Register(nullptr); } vtkFoamLabelListListImpl(const vtkFoamLabelListListImpl& rhs) : Offsets(rhs.Offsets) , Data(rhs.Data) { this->Offsets->Register(nullptr); // ref count the copy this->Data->Register(nullptr); } void operator=(const vtkFoamLabelListListImpl&) = delete; // Destructor ~vtkFoamLabelListListImpl() override { this->Offsets->Delete(); this->Data->Delete(); } size_t GetLabelSize() const override { return sizeof(LabelType); } vtkIdType GetNumberOfElements() const override { return this->Offsets->GetNumberOfTuples() - 1; } vtkDataArray* GetOffsetsArray() override { return this->Offsets; } vtkDataArray* GetDataArray() override { return this->Data; } void ResizeExact(vtkIdType numElem, vtkIdType numValues) override { this->Offsets->SetNumberOfValues(numElem + 1); this->Data->SetNumberOfValues(numValues); this->Offsets->SetValue(0, 0); } void ResizeData(vtkIdType numValues) override { this->Data->Resize(numValues); } void ResetOffsets() override { this->Offsets->FillValue(0); } vtkTypeInt64 GetBeginOffset(vtkIdType i) const override { return this->Offsets->GetValue(i); } vtkTypeInt64 GetEndOffset(vtkIdType i) const override { return this->Offsets->GetValue(i + 1); } vtkIdType GetSize(vtkIdType i) const override { return this->Offsets->GetValue(i + 1) - this->Offsets->GetValue(i); } void SetOffset(vtkIdType i, vtkIdType val) override { this->Offsets->SetValue(i, static_cast(val)); } void IncrementOffset(vtkIdType i) override { this->Offsets->SetValue(i, this->Offsets->GetValue(i) + 1); } void* WritePointer(vtkIdType cellId, vtkIdType dataOffset, vtkIdType subLength) override { return this->Data->WritePointer(*(this->Offsets->GetPointer(cellId)) = dataOffset, subLength); } vtkTypeInt64 GetValue(vtkIdType bodyIndex) const override { return this->Data->GetValue(bodyIndex); } void SetValue(vtkIdType bodyIndex, vtkTypeInt64 value) override { this->Data->SetValue(bodyIndex, static_cast(value)); } vtkTypeInt64 GetValue(vtkIdType cellId, vtkIdType subIndex) const override { return this->Data->GetValue(this->Offsets->GetValue(cellId) + subIndex); } void SetValue(vtkIdType cellId, vtkIdType subIndex, vtkTypeInt64 value) override { this->Data->SetValue(this->Offsets->GetValue(cellId) + subIndex, static_cast(value)); } void InsertValue(vtkIdType bodyIndex, vtkTypeInt64 value) override { this->Data->InsertValue(bodyIndex, value); } void GetCell(vtkIdType i, CellType& cell) const override { auto idx = this->Offsets->GetValue(i); const auto last = this->Offsets->GetValue(i + 1); cell.fast_resize(last - idx); auto outIter = cell.begin(); while (idx != last) { *outIter = this->Data->GetValue(idx); ++outIter; ++idx; } } }; // Forward Declarations typedef vtkFoamLabelListListImpl vtkFoamLabelListList32; typedef vtkFoamLabelListListImpl vtkFoamLabelListList64; //------------------------------------------------------------------------------ // struct vtkFoamPatch // A simple struct to hold OpenFOAM boundary patch information extracted // from polyMesh/boundary. Similar to Foam::polyPatch struct vtkFoamPatch { // General patch types (fits as vtkTypeInt8) enum patchType { GEOMETRICAL = 0, // symmetryPlane, wedge, cyclic, empty, etc. PHYSICAL = 1, // patch, wall PROCESSOR = 2 // processor }; std::string name_; vtkIdType index_ = 0; vtkIdType start_ = 0; vtkIdType size_ = 0; vtkIdType offset_ = 0; // The start-face offset into all boundaries patchType type_ = patchType::GEOMETRICAL; bool owner_ = true; // Patch owner (processor patch) // The first patch face vtkIdType startFace() const noexcept { return (this->start_); } // One beyond the last patch face vtkIdType endFace() const noexcept { return (this->start_ + this->size_); } // The patch local face (as per OpenFOAM polyPatch) vtkIdType whichFace(vtkIdType meshFacei) const { return (meshFacei - this->start_); } }; //------------------------------------------------------------------------------ // struct vtkFoamBoundaries // A collection of boundary patches with additional grouping and selection information struct vtkFoamBoundaries : public std::vector { // Collect and forwarding of errors (cannot use vtkErrorMacro here) vtkFoamError error_; // Patch groups, according to the inGroups keyword std::map> groups; // Active patch groups std::unordered_set groupActive; // Active patch indices, selected directly std::unordered_set patchActive; // Active patch indices, selected by group std::unordered_set patchActiveByGroup; // Reset group and patch selections void clearSelections() { groupActive.clear(); patchActive.clear(); patchActiveByGroup.clear(); } // Reset storage and errors, leaves timeName intact void clearAll() { this->clear(); error_.clear(); groups.clear(); this->clearSelections(); } const vtkFoamError& error() const noexcept { return error_; } vtkFoamError& error() noexcept { return error_; } // The start label of boundary faces in the polyMesh face list. // Same as mesh nInternalFaces() if boundaries exist vtkIdType startFace() const { return this->empty() ? 0 : this->front().startFace(); } // One beyond the last boundary face vtkIdType endFace() const { return this->empty() ? 0 : this->back().endFace(); } void enablePatch(vtkIdType patchIndex) { patchActive.emplace(patchIndex); } void enableGroup(const std::string& groupName) { auto citer = groups.find(groupName); if (citer != groups.end()) { const std::vector& patchIndices = citer->second; for (const vtkIdType patchIndex : patchIndices) { patchActiveByGroup.emplace(patchIndex); } } } // True if given patch index is active bool isActive(vtkIdType patchIndex) const { return (patchActive.find(patchIndex) != patchActive.end()) || (patchActiveByGroup.find(patchIndex) != patchActiveByGroup.end()); } // Set contents from dictionary // Return false on errors bool update(const vtkFoamDict& dict); // The patch index for a given face label, -1 for internal face or out-of-bounds vtkIdType whichPatch(vtkIdType faceIndex) const; }; //------------------------------------------------------------------------------ // struct vtkFoamZones // A collection of names id-lists, used for OpenFOAM zones or sets. // Stored as an unordered map instead of being ordered or a vector of items, // since any ordering (like zones) will seen during input and managed with the // VTK block structure. // // The idea is to maintain a list of ids (cell,face,point) in a cache that is // separate from the mesh to allow flexible usage later. // Also, it does not make any sense to have an entry like "CellId" in the CellData // since that not only complicates handling, but is also quite misleading when local // cell ids have been assembled from different processors. struct vtkFoamZones { // Representation for the zone or set type enum zoneType { UNKNOWN = 0, // placeholder POINT = 1, // pointZone FACE = 2, // faceZone CELL = 3 // cellZone }; // Collect and forwarding of errors (cannot use vtkErrorMacro here) vtkFoamError error_; // The {cell,face,point}Labels per zone std::unordered_map> zones_; // The zone type zoneType type_ = zoneType::UNKNOWN; // If zone map ids have content bool empty() const { return zones_.empty(); } // Reset storage and errors void clearAll() { error_.clear(); zones_.clear(); } // Clear and reset the zone type void reset(enum zoneType ztype) { this->clearAll(); type_ = ztype; } const vtkFoamError& error() const noexcept { return error_; } vtkFoamError& error() noexcept { return error_; } // Find zone by name and return list of ids or nullptr on failure vtkIdList* findZone(const std::string& zoneName) { auto iter = zones_.find(zoneName); if (iter != zones_.end()) { return iter->second; } return nullptr; } }; //------------------------------------------------------------------------------ // Simple handling of common OpenFOAM data types struct vtkFoamTypes { // Primitive types, with nComponents encoded in lower 4 bits enum dataType { NO_TYPE = 0, SCALAR_TYPE = 1, VECTOR_TYPE = 3, SYMM_TENSOR_TYPE = 6, TENSOR_TYPE = 9, // Single-component types, but disambiguate from SCALAR_TYPE BOOL_TYPE = (0x10 | SCALAR_TYPE), LABEL_TYPE = (0x20 | SCALAR_TYPE), SPH_TENSOR_TYPE = (0x30 | SCALAR_TYPE) }; // The number of data components static int GetNumberOfComponents(const dataType dtype) noexcept { return (dtype & 0xF); } static bool IsGood(dataType dtype) noexcept { return dtype != NO_TYPE; } static bool IsBool(dataType dtype) noexcept { return dtype == BOOL_TYPE; } static bool IsLabel(dataType dtype) noexcept { return dtype == LABEL_TYPE; } static bool IsScalar(dataType dtype) noexcept { return dtype == SCALAR_TYPE; } static bool IsNumeric(dataType dtype) noexcept { return IsLabel(dtype) || IsScalar(dtype); } // Is a VectorSpace type? static bool IsVectorSpace(dataType dtype) noexcept { return GetNumberOfComponents(dtype) > 1 || dtype == SPH_TENSOR_TYPE; } // Parse things like "scalarField" or "ScalarField" -> SCALAR_TYPE etc. // Ignore case on first letter (at pos), which makes it convenient for "volScalarField" too. static dataType FieldToEnum(const std::string& fieldTypeName, size_t pos = 0); // Handle "List" -> SCALAR_TYPE etc. static dataType ListToEnum(const std::string& listTypeName); private: // Implementation for FieldToEnum, ListToEnum static dataType ToEnumImpl(const std::string& str, size_t pos, size_t len, bool ignoreCase); }; //------------------------------------------------------------------------------ // class vtkOpenFOAMReaderPrivate // the reader core of vtkOpenFOAMReader class vtkOpenFOAMReaderPrivate : public vtkObject { public: // Use sparingly friend class vtkOpenFOAMReader; static vtkOpenFOAMReaderPrivate* New(); vtkTypeMacro(vtkOpenFOAMReaderPrivate, vtkObject); vtkGetMacro(TimeStep, int); vtkSetMacro(TimeStep, int); double GetTimeValue() const; void SetTimeValue(double requestedTime); vtkStringArray* GetTimeNames() { return this->TimeNames; } vtkDoubleArray* GetTimeValues() { return this->TimeValues; } // Print some time information (names, current time-step) void PrintTimes(std::ostream& os, vtkIndent indent, bool full = false) const; bool HasPolyMesh() const noexcept { return !this->PolyMeshTimeIndexFaces.empty(); } const std::string& GetRegionName() const noexcept { return this->RegionName; } vtkStringArray* GetLagrangianPaths() { return this->LagrangianPaths; } // Read mesh/fields and create dataset int RequestData(vtkMultiBlockDataSet* output); int MakeMetaDataAtTimeStep(vtkStringArray*, vtkStringArray*, vtkStringArray*, bool); // Gather time instances information and create cache for mesh times bool MakeInformationVector(const std::string& casePath, const std::string& controlDictPath, const std::string& procName, vtkOpenFOAMReader* parent, bool requirePolyMesh = true); // Use given time instances information and create cache for mesh times bool MakeInformationVector(const std::string& casePath, const std::string& procName, vtkOpenFOAMReader* parent, vtkStringArray* timeNames, vtkDoubleArray* timeValues, bool requirePolyMesh = true); // Copy time instances information and create cache for mesh times void SetupInformation(const std::string& casePath, const std::string& regionName, const std::string& procName, vtkOpenFOAMReaderPrivate* master, bool requirePolyMesh = true); private: vtkOpenFOAMReader* Parent; std::string CasePath; // The full path to the case - includes trailing '/' std::string RegionName; // Region name. Empty for default region std::string ProcessorName; // Processor subdirectory. Empty for serial case // Time information vtkDoubleArray* TimeValues; // Time values vtkStringArray* TimeNames; // Directory names // Topology indices into TimeValues, TimeName std::vector PolyMeshTimeIndexPoints; std::vector PolyMeshTimeIndexFaces; // Indices into TimeValues, TimeName int TimeStep; int TimeStepOld; // Topology time index, driven by PolyMeshTimeIndexFaces int TopologyTimeIndex; int InternalMeshSelectionStatus; int InternalMeshSelectionStatusOld; // filenames / directories vtkStringArray* VolFieldFiles; vtkStringArray* DimFieldFiles; vtkStringArray* AreaFieldFiles; vtkStringArray* PointFieldFiles; vtkStringArray* LagrangianFieldFiles; // The cloud paths (region-local) vtkNew LagrangianPaths; // Mesh dimensions and construction information vtkIdType NumPoints; vtkIdType NumInternalFaces; vtkIdType NumFaces; vtkIdType NumCells; // The face owner, neighbour (labelList) vtkDataArray* FaceOwner; vtkDataArray* FaceNeigh; // For cell-to-point interpolation vtkPolyData* AllBoundaries; vtkDataArray* AllBoundariesPointMap; vtkDataArray* InternalPoints; // For caching mesh vtkUnstructuredGrid* InternalMesh; vtkMultiBlockDataSet* BoundaryMesh; vtkFoamLabelArrayVector* BoundaryPointMap; vtkFoamBoundaries BoundaryDict; // Zones vtkFoamZones cellZoneMap; vtkFoamZones faceZoneMap; vtkFoamZones pointZoneMap; vtkMultiBlockDataSet* CellZoneMesh; vtkMultiBlockDataSet* FaceZoneMesh; vtkMultiBlockDataSet* PointZoneMesh; #if VTK_FOAMFILE_DECOMPOSE_POLYHEDRA // For polyhedral decomposition vtkIdType NumTotalAdditionalCells; vtkIdTypeArray* AdditionalCellIds; vtkIntArray* NumAdditionalCells; vtkFoamLabelArrayVector* AdditionalCellPoints; #endif #if VTK_FOAMFILE_FINITE_AREA vtkFoamZones areaMeshMap; vtkPolyData* AreaMesh; #endif // Constructor and destructor are kept private vtkOpenFOAMReaderPrivate(); ~vtkOpenFOAMReaderPrivate() override; vtkOpenFOAMReaderPrivate(const vtkOpenFOAMReaderPrivate&) = delete; void operator=(const vtkOpenFOAMReaderPrivate&) = delete; // Clear mesh construction void ClearInternalMeshes(); void ClearBoundaryMeshes(); void ClearZoneMeshes(); void ClearAreaMeshes(); void ClearMeshes(); // The subdirectory for a region. Eg, "/solid". Empty for default region. std::string RegionPath() const { if (this->RegionName.empty()) { return ""; } return ("/" + this->RegionName); } // Prefix display qualifier for a region. Eg, "/solid/". Empty for default region. std::string RegionPrefix() const { if (this->RegionName.empty()) { return ""; } return ("/" + this->RegionName + "/"); } // Test if display (selection) name matches the current region. // See RegionPrefix() comments bool IsDisplayRegion(const std::string& displayName) const { if (this->RegionName.empty()) { return (displayName[0] != '/'); } else if (displayName[0] != '/') { return false; } // Match "/regionName/..." const auto slash1 = displayName.find('/', 1); return (slash1 != std::string::npos) && (displayName.compare(1, slash1 - 1, this->RegionName) == 0); } // The timeName for the given index, with special handling for "constant" time directory std::string TimePath(int timeIndex) const { if (timeIndex < 0) { return this->CasePath + "constant"; } return this->CasePath + this->TimeNames->GetValue(timeIndex); } std::string CurrentTimePath() const { return this->TimePath(this->TimeStep); } // TimePath + region std::string TimeRegionPath(int timeIndex) const { return this->TimePath(timeIndex) + this->RegionPath(); } std::string CurrentTimeRegionPath() const { return this->TimeRegionPath(this->TimeStep); } std::string CurrentTimeRegionPath(const std::vector& indexer) const { return this->TimeRegionPath(indexer[this->TimeStep]); } #if VTK_FOAMFILE_DEBUG void PrintMeshTimes(const char* name, const std::vector&) const; // For debugging #endif // Search time directories for mesh void PopulateMeshTimeIndices(); void AddFieldName( const std::string& fieldName, const std::string& fieldType, bool isLagrangian = false); // Search a time directory for field objects void GetFieldNames(const std::string&, bool isLagrangian = false); void SortFieldFiles(vtkStringArray* selections, vtkStringArray* files); void LocateLagrangianClouds(const std::string& timePath); #if VTK_FOAMFILE_LIST_TIMEDIRS_BY_CONTROLDICT // List time directories according to system/controlDict vtkFoamError ListTimeDirectoriesByControlDict(const std::string& controlDictPath); #endif // List time directories by searching in a case directory bool ListTimeDirectoriesByInstances(); // Read polyMesh/points (vectorField) vtkSmartPointer ReadPointsFile(const std::string& timeRegionDir); // Read polyMesh/faces (faceCompactList or faceList) std::unique_ptr ReadFacesFile(const std::string& timeRegionDir); // Read polyMesh/{owner,neighbour}, check overall number of faces. bool ReadOwnerNeighbourFiles(const std::string& timeRegionDir); // Create meshCells from owner/neighbour information std::unique_ptr CreateCellFaces(); bool CheckFaceList(const vtkFoamLabelListList& faces); // Create volume mesh void InsertCellsToGrid(vtkUnstructuredGrid*, std::unique_ptr& meshCellsPtr, const vtkFoamLabelListList& meshFaces, vtkIdList* cellLabels = nullptr #if VTK_FOAMFILE_DECOMPOSE_POLYHEDRA , vtkIdTypeArray* additionalCellIds = nullptr, vtkFloatArray* pointArray = nullptr #endif ); vtkUnstructuredGrid* MakeInternalMesh(std::unique_ptr& meshCellsPtr, const vtkFoamLabelListList& meshFaces, vtkFloatArray* pointArray); void InsertFacesToGrid(vtkPolyData*, const vtkFoamLabelListList& meshFaces, vtkIdType startFace, vtkIdType endFace, vtkIdList* faceLabels = nullptr, vtkDataArray* pointMap = nullptr, bool isLookupValue = false); vtkMultiBlockDataSet* MakeBoundaryMesh( const vtkFoamLabelListList& meshFaces, vtkFloatArray* pointArray); // Move additional points for decomposed cells bool MoveInternalMesh(vtkUnstructuredGrid*, vtkFloatArray*); bool MoveBoundaryMesh(vtkMultiBlockDataSet*, vtkFloatArray*); // cell-to-point interpolator void InterpolateCellToPoint( vtkFloatArray*, vtkFloatArray*, vtkPointSet*, vtkDataArray*, vtkTypeInt64); // Convert OpenFOAM dimension array to string std::string ConstructDimensions(const vtkFoamDict& dict) const; // read and create cell/point fields bool ReadFieldFile(vtkFoamIOobject& io, vtkFoamDict& dict, const std::string& varName, const vtkDataArraySelection* selection); vtkSmartPointer FillField(vtkFoamEntry& entry, vtkIdType nElements, const vtkFoamIOobject& io, vtkFoamTypes::dataType fieldDataType); void GetVolFieldAtTimeStep(const std::string& varName, bool isInternalField = false); void GetPointFieldAtTimeStep(const std::string& varName); #if VTK_FOAMFILE_FINITE_AREA void GetAreaFieldAtTimeStep(const std::string& varName); #endif // Create lagrangian mesh/fields vtkMultiBlockDataSet* MakeLagrangianMesh(); // Read specified file (typeName) from polyMesh directory, using faces instance std::unique_ptr GetPolyMeshFile(const std::string& typeName, bool mandatory); // Create (cell|face|point) zones bool GetCellZoneMesh(vtkMultiBlockDataSet* zoneMesh, std::unique_ptr& meshCellsPtr, const vtkFoamLabelListList& meshFaces, vtkPoints*); bool GetFaceZoneMesh( vtkMultiBlockDataSet* zoneMesh, const vtkFoamLabelListList& meshFaces, vtkPoints*); bool GetPointZoneMesh(vtkMultiBlockDataSet* zoneMesh, vtkPoints*); #if VTK_FOAMFILE_FINITE_AREA // Mechanism for finiteArea mesh is similar to faceZone bool GetAreaMesh(vtkPolyData* areaMesh, const vtkFoamLabelListList& meshFaces, vtkPoints*); #endif }; vtkStandardNewMacro(vtkOpenFOAMReaderPrivate); //------------------------------------------------------------------------------ // Local Functions namespace { // Set named block void SetBlock(vtkMultiBlockDataSet* parent, unsigned int blockIndex, vtkDataObject* block, const std::string& name) { parent->SetBlock(blockIndex, block); parent->GetMetaData(blockIndex)->Set(vtkCompositeDataSet::NAME(), name.c_str()); } // Append named block void AppendBlock(vtkMultiBlockDataSet* parent, vtkDataObject* block, const std::string& name) { ::SetBlock(parent, parent->GetNumberOfBlocks(), block, name); } // Set array name and fieldData attributes // The optional suffix is for dimensions etc void AddArrayToFieldData(vtkDataSetAttributes* fieldData, vtkDataArray* array, const std::string& name, const std::string& suffix = "") { if (suffix.empty()) { array->SetName(name.c_str()); } else { array->SetName((name + suffix).c_str()); } if (array->GetNumberOfComponents() == 1 && name == "p") { fieldData->SetScalars(array); } else if (array->GetNumberOfComponents() == 3 && name == "U") { fieldData->SetVectors(array); } else { fieldData->AddArray(array); } } } // End anonymous namespace //------------------------------------------------------------------------------ // Simple handling of common OpenFOAM data types // Low-level implementation vtkFoamTypes::dataType vtkFoamTypes::ToEnumImpl( const std::string& str, size_t pos, size_t last, bool ignoreCase) { vtkFoamTypes::dataType dtype(vtkFoamTypes::NO_TYPE); char firstChar = str[pos]; if (ignoreCase) { firstChar = std::tolower(firstChar); } ++pos; // First character handled separately (for ignoring case) size_t len = std::string::npos; if (last != std::string::npos) { if (last > pos) { len = last - pos; } else { // Caught bad input firstChar = '\0'; } } switch (firstChar) { case '\0': { break; } case 'b': { if (str.compare(pos, len, "ool") == 0) { // (Bool | bool) dtype = vtkFoamTypes::BOOL_TYPE; } break; } case 'l': { if (str.compare(pos, len, "abel") == 0) { // (Label | label) dtype = vtkFoamTypes::LABEL_TYPE; } break; } case 's': { if (str.compare(pos, len, "calar") == 0) { // (Scalar | scalar) dtype = vtkFoamTypes::SCALAR_TYPE; } else if (str.compare(pos, len, "phericalTensor") == 0) { // (SphericalTensor | sphericalTensor) dtype = vtkFoamTypes::SPH_TENSOR_TYPE; } else if (str.compare(pos, len, "ymmTensor") == 0) { // (SymmTensor | symmTensor) dtype = vtkFoamTypes::SYMM_TENSOR_TYPE; } break; } case 't': { if (str.compare(pos, len, "ensor") == 0) { // (Tensor | tensor) dtype = vtkFoamTypes::TENSOR_TYPE; } break; } case 'v': { if (str.compare(pos, len, "ector") == 0) { // (Vector | vector) dtype = vtkFoamTypes::VECTOR_TYPE; } break; } } return dtype; } // Fields: expects scalarField, volScalarField etc. vtkFoamTypes::dataType vtkFoamTypes::FieldToEnum(const std::string& fieldTypeName, size_t pos) { // With ignoreCase return vtkFoamTypes::ToEnumImpl(fieldTypeName, pos, fieldTypeName.find("Field", pos), true); } // Lists: expects "List", "List" etc. vtkFoamTypes::dataType vtkFoamTypes::ListToEnum(const std::string& listTypeName) { const auto endp = listTypeName.find('>'); if ((endp != std::string::npos) && (endp + 1 == listTypeName.length()) && listTypeName.compare(0, 5, "List<") == 0) { // Without ignoreCase return vtkFoamTypes::ToEnumImpl(listTypeName, 5, endp, false); } return vtkFoamTypes::NO_TYPE; } //------------------------------------------------------------------------------ // class vtkFoamStreamOption // Some elements from Foam::IOstreamOption and from Foam::IOstream // - format (ASCII | BINARY) // - label, scalar sizes // // Note: all enums pack into 32-bits, so we can use them in vtkFoamToken, vtkFoamFile etc. // without adversely affecting the size of the structures struct vtkFoamStreamOption { public: // The OpenFOAM input stream format is ASCII or BINARY enum fileFormat : unsigned char { ASCII = 0, // ASCII unless otherwise specified BINARY }; // Bitwidth of an OpenFOAM label (integer type). // Corresponds to WM_LABEL_SIZE (32|64) enum labelType : unsigned char { INT32, INT64 }; // Bitwidth of an OpenFOAM scalar (floating-point type) // Corresponds to WM_PRECISION_OPTION (SP|DP|SPDP) enum scalarType : unsigned char { FLOAT32, FLOAT64 }; private: fileFormat Format = fileFormat::ASCII; labelType LabelType = labelType::INT32; scalarType ScalarType = scalarType::FLOAT64; public: // Default construct. ASCII, Int32, double precision vtkFoamStreamOption() = default; // Construct with specified handling for labels/floats vtkFoamStreamOption(const bool use64BitLabels, const bool use64BitFloats) { this->SetLabel64(use64BitLabels); this->SetFloat64(use64BitFloats); } bool IsAsciiFormat() const noexcept { return this->Format == fileFormat::ASCII; } bool IsLabel64() const noexcept { return this->LabelType == labelType::INT64; } bool IsFloat64() const noexcept { return this->ScalarType == scalarType::FLOAT64; } void SetBinaryFormat(const bool on) { this->Format = (on ? fileFormat::BINARY : fileFormat::ASCII); } void SetLabel64(const bool on) noexcept { this->LabelType = (on ? labelType::INT64 : labelType::INT32); } void SetFloat64(const bool on) noexcept { this->ScalarType = (on ? scalarType::FLOAT64 : scalarType::FLOAT32); } const vtkFoamStreamOption& GetStreamOption() const noexcept { return static_cast(*this); } void SetStreamOption(const vtkFoamStreamOption& opt) noexcept { static_cast(*this) = opt; } }; //------------------------------------------------------------------------------ // class vtkFoamToken // token class which also works as container for list types // - a word token is treated as a string token for simplicity // - handles only atomic types. Handling of list types are left to the // derived classes. struct vtkFoamToken : public vtkFoamStreamOption { public: enum tokenType { // Undefined type UNDEFINED = 0, // atomic types PUNCTUATION, LABEL, SCALAR, STRING, IDENTIFIER, // List types (vtkObject-derived) BOOLLIST, LABELLIST, SCALARLIST, VECTORLIST, STRINGLIST, // List types (non-vtkObject) LABELLISTLIST, ENTRYVALUELIST, EMPTYLIST, DICTIONARY, // error state TOKEN_ERROR }; protected: tokenType Type = tokenType::UNDEFINED; union { char Char; vtkTypeInt64 Int; double Double; // Any/all pointer types void* AnyPointer; std::string* StringPtr; // List types (vtkObject-derived) vtkObjectBase* VtkObjectPtr; vtkTypeInt8Array* BoolListPtr; vtkDataArray* LabelListPtr; vtkFloatArray* ScalarListPtr; vtkFloatArray* VectorListPtr; vtkStringArray* StringListPtr; // List types (non-vtkObject) vtkFoamLabelListList* LabelListListPtr; vtkFoamPtrList* EntryValuePtrs; vtkFoamDict* DictPtr; }; void Clear() { if (this->Type == STRING || this->Type == IDENTIFIER) // IsStringType { delete this->StringPtr; } } void AssignData(const vtkFoamToken& tok) { switch (tok.Type) { case PUNCTUATION: this->Char = tok.Char; break; case LABEL: this->Int = tok.Int; break; case SCALAR: this->Double = tok.Double; break; case STRING: case IDENTIFIER: this->StringPtr = new std::string(*tok.StringPtr); break; default: break; } } public: // Default construct vtkFoamToken() = default; vtkFoamToken(const vtkFoamToken& tok) : vtkFoamStreamOption(tok) , Type(tok.Type) { this->AssignData(tok); } ~vtkFoamToken() { this->Clear(); } tokenType GetType() const { return this->Type; } template bool Is() const; template T To() const; #if defined(_MSC_VER) // workaround for Win32-64ids-nmake70 template <> bool Is() const; template <> bool Is() const; template <> bool Is() const; template <> bool Is() const; template <> vtkTypeInt32 To() const; template <> vtkTypeInt64 To() const; template <> float To() const; template <> double To() const; #endif // Token represents PUNCTUATION bool IsPunctuation() const noexcept { return this->Type == PUNCTUATION; } // Token is PUNCTUATION and equal to parameter bool IsPunctuation(const char c) const noexcept { return this->Type == PUNCTUATION && c == this->Char; } // Token represents an LABEL (integer) value bool IsLabel() const noexcept { return this->Type == LABEL; } // Token is LABEL (integer) value and equal to parameter bool IsLabel(const vtkTypeInt64 val) const noexcept { return this->Type == LABEL && val == this->Int; } // Token represents a SCALAR (floating-point) value bool IsScalar() const noexcept { return this->Type == SCALAR; } // Token represents a numerical value bool IsNumeric() const noexcept { return this->Type == LABEL || this->Type == SCALAR; } // Token is STRING bool IsString() const noexcept { return this->Type == STRING; } // Token is STRING and equal to parameter bool IsString(const std::string& str) const { return this->Type == STRING && str == *this->StringPtr; } // Token represents string content bool IsStringType() const noexcept { return this->Type == STRING || this->Type == IDENTIFIER; } // Integer value from LABEL token without checks vtkTypeInt64 ToInt() const noexcept { return this->Int; } // Mostly the same as To, with additional check float ToFloat() const noexcept { return this->Type == LABEL ? static_cast(this->Int) : this->Type == SCALAR ? static_cast(this->Double) : 0.0F; } // Mostly the same as To, with additional check double ToDouble() const noexcept { return this->Type == LABEL ? static_cast(this->Int) : this->Type == SCALAR ? this->Double : 0.0; } std::string ToString() const { return *this->StringPtr; } std::string ToIdentifier() const { return *this->StringPtr; } // Clear token and set to be ERROR. void SetBad() { this->Clear(); this->Type = TOKEN_ERROR; } void SetIdentifier(const std::string& idString) { this->operator=(idString); this->Type = IDENTIFIER; } void operator=(const char c) { this->Clear(); this->Type = PUNCTUATION; this->Char = c; } void operator=(const vtkTypeInt32 val) { this->Clear(); this->Type = LABEL; this->Int = static_cast(val); if (this->IsLabel64()) { vtkGenericWarningMacro("Assigned int32 to int64 label"); } } void operator=(const vtkTypeInt64 val) { this->Clear(); this->Type = LABEL; this->Int = val; if (!this->IsLabel64()) { vtkGenericWarningMacro("Assigned int64 to int32 label - may lose precision"); } } void operator=(const double val) { this->Clear(); this->Type = SCALAR; this->Double = val; } void operator=(const char* str) { this->Clear(); this->Type = STRING; this->StringPtr = new std::string(str); } void operator=(const std::string& str) { this->Clear(); this->Type = STRING; this->StringPtr = new std::string(str); } vtkFoamToken& operator=(const vtkFoamToken& tok) { this->Clear(); this->SetStreamOption(tok); this->Type = tok.Type; this->AssignData(tok); return *this; } bool operator==(const char c) const noexcept { return this->IsPunctuation(c); } bool operator!=(const char c) const noexcept { return !this->IsPunctuation(c); } bool operator==(const vtkTypeInt32 val) const { return this->IsLabel(val); } bool operator==(const vtkTypeInt64 val) const { return this->IsLabel(val); } bool operator==(const std::string& str) const { return this->IsString(str); } bool operator!=(const std::string& str) const { return !this->IsString(str); } friend std::ostringstream& operator<<(std::ostringstream& os, const vtkFoamToken& tok) { switch (tok.GetType()) { case TOKEN_ERROR: os << "badToken (an unexpected EOF?)"; break; case PUNCTUATION: os << tok.Char; break; case LABEL: if (tok.IsLabel64()) { os << tok.Int; } else { os << static_cast(tok.Int); } break; case SCALAR: os << tok.Double; break; case STRING: case IDENTIFIER: os << *(tok.StringPtr); break; default: break; } return os; } }; //------------------------------------------------------------------------------ // Specializations for vtkFoamToken template <> inline bool vtkFoamToken::Is() const { // masquerade for bool return this->Type == LABEL; } template <> inline bool vtkFoamToken::Is() const { return this->Type == LABEL && !(this->IsLabel64()); } template <> inline bool vtkFoamToken::Is() const { return this->Type == LABEL; } template <> inline bool vtkFoamToken::Is() const { return this->Type == LABEL || this->Type == SCALAR; } template <> inline bool vtkFoamToken::Is() const { return this->Type == SCALAR; } // ie, a bool value template <> inline vtkTypeInt8 vtkFoamToken::To() const { return static_cast(this->Int); } template <> inline vtkTypeInt32 vtkFoamToken::To() const { if (this->IsLabel64()) { vtkGenericWarningMacro("Casting int64 label to int32 - may lose precision"); } return static_cast(this->Int); } template <> inline vtkTypeInt64 vtkFoamToken::To() const { return this->Int; } template <> inline float vtkFoamToken::To() const { return this->Type == LABEL ? static_cast(this->Int) : static_cast(this->Double); } template <> inline double vtkFoamToken::To() const { return this->Type == LABEL ? static_cast(this->Int) : this->Double; } //------------------------------------------------------------------------------ // class vtkFoamFileStack // list of variables that have to be saved when a file is included. struct vtkFoamFileStack { protected: vtkOpenFOAMReader* Reader; // GUI preference std::string FileName; FILE* File; z_stream Z; int ZStatus; int LineNumber; bool IsCompressed; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD bool WasNewline; #endif // buffer pointers. using raw pointers for performance reason. unsigned char* Inbuf; unsigned char* Outbuf; unsigned char* BufPtr; unsigned char* BufEndPtr; vtkFoamFileStack(vtkOpenFOAMReader* reader) : Reader(reader) , File(nullptr) , ZStatus(Z_OK) , LineNumber(0) , IsCompressed(false) #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD , WasNewline(true) #endif , Inbuf(nullptr) , Outbuf(nullptr) , BufPtr(nullptr) , BufEndPtr(nullptr) { this->Z.zalloc = Z_NULL; this->Z.zfree = Z_NULL; this->Z.opaque = Z_NULL; } void Reset() { // this->FileName = ""; this->File = nullptr; // this->ZStatus = Z_OK; this->Z.zalloc = Z_NULL; this->Z.zfree = Z_NULL; this->Z.opaque = Z_NULL; // this->LineNumber = 0; this->IsCompressed = false; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->WasNewline = true; #endif this->Inbuf = nullptr; this->Outbuf = nullptr; // this->BufPtr = nullptr; // this->BufEndPtr = nullptr; } public: const std::string& GetFileName() const noexcept { return this->FileName; } int GetLineNumber() const noexcept { return this->LineNumber; } // Try to open the file. Return non-empty error string on failure vtkFoamError TryOpen(const std::string& fileName) { vtkFoamError errors; do { // Line number 0 to indicate beginning of file when an exception is thrown this->LineNumber = 0; this->FileName = fileName; if (this->File) { errors << "File already opened within this object"; break; } this->File = vtksys::SystemTools::Fopen(this->FileName, "rb"); if (this->File == nullptr) { errors << "Cannot open file for reading"; break; } unsigned char zMagic[2]; if (fread(zMagic, 1, 2, this->File) == 2 && zMagic[0] == 0x1f && zMagic[1] == 0x8b) { // gzip-compressed format this->Z.avail_in = 0; this->Z.next_in = Z_NULL; // + 32 to automatically recognize gzip format if (inflateInit2(&this->Z, 15 + 32) == Z_OK) { this->IsCompressed = true; this->Inbuf = new unsigned char[VTK_FOAMFILE_INBUFSIZE]; } else { fclose(this->File); this->File = nullptr; errors << "Cannot init zstream"; if (this->Z.msg) { errors << " " << this->Z.msg; } break; } } else { this->IsCompressed = false; } rewind(this->File); this->ZStatus = Z_OK; this->Outbuf = new unsigned char[VTK_FOAMFILE_OUTBUFSIZE + 1]; this->BufPtr = this->Outbuf + 1; this->BufEndPtr = this->BufPtr; this->LineNumber = 1; } while (false); return errors; } void CloseCurrentFile() { if (this->IsCompressed) { inflateEnd(&this->Z); } delete[] this->Inbuf; delete[] this->Outbuf; this->Inbuf = this->Outbuf = nullptr; if (this->File) { fclose(this->File); this->File = nullptr; } // don't reset the line number so that the last line number is // retained after close // lineNumber_ = 0; } }; //------------------------------------------------------------------------------ // class vtkFoamFile // Read and tokenize the input. Retains format and label/scalar size informatio struct vtkFoamFile : public vtkFoamStreamOption , public vtkFoamFileStack { private: typedef vtkFoamFileStack Superclass; // Find last slash (os-specific) static size_t rfind_slash(const std::string& str, size_t pos = std::string::npos) noexcept { #if defined(_WIN32) return str.find_last_of("/\\", pos); #else return str.find_last_of('/', pos); #endif } // String equivalent cwd (os-specific) static std::string cwd_string() noexcept { #if defined(_WIN32) return std::string(".\\"); #else return std::string("./"); #endif } public: // The dictionary #inputMode values enum inputMode { INPUT_MODE_MERGE, INPUT_MODE_OVERWRITE, INPUT_MODE_PROTECT, INPUT_MODE_WARN, INPUT_MODE_ERROR }; // Generic exception throwing with stack trace void ThrowStackTrace(const std::string& msg); private: std::string CasePath; // The full path to the case - includes trailing '/' // The current input mode inputMode InputMode; // Handling include files vtkFoamFileStack* Stack[VTK_FOAMFILE_INCLUDE_STACK_SIZE]; int StackI; bool InflateNext(unsigned char* buf, size_t requestSize, vtkTypeInt64* readSize = nullptr); int NextTokenHead(); // Keep exception throwing / recursive codes out-of-line to make // putBack(), getc() and readExpecting() inline expandable void ThrowDuplicatedPutBackException(); void ThrowUnexpectedEOFException(); void ThrowUnexpectedNondigitException(int c); void ThrowUnexpectedTokenException(char, int c); int ReadNext(); void PutBack(const int c) { if (--this->Superclass::BufPtr < this->Superclass::Outbuf) { this->ThrowDuplicatedPutBackException(); } *this->Superclass::BufPtr = static_cast(c); } // get a character int Getc() { return this->Superclass::BufPtr == this->Superclass::BufEndPtr ? this->ReadNext() : *this->Superclass::BufPtr++; } vtkFoamError StackString() { vtkFoamError err; if (this->StackI > 0) { err << "\n included"; for (int stackI = this->StackI - 1; stackI >= 0; stackI--) { err << " from line " << this->Stack[stackI]->GetLineNumber() << " of " << this->Stack[stackI]->GetFileName() << "\n"; } err << ": "; } return err; } bool CloseIncludedFile() { if (this->StackI == 0) { return false; } this->StackI--; this->Superclass::CloseCurrentFile(); // use the default bitwise assignment operator this->Superclass::operator=(*this->Stack[this->StackI]); delete this->Stack[this->StackI]; return true; } public: // No default construct, copy or assignment vtkFoamFile() = delete; vtkFoamFile(const vtkFoamFile&) = delete; void operator=(const vtkFoamFile&) = delete; vtkFoamFile(const std::string& casePath, vtkOpenFOAMReader* reader) : vtkFoamStreamOption(reader->GetUse64BitLabels(), reader->GetUse64BitFloats()) , vtkFoamFileStack(reader) , CasePath(casePath) , InputMode(INPUT_MODE_MERGE) , StackI(0) { } ~vtkFoamFile() { this->Close(); } std::string GetCasePath() const noexcept { return this->CasePath; } std::string GetFilePath() const { return vtkFoamFile::ExtractPath(this->FileName); } inputMode GetInputMode() const noexcept { return this->InputMode; } void Open(const std::string& fileName) { vtkFoamError err = this->Superclass::TryOpen(fileName); if (!err.empty()) { this->ThrowStackTrace(err); } } void Close() { while (this->CloseIncludedFile()) ; this->CloseCurrentFile(); // Reinstate values from reader (eg, GUI) auto& streamOpt = static_cast(*this); streamOpt.SetLabel64(this->Reader->GetUse64BitLabels()); streamOpt.SetFloat64(this->Reader->GetUse64BitFloats()); } // Static File Functions // Check for existence of specified file static bool IsFile(const std::string& file, bool checkGzip = true) { return (vtksys::SystemTools::FileExists(file, true) || (checkGzip && vtksys::SystemTools::FileExists(file + ".gz", true))); } //! Return file name (part beyond last /) static std::string ExtractName(const std::string& path) { auto pos = vtkFoamFile::rfind_slash(path); if (pos == std::string::npos) { // no slash return path; } else if (pos + 1 == path.length()) { // trailing slash const auto endPos = pos; pos = vtkFoamFile::rfind_slash(path, pos - 1); if (pos == std::string::npos) { // no further slash return path.substr(0, endPos); } else { return path.substr(pos + 1, endPos - pos - 1); } } else { return path.substr(pos + 1); } } //! Return directory path name (part before last /). Return includes trailing slash! static std::string ExtractPath(const std::string& path) { const auto pos = vtkFoamFile::rfind_slash(path); return pos == std::string::npos ? vtkFoamFile::cwd_string() : path.substr(0, pos + 1); } // Member Functions std::string ExpandPath(const std::string& pathIn, const std::string& defaultPath) { std::string expandedPath; bool isExpanded = false, wasPathSeparator = true; size_t charI = 0; const size_t nChars = pathIn.length(); std::string::size_type delim = 0; if ('<' == pathIn[0] && (delim = pathIn.find(">/")) != std::string::npos) { // Expand a leading / // Convenient for frequently used directories - see OpenFOAM stringOps.C // // Handle // / => FOAM_CASE directory // / => FOAM_CASE/constant directory // / => FOAM_CASE/system directory // / => not handled const std::string tag(pathIn, 1, delim - 2); if (tag == "case") { expandedPath = this->CasePath + '/'; isExpanded = true; wasPathSeparator = false; } else if (tag == "constant" || tag == "system") { expandedPath = this->CasePath + '/' + tag + '/'; isExpanded = true; wasPathSeparator = false; } // in not handled if (isExpanded) { charI = delim + 2; } } while (charI < nChars) { const char c = pathIn[charI]; switch (c) { case '$': // $-variable expansion { std::string variable; while (++charI < nChars && (isalnum(pathIn[charI]) || pathIn[charI] == '_')) { variable += pathIn[charI]; } if (variable == "FOAM_CASE") // discard path until the variable { expandedPath = this->CasePath; wasPathSeparator = true; isExpanded = true; } else if (variable == "FOAM_CASENAME") { // FOAM_CASENAME is the final directory name from CasePath expandedPath += vtkFoamFile::ExtractName(this->CasePath); wasPathSeparator = false; isExpanded = true; } else { std::string value; if (vtksys::SystemTools::GetEnv(variable, value)) { expandedPath += value; } const auto len = expandedPath.length(); if (len > 0) { const char c2 = expandedPath[len - 1]; wasPathSeparator = (c2 == '/' || c2 == '\\'); } else { wasPathSeparator = false; } } } break; case '~': // home directory expansion // not using vtksys::SystemTools::ConvertToUnixSlashes() for // a bit better handling of "~" if (wasPathSeparator) { std::string userName; while (++charI < nChars && (pathIn[charI] != '/' && pathIn[charI] != '\\') && pathIn[charI] != '$') { userName += pathIn[charI]; } std::string homeDir; if (userName.empty()) { if (!vtksys::SystemTools::GetEnv("HOME", homeDir) || homeDir.empty()) { #if defined(_WIN32) && !defined(__CYGWIN__) || defined(__LIBCATAMOUNT__) // No fallback homeDir.clear(); #else const struct passwd* pwentry = getpwuid(getuid()); if (pwentry == nullptr) { this->ThrowStackTrace("Home directory path not found"); } homeDir = pwentry->pw_dir; #endif } expandedPath = homeDir; } else if (userName == "OpenFOAM") { // So far only "~/.OpenFOAM" expansion is supported if (!vtksys::SystemTools::GetEnv("HOME", homeDir) || homeDir.empty()) { #if defined(_WIN32) && !defined(__CYGWIN__) || defined(__LIBCATAMOUNT__) // No fallback homeDir.clear(); #else const struct passwd* pwentry = getpwuid(getuid()); if (pwentry == nullptr) { this->ThrowStackTrace("Home directory path not found"); } homeDir = pwentry->pw_dir; #endif } if (homeDir.empty()) { expandedPath = homeDir; } else { expandedPath = homeDir + "/.OpenFOAM"; } } else { #if defined(_WIN32) && !defined(__CYGWIN__) || defined(__LIBCATAMOUNT__) if (!vtksys::SystemTools::GetEnv("HOME", homeDir)) { // No fallback homeDir.clear(); } expandedPath = vtkFoamFile::ExtractPath(homeDir) + userName; #else const struct passwd* pwentry = getpwnam(userName.c_str()); if (pwentry == nullptr) { this->ThrowStackTrace("No home directory for user " + userName); } expandedPath = pwentry->pw_dir; #endif } wasPathSeparator = false; isExpanded = true; break; } VTK_FALLTHROUGH; default: wasPathSeparator = (c == '/' || c == '\\'); expandedPath += c; charI++; } } if (isExpanded || expandedPath[0] == '/' || expandedPath[0] == '\\') { return expandedPath; } else { return defaultPath + expandedPath; } } void IncludeFile(const std::string& includedFileName, const std::string& defaultPath) { if (this->StackI >= VTK_FOAMFILE_INCLUDE_STACK_SIZE) { throw this->StackString() << "Exceeded maximum #include recursions of " << VTK_FOAMFILE_INCLUDE_STACK_SIZE; } // use the default bitwise copy constructor this->Stack[this->StackI++] = new vtkFoamFileStack(*this); this->Superclass::Reset(); this->Open(this->ExpandPath(includedFileName, defaultPath)); } // the tokenizer // returns true if success, false if encountered EOF bool Read(vtkFoamToken& token) { token.SetStreamOption(this->GetStreamOption()); const bool use64BitLabels = this->IsLabel64(); // expanded the outermost loop in nextTokenHead() for performance int c; while (isspace(c = this->Getc())) // isspace() accepts -1 as EOF { if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (c == '/') { this->PutBack(c); c = this->NextTokenHead(); } #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD if (c != '#') { this->Superclass::WasNewline = false; } #endif constexpr int MAXLEN = 1024; char buf[MAXLEN + 1]; int charI = 0; switch (c) { case '(': case ')': // high-priority punctuation token token = static_cast(c); return true; case '1': case '2': case '3': case '4': case '5': case '6': case '7': case '8': case '9': case '0': case '-': // undetermined number token do { buf[charI++] = static_cast(c); } while (isdigit(c = this->Getc()) && charI < MAXLEN); if (c != '.' && c != 'e' && c != 'E' && charI < MAXLEN && c != EOF) { // label token buf[charI] = '\0'; if (use64BitLabels) { token = static_cast(strtoll(buf, nullptr, 10)); } else { token = static_cast(strtol(buf, nullptr, 10)); } this->PutBack(c); return true; } VTK_FALLTHROUGH; case '.': // scalar token if (c == '.' && charI < MAXLEN) { // read decimal fraction part buf[charI++] = static_cast(c); while (isdigit(c = this->Getc()) && charI < MAXLEN) { buf[charI++] = static_cast(c); } } if ((c == 'e' || c == 'E') && charI < MAXLEN) { // read exponent part buf[charI++] = static_cast(c); if (((c = this->Getc()) == '+' || c == '-') && charI < MAXLEN) { buf[charI++] = static_cast(c); c = this->Getc(); } while (isdigit(c) && charI < MAXLEN) { buf[charI++] = static_cast(c); c = this->Getc(); } } if (charI == 1 && buf[0] == '-') { token = '-'; this->PutBack(c); return true; } buf[charI] = '\0'; token = strtod(buf, nullptr); this->PutBack(c); break; case ';': case '{': case '}': case '[': case ']': case ':': case ',': case '=': case '+': case '*': case '/': // low-priority punctuation token token = static_cast(c); return true; case '"': { // string token bool wasEscape = false; while ((c = this->Getc()) != EOF && charI < MAXLEN) { if (c == '\\' && !wasEscape) { wasEscape = true; continue; } else if (c == '"' && !wasEscape) { break; } else if (c == '\n') { ++this->Superclass::LineNumber; if (!wasEscape) { this->ThrowStackTrace("Unescaped newline in string constant"); } } buf[charI++] = static_cast(c); wasEscape = false; } buf[charI] = '\0'; token = buf; } break; case EOF: // end of file token.SetBad(); return false; case '$': { vtkFoamToken identifierToken; if (!this->Read(identifierToken)) { this->ThrowStackTrace("Unexpected EOF reading identifier"); } if (identifierToken.GetType() != vtkFoamToken::STRING) { throw this->StackString() << "Expected a word, found " << identifierToken; } token.SetIdentifier(identifierToken.ToString()); return true; } case '#': { #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD // the OpenFOAM #-directives can indeed be placed in the // middle of a line if (!this->Superclass::WasNewline) { this->ThrowStackTrace("Encountered #-directive in the middle of a line"); } this->Superclass::WasNewline = false; #endif // read directive vtkFoamToken directiveToken; if (!this->Read(directiveToken)) { this->ThrowStackTrace("Unexpected EOF reading directive"); } if (directiveToken == "include") { vtkFoamToken fileNameToken; if (!this->Read(fileNameToken)) { this->ThrowStackTrace("Unexpected EOF reading filename"); } this->IncludeFile(fileNameToken.ToString(), vtkFoamFile::ExtractPath(this->FileName)); } else if (directiveToken == "sinclude" || directiveToken == "includeIfPresent") { vtkFoamToken fileNameToken; if (!this->Read(fileNameToken)) { this->ThrowStackTrace("Unexpected EOF reading filename"); } // special treatment since the file is allowed to be missing const std::string fullName = this->ExpandPath(fileNameToken.ToString(), vtkFoamFile::ExtractPath(this->FileName)); FILE* fh = vtksys::SystemTools::Fopen(fullName, "rb"); if (fh) { fclose(fh); this->IncludeFile(fileNameToken.ToString(), vtkFoamFile::ExtractPath(this->FileName)); } } else if (directiveToken == "inputMode") { vtkFoamToken modeToken; if (!this->Read(modeToken)) { this->ThrowStackTrace("Unexpected EOF reading inputMode specifier"); } if (modeToken == "merge" || modeToken == "default") { this->InputMode = INPUT_MODE_MERGE; } else if (modeToken == "overwrite") { this->InputMode = INPUT_MODE_OVERWRITE; } else if (modeToken == "protect") { // not properly supported - treat like "merge" for now // this->InputMode = INPUT_MODE_PROTECT; this->InputMode = INPUT_MODE_MERGE; } else if (modeToken == "warn") { // not properly supported - treat like "error" for now // this->InputMode = INPUT_MODE_WARN; this->InputMode = INPUT_MODE_ERROR; } else if (modeToken == "error") { this->InputMode = INPUT_MODE_ERROR; } else { throw this->StackString() << "Expected one of inputMode specifiers " "(merge, overwrite, protect, warn, error, default), found " << modeToken; } } else if (directiveToken == '{') { // '#{' verbatim/code block. swallow everything until a closing '#}' // This hopefully matches the first one... while (true) { c = this->NextTokenHead(); if (c == EOF) { this->ThrowStackTrace("Unexpected EOF while skipping over #{ directive"); } else if (c == '#') { c = this->Getc(); if (c == '/') { this->PutBack(c); } else if (c == '}') { break; } } } } else { throw this->StackString() << "Unsupported directive " << directiveToken; } return this->Read(token); } default: { // parses as a word token, but gives the STRING type for simplicity int inBrace = 0; do { if (c == '(') { inBrace++; } else if (c == ')' && --inBrace == -1) { break; } buf[charI++] = static_cast(c); // valid characters that constitutes a word // cf. src/OpenFOAM/primitives/strings/word/wordI.H } while ((c = this->Getc()) != EOF && !isspace(c) && c != '"' && c != '/' && c != ';' && c != '{' && c != '}' && charI < MAXLEN); buf[charI] = '\0'; token = buf; this->PutBack(c); } } if (c == EOF) { this->ThrowUnexpectedEOFException(); } if (charI == MAXLEN) { throw this->StackString() << "Exceeded maximum allowed length of " << MAXLEN; } return true; } // fread or gzread with buffering handling vtkTypeInt64 Read(unsigned char* buf, const vtkTypeInt64 len) { const size_t buflen = (this->Superclass::BufEndPtr - this->Superclass::BufPtr); vtkTypeInt64 readlen; if (static_cast(len) > buflen) { memcpy(buf, this->Superclass::BufPtr, buflen); this->InflateNext(buf + buflen, len - buflen, &readlen); if (readlen >= 0) { readlen += buflen; } else { if (buflen == 0) // return EOF { readlen = -1; } else { readlen = buflen; } } this->Superclass::BufPtr = this->Superclass::BufEndPtr; } else { memcpy(buf, this->Superclass::BufPtr, len); this->Superclass::BufPtr += len; readlen = len; } for (vtkTypeInt64 i = 0; i < readlen; ++i) { if (buf[i] == '\n') { ++this->Superclass::LineNumber; } } return readlen; } void ReadExpecting(const char expected) { // skip prepending invalid chars // expanded the outermost loop in nextTokenHead() for performance int c; while (isspace(c = this->Getc())) // isspace() accepts -1 as EOF { if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (c == '/') { this->PutBack(c); c = this->NextTokenHead(); } if (c != expected) { this->ThrowUnexpectedTokenException(expected, c); } } void ReadExpecting(const char* str) { vtkFoamToken tok; if (!this->Read(tok) || tok != str) { throw this->StackString() << "Expected string \"" << str << "\", found " << tok; } } // ASCII read of longest integer vtkTypeInt64 ReadIntegerValue(); // ASCII read of longest floating-point double ReadDoubleValue(); }; //------------------------------------------------------------------------------ // Code: vtkFoamFile int vtkFoamFile::ReadNext() { if (!this->InflateNext(this->Superclass::Outbuf + 1, VTK_FOAMFILE_OUTBUFSIZE)) { return this->CloseIncludedFile() ? this->Getc() : EOF; } return *this->Superclass::BufPtr++; } // specialized for reading an integer value. // not using the standard strtol() for speed reason. vtkTypeInt64 vtkFoamFile::ReadIntegerValue() { // skip prepending invalid chars // expanded the outermost loop in nextTokenHead() for performance int c; while (isspace(c = this->Getc())) // isspace() accepts -1 as EOF { if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (c == '/') { this->PutBack(c); c = this->NextTokenHead(); } // leading sign? const bool negNum = (c == '-'); if (negNum || c == '+') { c = this->Getc(); if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (!isdigit(c)) // isdigit() accepts -1 as EOF { if (c == EOF) { this->ThrowUnexpectedEOFException(); } else { this->ThrowUnexpectedNondigitException(c); } } vtkTypeInt64 num = c - '0'; while (isdigit(c = this->Getc())) { num = 10 * num + c - '0'; } if (c == EOF) { this->ThrowUnexpectedEOFException(); } this->PutBack(c); return negNum ? -num : num; } // extremely simplified high-performing string to floating point // conversion code based on // ParaView3/VTK/Utilities/vtksqlite/vtk_sqlite3.c double vtkFoamFile::ReadDoubleValue() { // skip prepending invalid chars // expanded the outermost loop in nextTokenHead() for performance int c; while (isspace(c = this->Getc())) // isspace() accepts -1 as EOF { if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (c == '/') { this->PutBack(c); c = this->NextTokenHead(); } // leading sign? const bool negNum = (c == '-'); if (negNum || c == '+') { c = this->Getc(); if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (!isdigit(c) && c != '.') // Attention: isdigit() accepts EOF { this->ThrowUnexpectedNondigitException(c); } double num = 0; // read integer part (before '.') if (c != '.') { num = c - '0'; while (isdigit(c = this->Getc())) { num = num * 10.0 + (c - '0'); } } // read decimal part (after '.') if (c == '.') { double divisor = 1.0; while (isdigit(c = this->Getc())) { num = num * 10.0 + (c - '0'); divisor *= 10.0; } num /= divisor; } // read exponent part if (c == 'E' || c == 'e') { int esign = 1; int eval = 0; double scale = 1.0; c = this->Getc(); if (c == '-') { esign = -1; c = this->Getc(); } else if (c == '+') { c = this->Getc(); } while (isdigit(c)) { eval = eval * 10 + (c - '0'); c = this->Getc(); } // fast exponent multiplication! while (eval >= 64) { scale *= 1.0e+64; eval -= 64; } while (eval >= 16) { scale *= 1.0e+16; eval -= 16; } while (eval >= 4) { scale *= 1.0e+4; eval -= 4; } while (eval >= 1) { scale *= 1.0e+1; eval -= 1; } if (esign < 0) { num /= scale; } else { num *= scale; } } if (c == EOF) { this->ThrowUnexpectedEOFException(); } this->PutBack(c); return negNum ? -num : num; } void vtkFoamFile::ThrowStackTrace(const std::string& msg) { throw this->StackString() << msg; } // hacks to keep exception throwing code out-of-line to make // putBack() and readExpecting() inline expandable void vtkFoamFile::ThrowUnexpectedEOFException() { this->ThrowStackTrace("Unexpected EOF"); } void vtkFoamFile::ThrowUnexpectedNondigitException(int c) { throw this->StackString() << "Expected a number, found a non-digit character " << static_cast(c); } void vtkFoamFile::ThrowUnexpectedTokenException(char expected, int c) { vtkFoamError err; err << this->StackString() << "Expected punctuation token '" << expected << "', found "; if (c == EOF) { err << "EOF"; } else { err << static_cast(c); } throw err; } void vtkFoamFile::ThrowDuplicatedPutBackException() { this->ThrowStackTrace("Attempted duplicated putBack()"); } bool vtkFoamFile::InflateNext(unsigned char* buf, size_t requestSize, vtkTypeInt64* readSize) { if (readSize) { *readSize = -1; // Set to an error state for early returns } size_t size; if (this->Superclass::IsCompressed) { if (this->Superclass::ZStatus != Z_OK) { return false; } this->Superclass::Z.next_out = buf; this->Superclass::Z.avail_out = static_cast(requestSize); do { if (this->Superclass::Z.avail_in == 0) { this->Superclass::Z.next_in = this->Superclass::Inbuf; this->Superclass::Z.avail_in = static_cast( fread(this->Superclass::Inbuf, 1, VTK_FOAMFILE_INBUFSIZE, this->Superclass::File)); if (ferror(this->Superclass::File)) { this->ThrowStackTrace("failed in fread()"); } } this->Superclass::ZStatus = inflate(&this->Superclass::Z, Z_NO_FLUSH); if (this->Superclass::ZStatus == Z_STREAM_END #if VTK_FOAMFILE_OMIT_CRCCHECK // the dummy CRC function causes data error when finalizing // so we have to proceed even when a data error is detected || this->Superclass::ZStatus == Z_DATA_ERROR #endif ) { break; } if (this->Superclass::ZStatus != Z_OK) { throw this->StackString() << "Inflation failed: " << (this->Superclass::Z.msg ? this->Superclass::Z.msg : ""); } } while (this->Superclass::Z.avail_out > 0); size = requestSize - this->Superclass::Z.avail_out; } else { // not compressed size = fread(buf, 1, requestSize, this->Superclass::File); } if (size <= 0) { // retain the current location bufPtr_ to the end of the buffer so that // getc() returns EOF again when called next time return false; } // size > 0 // reserve the first byte for getback char this->Superclass::BufPtr = this->Superclass::Outbuf + 1; this->Superclass::BufEndPtr = this->Superclass::BufPtr + size; if (readSize) { // Cast size_t to int64. Should be OK since requestSize came from OpenFOAM (signed integer) *readSize = static_cast(size); } return true; } // get next semantically valid character int vtkFoamFile::NextTokenHead() { while (true) { int c; while (isspace(c = this->Getc())) // isspace() accepts -1 as EOF { if (c == '\n') { ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } } if (c == '/') { if ((c = this->Getc()) == '/') { while ((c = this->Getc()) != EOF && c != '\n') ; if (c == EOF) { return c; } ++this->Superclass::LineNumber; #if VTK_FOAMFILE_RECOGNIZE_LINEHEAD this->Superclass::WasNewline = true; #endif } else if (c == '*') { while (true) { while ((c = this->Getc()) != EOF && c != '*') { if (c == '\n') { ++this->Superclass::LineNumber; } } if (c == EOF) { return c; } else if ((c = this->Getc()) == '/') { break; } this->PutBack(c); } } else { this->PutBack(c); // may be an EOF return '/'; } } else // may be an EOF { return c; } } #if defined(__hpux) return EOF; // this line should not be executed; workaround for HP-UXia64-aCC #endif } //------------------------------------------------------------------------------ // class vtkFoamDict // A class to holds a FoamFile data structure struct vtkFoamDict : public std::vector { private: typedef std::vector Superclass; vtkFoamToken Token; const vtkFoamDict* UpperDictPtr; vtkFoamDict(const vtkFoamDict&) = delete; public: // Default construct with given parent directory explicit vtkFoamDict(const vtkFoamDict* upperDictPtr = nullptr) : UpperDictPtr(upperDictPtr) { } // Copy construct with given parent directory vtkFoamDict(const vtkFoamDict& dict, const vtkFoamDict* upperDictPtr); // Destructor - delete list elements if held by the dictionary ~vtkFoamDict(); // Remove top element, deleting its pointer void remove_back(); void SetStreamOption(const vtkFoamStreamOption& opt) { this->Token.SetStreamOption(opt); } bool IsLabel64() const { return this->Token.IsLabel64(); } // convenience const vtkFoamToken& GetToken() const { return this->Token; } const vtkFoamDict* GetUpperDictPtr() const { return this->UpperDictPtr; } vtkFoamToken::tokenType GetType() const { return (this->Token.GetType() == vtkFoamToken::UNDEFINED ? vtkFoamToken::DICTIONARY : this->Token.GetType()); } // Return list of keywords - table of contents std::vector Toc() const; // Search dictionary for specified keyword. Return nullptr on failure. vtkFoamEntry* Lookup(const std::string& keyword, bool isPattern = false) const; // Reads a FoamFile or a subdictionary. // If the stream to be read is a subdictionary, // the preceding '{' is assumed to have already been discarded. bool Read( vtkFoamIOobject& io, bool isSubDict = false, const vtkFoamToken& firstKeyword = vtkFoamToken()); }; //------------------------------------------------------------------------------ // class vtkFoamIOobject // Extends vtkFoamFile with OpenFOAM class/object information struct vtkFoamIOobject : public vtkFoamFile { private: typedef vtkFoamFile Superclass; std::string objectName_; std::string headerClassName_; vtkFoamError error_; // Inform IO object that lagrangian/positions has extra data (OpenFOAM v1.4 - v2.4) const bool LagrangianPositionsExtraData_; // Reads OpenFOAM format/class/object information and handles "arch" information void ReadHeader(); // Attempt to open file (or file.gz) and read header bool OpenFile(const std::string& file, bool checkGzip = false) { this->ClearError(); // Discard any previous errors try { this->Superclass::Open(file); checkGzip = false; } catch (const vtkFoamError& err) { if (checkGzip) { // Avoid checking again if already ends_with(".gz") const auto len = file.length(); if (len > 3 && file.compare(len - 3, std::string::npos, ".gz") == 0) { checkGzip = false; } } if (!checkGzip) { this->SetError(err); return false; } } if (checkGzip) { try { this->Superclass::Open(file + ".gz"); } catch (const vtkFoamError& err) { this->SetError(err); return false; } } try { this->ReadHeader(); } catch (const vtkFoamError& err) { this->Superclass::Close(); this->SetError(err); return false; } return true; } void CloseFile() { this->Superclass::Close(); this->objectName_.clear(); this->headerClassName_.clear(); this->error_.clear(); } public: // No generated methods vtkFoamIOobject() = delete; vtkFoamIOobject(const vtkFoamIOobject&) = delete; void operator=(const vtkFoamIOobject&) = delete; // Construct for specified case -path vtkFoamIOobject(const std::string& casePath, vtkOpenFOAMReader* reader) : vtkFoamFile(casePath, reader) , LagrangianPositionsExtraData_(static_cast(!reader->GetPositionsIsIn13Format())) { } ~vtkFoamIOobject() { this->Close(); } // Attempt to open file (without gzip fallback) and read FoamFile header bool Open(const std::string& file) { return this->OpenFile(file, false); } // Attempt to open file (with gzip fallback) and read FoamFile header bool OpenOrGzip(const std::string& file) { return this->OpenFile(file, true); } // Attempt to open file relative to the case, and read FoamFile header bool OpenCaseFile(const std::string& file) { return this->OpenFile(this->GetCasePath() + file); } // Attempt to open file (or file.gz) relative to the case, and read FoamFile header bool OpenCaseFileOrGzip(const std::string& file) { return this->OpenFile(this->GetCasePath() + file); } void Close() { this->CloseFile(); } const std::string& GetClassName() const noexcept { return this->headerClassName_; } const std::string& GetObjectName() const noexcept { return this->objectName_; } const vtkFoamError& GetError() const noexcept { return this->error_; } void ClearError() noexcept { this->error_.clear(); } bool HasError() const noexcept { return !this->error_.empty(); } void SetError(const vtkFoamError& err) { this->error_ = err; } bool GetLagrangianPositionsExtraData() const { return this->LagrangianPositionsExtraData_; } }; //------------------------------------------------------------------------------ // ASCII read of primitive, with type casting template struct vtkFoamReadValue { static T ReadValue(vtkFoamIOobject& io); }; template <> inline vtkTypeInt8 vtkFoamReadValue::ReadValue(vtkFoamIOobject& io) { return static_cast(io.ReadIntegerValue()); } template <> inline vtkTypeInt32 vtkFoamReadValue::ReadValue(vtkFoamIOobject& io) { return static_cast(io.ReadIntegerValue()); } template <> inline vtkTypeInt64 vtkFoamReadValue::ReadValue(vtkFoamIOobject& io) { return io.ReadIntegerValue(); } template <> inline float vtkFoamReadValue::ReadValue(vtkFoamIOobject& io) { return static_cast(io.ReadDoubleValue()); } template <> inline double vtkFoamReadValue::ReadValue(vtkFoamIOobject& io) { return io.ReadDoubleValue(); } //------------------------------------------------------------------------------ // struct vtkFoamRead for reading primitives, lists etc. struct vtkFoamRead { // -------------------------------------------------------------------------- // Reading lists of primitives (int/float/...) template class listTraits { listT* Ptr; public: using ValueType = typename listT::ValueType; listTraits() : Ptr(listT::New()) { } // Get the contained pointer listT* GetPointer() const noexcept { return this->Ptr; } // De-reference pointer for operation listT* operator->() const noexcept { return this->Ptr; } void ReadValue(vtkFoamIOobject&, const vtkFoamToken& currToken) { if (!currToken.Is()) { throw vtkFoamError() << "Expected an integer or a (, found " << currToken; } this->Ptr->InsertNextValue(currToken.To()); } void ReadUniformValues(vtkFoamIOobject& io) { primitiveT value = vtkFoamReadValue::ReadValue(io); this->Ptr->FillValue(value); } void ReadAsciiList(vtkFoamIOobject& io) { const vtkIdType nTuples = this->Ptr->GetNumberOfTuples(); for (vtkIdType i = 0; i < nTuples; ++i) { this->Ptr->SetValue(i, vtkFoamReadValue::ReadValue(io)); } } void ReadBinaryList(vtkFoamIOobject& io) { // nComponents == 1 const vtkIdType nTuples = this->Ptr->GetNumberOfTuples(); const size_t nbytes = (nTuples * sizeof(primitiveT)); if (typeid(ValueType) == typeid(primitiveT)) { io.Read(reinterpret_cast(this->Ptr->GetPointer(0)), nbytes); } else { auto* fileData = vtkDataArray::CreateDataArray(vtkTypeTraits::VTKTypeID()); // nComponents == 1 fileData->SetNumberOfTuples(nTuples); io.Read(reinterpret_cast(fileData->GetVoidPointer(0)), nbytes); this->Ptr->DeepCopy(fileData); fileData->Delete(); } } }; // -------------------------------------------------------------------------- // Reads rank 1 lists of types vector, sphericalTensor, symmTensor and tensor. // If isPositions is true it reads Cloud type of data as // particle positions. cf. (the positions format) // src/lagrangian/basic/particle/particleIO.C - writePosition() template class vectorListTraits { listT* Ptr; // Items to skip for lagrangian/positions (class Cloud) after the x/y/z values. // xyz (3*scalar) + celli (label) // in OpenFOAM 1.4 -> 2.4 also had facei (label) and stepFraction (scalar) // ASCII only static void LagrangianPositionsSkip(vtkFoamIOobject& io) { (void)io.ReadIntegerValue(); // Skip celli (label) if (io.GetLagrangianPositionsExtraData()) { (void)io.ReadIntegerValue(); // Skip facei (label) (void)io.ReadDoubleValue(); // Skip stepFraction (scalar) } } public: using ValueType = typename listT::ValueType; vectorListTraits() : Ptr(listT::New()) { this->Ptr->SetNumberOfComponents(nComponents); } // Get the contained pointer listT* GetPointer() const noexcept { return this->Ptr; } // De-reference pointer for operation listT* operator->() const noexcept { return this->Ptr; } void ReadValue(vtkFoamIOobject& io, const vtkFoamToken& currToken) { if (currToken != '(') { throw vtkFoamError() << "Expected '(', found " << currToken; } primitiveT tuple[nComponents]; for (int cmpt = 0; cmpt < nComponents; ++cmpt) { tuple[cmpt] = vtkFoamReadValue::ReadValue(io); } ::remapFoamTuple(tuple); // For symmTensor io.ReadExpecting(')'); this->Ptr->InsertNextTuple(tuple); } void ReadUniformValues(vtkFoamIOobject& io) { const vtkIdType nTuples = this->Ptr->GetNumberOfTuples(); io.ReadExpecting('('); primitiveT tuple[nComponents]; for (int cmpt = 0; cmpt < nComponents; ++cmpt) { tuple[cmpt] = vtkFoamReadValue::ReadValue(io); } ::remapFoamTuple(tuple); // For symmTensor io.ReadExpecting(')'); if (isPositions) { this->LagrangianPositionsSkip(io); } for (vtkIdType i = 0; i < nTuples; ++i) { this->Ptr->SetTuple(i, tuple); } } void ReadAsciiList(vtkFoamIOobject& io) { const vtkIdType nTuples = this->Ptr->GetNumberOfTuples(); for (vtkIdType i = 0; i < nTuples; ++i) { io.ReadExpecting('('); ValueType* tuple = this->Ptr->GetPointer(nComponents * i); for (int cmpt = 0; cmpt < nComponents; ++cmpt) { tuple[cmpt] = static_cast(vtkFoamReadValue::ReadValue(io)); } ::remapFoamTuple(tuple); // For symmTensor io.ReadExpecting(')'); if (isPositions) { this->LagrangianPositionsSkip(io); } } } void ReadBinaryList(vtkFoamIOobject& io) { const vtkTypeInt64 nTuples = this->Ptr->GetNumberOfTuples(); if (isPositions) // lagrangian/positions (class Cloud) { // xyz (3*scalar) + celli (label) // in OpenFOAM 1.4 -> 2.4 also had facei (label) and stepFraction (scalar) const unsigned labelWidth = (io.IsLabel64() ? 8 : 4); const unsigned tupleLength = (sizeof(primitiveT) * nComponents + labelWidth + (io.GetLagrangianPositionsExtraData() ? (labelWidth + sizeof(primitiveT)) : 0)); // Variable-sized stack arrays are non-standard, so use our own version // Have a max of 6 double/int64 elements = 48 bytes. Allocate 64 bytes for good measure vtkFoamStackVector buffer(tupleLength); primitiveT* tuple = reinterpret_cast(buffer.data()); for (vtkTypeInt64 i = 0; i < nTuples; ++i) { io.ReadExpecting('('); io.Read(reinterpret_cast(tuple), tupleLength); io.ReadExpecting(')'); this->Ptr->SetTuple(i, tuple); } } else { // Compiler hint for better unrolling: VTK_ASSUME(this->Ptr->GetNumberOfComponents() == nComponents); const unsigned tupleLength = (sizeof(primitiveT) * nComponents); primitiveT tuple[nComponents]; for (vtkTypeInt64 i = 0; i < nTuples; ++i) { const int readLength = io.Read(reinterpret_cast(tuple), tupleLength); if (readLength != tupleLength) { throw vtkFoamError() << "Failed to read tuple " << i << '/' << nTuples << ": Expected " << tupleLength << " bytes, got " << readLength << " bytes."; } ::remapFoamTuple(tuple); // For symmTensor for (int cmpt = 0; cmpt < nComponents; ++cmpt) { this->Ptr->SetTypedComponent(i, cmpt, static_cast(tuple[cmpt])); } } } } }; }; //------------------------------------------------------------------------------ // class vtkFoamEntryValue // a class that represents a value of a dictionary entry that corresponds to // its keyword. note that an entry can have more than one value. struct vtkFoamEntryValue : public vtkFoamToken { private: typedef vtkFoamToken Superclass; bool IsUniformEntry; bool Managed; const vtkFoamEntry* UpperEntryPtr; vtkFoamEntryValue() = delete; vtkObjectBase* ToVTKObject() { return this->Superclass::VtkObjectPtr; } // Delete if managed void Clear(); void ReadList(vtkFoamIOobject& io); public: // Construct empty value with given parent explicit vtkFoamEntryValue(const vtkFoamEntry* parent) : IsUniformEntry(false) , Managed(true) , UpperEntryPtr(parent) { } // Copy construct vtkFoamEntryValue(vtkFoamEntryValue& val, const vtkFoamEntry* parent); ~vtkFoamEntryValue() { this->Clear(); } // Member Functions void SetEmptyList() { this->Clear(); this->IsUniformEntry = false; this->Superclass::Type = vtkFoamToken::EMPTYLIST; } bool IsUniform() const noexcept { return this->IsUniformEntry; } bool Read(vtkFoamIOobject& io); void ReadDictionary(vtkFoamIOobject& io, const vtkFoamToken& firstKeyword); const vtkDataArray& LabelList() const { return *this->Superclass::LabelListPtr; } vtkDataArray& LabelList() { return *this->Superclass::LabelListPtr; } const vtkFoamLabelListList& LabelListList() const { return *this->Superclass::LabelListListPtr; } const vtkFloatArray& ScalarList() const { return *this->Superclass::ScalarListPtr; } vtkFloatArray& ScalarList() { return *this->Superclass::ScalarListPtr; } const vtkFloatArray& VectorList() const { return *this->Superclass::VectorListPtr; } vtkFloatArray& VectorList() { return *this->Superclass::VectorListPtr; } const vtkStringArray& StringList() const { return *this->Superclass::StringListPtr; } vtkStringArray& StringList() { return *this->Superclass::StringListPtr; } const vtkFoamDict& Dictionary() const { return *this->Superclass::DictPtr; } vtkFoamDict& Dictionary() { return *this->Superclass::DictPtr; } // Release the managed pointer, cast to specified type template DataType* ReleasePtr() { // Not managed = do not delete pointer in destructor this->Managed = false; return static_cast(this->Superclass::AnyPointer); } std::string ToString() const { return this->Superclass::Type == STRING ? this->Superclass::ToString() : std::string(); } float ToFloat() const { return this->Superclass::IsNumeric() ? this->Superclass::To() : 0.0F; } double ToDouble() const { return this->Superclass::IsNumeric() ? this->Superclass::To() : 0.0; } // TODO is it ok to always use a 64bit int here? vtkTypeInt64 ToInt() const { return this->Superclass::Type == LABEL ? this->Superclass::To() : 0; } void MakeLabelList(const vtkIdType len, const vtkTypeInt64 val = 0) { this->Superclass::Type = vtkFoamToken::LABELLIST; if (this->IsLabel64()) { auto* array = vtkTypeInt64Array::New(); array->SetNumberOfValues(len); array->FillValue(val); this->Superclass::LabelListPtr = array; } else { auto* array = vtkTypeInt32Array::New(); array->SetNumberOfValues(len); array->FillValue(static_cast(val)); this->Superclass::LabelListPtr = array; } } void MakeScalarList(const vtkIdType len, const float val = 0.0f) { this->Superclass::Type = vtkFoamToken::SCALARLIST; this->Superclass::ScalarListPtr = vtkFloatArray::New(); this->Superclass::ScalarListPtr->SetNumberOfValues(len); this->Superclass::ScalarListPtr->FillValue(val); } template void ReadNonUniformList(vtkFoamIOobject& io); // Dispatch reading of uniform list based on the field data type (scalar, vector etc). // Return false if could not be dispatched bool ReadNonUniformList(vtkFoamIOobject& io, vtkFoamTypes::dataType fieldDataType); bool ReadField(vtkFoamIOobject& io); // Read a list of labelLists. requires size prefix of the listList // to be present. size of each sublist must also be present in the // stream if the format is binary. void ReadLabelListList(vtkFoamIOobject& io); // Read compact labelListList which has offsets and data void ReadCompactLabelListList(vtkFoamIOobject& io); // Read dimensions set (always ASCII). The leading '[' has already been removed before calling. // - can be integer or floating point // - user-generated files may have only the first five dimensions. // Note // - may even have "human-readable" values such as [kg m^-1 s^-2] but they are very rare // and we silently skip these void ReadDimensionSet(vtkFoamIOobject& io); }; //------------------------------------------------------------------------------ // Code: vtkFoamEntryValue // generic reader for nonuniform lists. requires size prefix of the // list to be present in the stream if the format is binary. template void vtkFoamEntryValue::ReadNonUniformList(vtkFoamIOobject& io) { this->SetStreamOption(io); vtkFoamToken currToken; currToken.SetStreamOption(io); if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } traitsT list; this->Superclass::Type = listType; this->Superclass::VtkObjectPtr = list.GetPointer(); if (currToken.Is()) { const vtkTypeInt64 size = currToken.To(); if (size < 0) { throw vtkFoamError() << "List size must not be negative: size = " << size; } list->SetNumberOfTuples(size); if (io.IsAsciiFormat()) { if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } // some objects have lists with only one element enclosed by {} // e. g. simpleFoam/pitzDaily3Blocks/constant/polyMesh/faceZones if (currToken == '{') { list.ReadUniformValues(io); io.ReadExpecting('}'); return; } else if (currToken != '(') { throw vtkFoamError() << "Expected '(', found " << currToken; } list.ReadAsciiList(io); io.ReadExpecting(')'); } else { if (size > 0) { // Non-empty (binary) list - only read parentheses only when size > 0 io.ReadExpecting('('); list.ReadBinaryList(io); io.ReadExpecting(')'); } } } else if (currToken == '(') { while (io.Read(currToken) && currToken != ')') { list.ReadValue(io, currToken); } list->Squeeze(); } else { throw vtkFoamError() << "Expected integer or '(', found " << currToken; } } // Dispatch known field/list types for non-uniform reading bool vtkFoamEntryValue::ReadNonUniformList(vtkFoamIOobject& io, vtkFoamTypes::dataType listDataType) { bool handled = true; switch (listDataType) { case vtkFoamTypes::BOOL_TYPE: { // List is read as a list of bytes (binary) or ints (ascii) // - primary location is the flipMap entry in faceZones this->ReadNonUniformList>(io); break; } case vtkFoamTypes::LABEL_TYPE: { if (io.IsLabel64()) { this->ReadNonUniformList>(io); } else { this->ReadNonUniformList>(io); } break; } case vtkFoamTypes::SCALAR_TYPE: { if (io.IsFloat64()) { this->ReadNonUniformList>(io); } else { this->ReadNonUniformList>(io); } break; } case vtkFoamTypes::SPH_TENSOR_TYPE: { if (io.IsFloat64()) { this->ReadNonUniformList>(io); } else { this->ReadNonUniformList>(io); } break; } case vtkFoamTypes::VECTOR_TYPE: { if (io.IsFloat64()) { this->ReadNonUniformList>(io); } else { this->ReadNonUniformList>(io); } break; } case vtkFoamTypes::SYMM_TENSOR_TYPE: { if (io.IsFloat64()) { this->ReadNonUniformList>(io); } else { this->ReadNonUniformList>(io); } break; } case vtkFoamTypes::TENSOR_TYPE: { if (io.IsFloat64()) { this->ReadNonUniformList>(io); } else { this->ReadNonUniformList>(io); } break; } default: { handled = false; break; } } return handled; } bool vtkFoamEntryValue::ReadField(vtkFoamIOobject& io) { this->SetStreamOption(io); // Basic field types: "boolField", "labelField", "scalarField" ... vtkFoamTypes::dataType listDataType(vtkFoamTypes::FieldToEnum(io.GetClassName())); try { if (vtkFoamTypes::IsGood(listDataType)) { this->ReadNonUniformList(io, listDataType); } else { throw vtkFoamError() << "Unsupported field type " << io.GetClassName(); } } catch (const vtkFoamError& err) { io.SetError(err); return false; } return true; } // Read a list of labelLists. requires size prefix of the listList // to be present. size of each sublist must also be present in the // stream if the format is binary. void vtkFoamEntryValue::ReadLabelListList(vtkFoamIOobject& io) { // NOTE: // when OpenFOAM writes a "faceCompactList" it automatically switches to ASCII // if it detects that the offsets will overflow 32bits. // // We risk the same overflow potential when constructing a compact labelListList. // Thus assume the worst and use 64bit sizing when reading ASCII. const bool use64BitLabels = (io.IsLabel64() || io.IsAsciiFormat()); vtkFoamToken currToken; currToken.SetStreamOption(io); if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } if (currToken.IsLabel()) { const vtkTypeInt64 listLen = currToken.To(); if (listLen < 0) { throw vtkFoamError() << "Illegal negative list length: " << listLen; } if (use64BitLabels) { this->LabelListListPtr = new vtkFoamLabelListList64; } else { this->LabelListListPtr = new vtkFoamLabelListList32; } // Initial guess for list length this->LabelListListPtr->ResizeExact(listLen, 4 * listLen); this->Superclass::Type = vtkFoamToken::LABELLISTLIST; io.ReadExpecting('('); vtkIdType nTotalElems = 0; for (vtkTypeInt64 idx = 0; idx < listLen; ++idx) { if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } if (currToken.IsLabel()) { const vtkTypeInt64 sublistLen = currToken.To(); if (sublistLen < 0) { throw vtkFoamError() << "Illegal negative list length: " << sublistLen; } // LabelListListPtr->SetOffset(idx, nTotalElems); void* sublist = this->LabelListListPtr->WritePointer(idx, nTotalElems, sublistLen); if (io.IsAsciiFormat()) { io.ReadExpecting('('); for (vtkTypeInt64 subIdx = 0; subIdx < sublistLen; ++subIdx) { vtkTypeInt64 value(vtkFoamReadValue::ReadValue(io)); this->LabelListListPtr->SetValue(idx, subIdx, value); } io.ReadExpecting(')'); } else if (sublistLen > 0) { // Non-empty (binary) list - only read parentheses only when size > 0 const size_t nbytes = static_cast(sublistLen * this->LabelListListPtr->GetLabelSize()); io.ReadExpecting('('); io.Read(reinterpret_cast(sublist), nbytes); io.ReadExpecting(')'); } nTotalElems += sublistLen; } else if (currToken == '(') { this->Superclass::LabelListListPtr->SetOffset(idx, nTotalElems); while (io.Read(currToken) && currToken != ')') { if (!currToken.IsLabel()) { throw vtkFoamError() << "Expected an integer, found " << currToken; } this->Superclass::LabelListListPtr->InsertValue(nTotalElems++, currToken.To()); ++nTotalElems; } } else { throw vtkFoamError() << "Expected integer or '(', found " << currToken; } } // Set the final offset this->Superclass::LabelListListPtr->SetOffset(listLen, nTotalElems); // Shrink to the actually used size this->Superclass::LabelListListPtr->ResizeData(nTotalElems); io.ReadExpecting(')'); } else { throw vtkFoamError() << "Expected integer, found " << currToken; } } // Read compact labelListList which has offsets and data void vtkFoamEntryValue::ReadCompactLabelListList(vtkFoamIOobject& io) { if (io.IsAsciiFormat()) { this->ReadLabelListList(io); return; } this->SetStreamOption(io); const bool use64BitLabels = io.IsLabel64(); if (use64BitLabels) { this->LabelListListPtr = new vtkFoamLabelListList64; } else { this->LabelListListPtr = new vtkFoamLabelListList32; } this->Superclass::Type = vtkFoamToken::LABELLISTLIST; for (int arrayI = 0; arrayI < 2; arrayI++) { vtkFoamToken currToken; currToken.SetStreamOption(io); if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } if (currToken.IsLabel()) { vtkTypeInt64 listLen = currToken.To(); if (listLen < 0) { throw vtkFoamError() << "Illegal negative list length: " << listLen; } vtkDataArray* array = (arrayI == 0 ? this->Superclass::LabelListListPtr->GetOffsetsArray() : this->Superclass::LabelListListPtr->GetDataArray()); array->SetNumberOfValues(static_cast(listLen)); if (listLen > 0) { // Non-empty (binary) list - only read parentheses only when size > 0 io.ReadExpecting('('); // Begin list io.Read(reinterpret_cast(array->GetVoidPointer(0)), static_cast(listLen * array->GetDataTypeSize())); io.ReadExpecting(')'); // End list } } else { throw vtkFoamError() << "Expected integer, found " << currToken; } } } // Read dimensions set (always ASCII). The leading '[' has already been removed before calling. // - can be integer or floating point // - user-generated files may have only the first five dimensions. // Note // - may even have "human-readable" values such as [kg m^-1 s^-2] but they are very rare // and we silently skip these void vtkFoamEntryValue::ReadDimensionSet(vtkFoamIOobject& io) { const int nDimensions = 7; // There are 7 base dimensions this->MakeScalarList(nDimensions, 0.0); vtkFloatArray& dims = *(this->Superclass::ScalarListPtr); // Read using tokenizer to handle scalar/label, variable lengths, and ignore human-readable vtkFoamToken tok; char expectEnding = ']'; bool goodInput = true; for (int ndims = 0; ndims < nDimensions && goodInput && expectEnding; ++ndims) { if (!io.Read(tok)) { goodInput = false; } else if (tok.IsNumeric()) { dims.SetValue(ndims, tok.ToFloat()); } else if (tok.IsPunctuation()) { if (tok == expectEnding) { expectEnding = '\0'; // Already got the closing ']' } else { goodInput = false; } } else { // Some unknown token type (eg, encountered human-readable units) // - skip until ']' while ((goodInput = io.Read(tok))) { if (tok.IsPunctuation() && (tok == expectEnding)) { expectEnding = '\0'; // Already got the closing ']' break; } } break; } } if (!goodInput) { io.ThrowStackTrace("Unexpected input while parsing dimensions array"); } else if (expectEnding) { io.ReadExpecting(expectEnding); } } //------------------------------------------------------------------------------ // class vtkFoamEntry // a class that represents an entry of a dictionary. note that an // entry can have more than one value. struct vtkFoamEntry : public vtkFoamPtrList { private: typedef vtkFoamPtrList Superclass; std::string Keyword; vtkFoamDict* UpperDictPtr; vtkFoamEntry() = delete; public: vtkFoamEntry(vtkFoamDict* upperDictPtr) : UpperDictPtr(upperDictPtr) { } vtkFoamEntry(const vtkFoamEntry& entry, vtkFoamDict* upperDictPtr) : Superclass(entry.size()) , Keyword(entry.GetKeyword()) , UpperDictPtr(upperDictPtr) { for (size_t i = 0; i < entry.size(); ++i) { (*this)[i] = new vtkFoamEntryValue(*entry[i], this); } } ~vtkFoamEntry() = default; void Clear() { this->Superclass::clear(); } const std::string& GetKeyword() const { return this->Keyword; } void SetKeyword(const std::string& keyword) { this->Keyword = keyword; } const vtkFoamEntryValue& FirstValue() const { return *this->Superclass::operator[](0); } vtkFoamEntryValue& FirstValue() { return *this->Superclass::operator[](0); } const vtkDataArray& LabelList() const { return this->FirstValue().LabelList(); } vtkDataArray& LabelList() { return this->FirstValue().LabelList(); } const vtkFoamLabelListList& LabelListList() const { return this->FirstValue().LabelListList(); } const vtkFloatArray& ScalarList() const { return this->FirstValue().ScalarList(); } vtkFloatArray& ScalarList() { return this->FirstValue().ScalarList(); } const vtkFloatArray& VectorList() const { return this->FirstValue().VectorList(); } vtkFloatArray& VectorList() { return this->FirstValue().VectorList(); } const vtkFoamDict& Dictionary() const { return this->FirstValue().Dictionary(); } vtkFoamDict& Dictionary() { return this->FirstValue().Dictionary(); } const vtkFoamDict* GetUpperDictPtr() const { return this->UpperDictPtr; } template DataType* ReleasePtr() { return this->FirstValue().ReleasePtr(); } std::string ToString() const { return this->empty() ? std::string() : this->FirstValue().ToString(); } float ToFloat() const { return this->empty() ? 0.0F : this->FirstValue().ToFloat(); } double ToDouble() const { return this->empty() ? 0.0 : this->FirstValue().ToDouble(); } vtkTypeInt64 ToInt() const { return this->empty() ? 0 : this->FirstValue().ToInt(); } void ReadDictionary(vtkFoamIOobject& io) { this->Superclass::push_back(new vtkFoamEntryValue(this)); this->Superclass::back()->ReadDictionary(io, vtkFoamToken()); } // read values of an entry void Read(vtkFoamIOobject& io); }; //------------------------------------------------------------------------------ // Code: vtkFoamDict // Copy construct vtkFoamDict::vtkFoamDict(const vtkFoamDict& dict, const vtkFoamDict* upperDictPtr) : Superclass(dict.size()) , UpperDictPtr(upperDictPtr) { if (dict.GetType() == vtkFoamToken::DICTIONARY) { for (size_t i = 0; i < dict.size(); ++i) { (*this)[i] = new vtkFoamEntry(*dict[i], this); } } else { Superclass::assign(dict.size(), nullptr); } } // Destructor vtkFoamDict::~vtkFoamDict() { if (this->Token.GetType() == vtkFoamToken::UNDEFINED) { for (auto* ptr : *this) { delete ptr; } } } // Remove top element, deleting its pointer void vtkFoamDict::remove_back() { if (!Superclass::empty()) { delete Superclass::back(); Superclass::pop_back(); } } // Return list of keywords - table of contents std::vector vtkFoamDict::Toc() const { std::vector list; list.reserve(this->size()); for (const vtkFoamEntry* eptr : *this) { const std::string& key = eptr->GetKeyword(); if (!key.empty()) // should not really happen anyhow { list.push_back(key); } } return list; } // Search dictionary for specified keyword. Return nullptr on failure. vtkFoamEntry* vtkFoamDict::Lookup(const std::string& keyword, const bool isPattern) const { if (this->Token.GetType() == vtkFoamToken::UNDEFINED) { int lastMatch = -1; for (size_t i = 0; i < this->Superclass::size(); i++) { const std::string& key = this->operator[](i)->GetKeyword(); vtksys::RegularExpression rex; if (key == keyword) // found { return this->operator[](i); } else if (isPattern && rex.compile(key) && rex.find(keyword) && rex.start(0) == 0 && rex.end(0) == keyword.size()) { // regular expression matches full keyword lastMatch = static_cast(i); } } if (lastMatch >= 0) { return this->operator[](lastMatch); } } // Not found return nullptr; } // Reads a FoamFile or a subdictionary. // If the stream to be read is a subdictionary, // the preceding '{' is assumed to have already been discarded. bool vtkFoamDict::Read(vtkFoamIOobject& io, const bool isSubDict, const vtkFoamToken& firstKeyword) { try { vtkFoamToken currToken; currToken.SetStreamOption(io); if (firstKeyword.GetType() == vtkFoamToken::UNDEFINED) { // read the first token if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } if (isSubDict) { // the following if clause is for an exceptional expression // of `LABEL{numeric}' without type prefix // (e. g. `2{-0}' in mixedRhoE B.C. in // rhopSonicFoam/shockTube) if (currToken.IsNumeric()) { this->Token = currToken; io.ReadExpecting('}'); return true; } // return as empty dictionary else if (currToken == '}') { return true; } } else { // list of dictionaries is read as a usual dictionary // polyMesh/boundary, point/face/cell-Zones if (currToken.IsLabel()) { io.ReadExpecting('('); if (currToken.To() > 0) { if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } // continue to read as a usual dictionary } else // return as empty dictionary { io.ReadExpecting(')'); return true; } } // some boundary files does not have the number of boundary // patches (e.g. settlingFoam/tank3D). in this case we need to // explicitly read the file as a dictionary. else if (currToken == '(' && io.GetClassName() == "polyBoundaryMesh") // polyMesh/boundary { if (!io.Read(currToken)) // read the first keyword { throw vtkFoamError() << "Unexpected EOF"; } if (currToken == ')') // return as empty dictionary { return true; } } } } // if firstKeyword is given as string read the following stream as // subdictionary else if (firstKeyword.GetType() == vtkFoamToken::STRING) { this->Superclass::push_back(new vtkFoamEntry(this)); this->Superclass::back()->SetKeyword(firstKeyword.ToString()); this->Superclass::back()->ReadDictionary(io); if (!io.Read(currToken) || currToken == '}' || currToken == ')') { return true; } } else // quite likely an identifier { currToken = firstKeyword; } if (currToken == ';' || currToken.IsStringType()) { // general dictionary do { if (currToken.GetType() == vtkFoamToken::STRING) { vtkFoamEntry* previousEntry = this->Lookup(currToken.ToString()); if (previousEntry != nullptr) { if (io.GetInputMode() == vtkFoamFile::INPUT_MODE_MERGE) { if (previousEntry->FirstValue().GetType() == vtkFoamToken::DICTIONARY) { io.ReadExpecting('{'); previousEntry->FirstValue().Dictionary().Read(io, true); // Read as sub-dict } else { previousEntry->Clear(); previousEntry->Read(io); } } else if (io.GetInputMode() == vtkFoamFile::INPUT_MODE_OVERWRITE) { previousEntry->Clear(); previousEntry->Read(io); } else // INPUT_MODE_ERROR { throw vtkFoamError() << "Found duplicated entries with keyword " << currToken.ToString(); } } else { this->Superclass::push_back(new vtkFoamEntry(this)); this->Superclass::back()->SetKeyword(currToken.ToString()); this->Superclass::back()->Read(io); } if (currToken == "FoamFile") { // Drop the FoamFile header subdictionary entry this->remove_back(); } } else if (currToken.GetType() == vtkFoamToken::IDENTIFIER) { // substitute identifier const std::string identifier(currToken.ToIdentifier()); for (const vtkFoamDict* uDictPtr = this;;) { const vtkFoamEntry* identifiedEntry = uDictPtr->Lookup(identifier); if (identifiedEntry != nullptr) { if (identifiedEntry->FirstValue().GetType() != vtkFoamToken::DICTIONARY) { throw vtkFoamError() << "Expected dictionary for substituting entry " << identifier; } const vtkFoamDict& identifiedDict = identifiedEntry->FirstValue().Dictionary(); for (size_t entryI = 0; entryI < identifiedDict.size(); entryI++) { // I think #inputMode handling should be done here // as well, but the genuine FoamFile parser for OF // 1.5 does not seem to be doing it. this->Superclass::push_back(new vtkFoamEntry(*identifiedDict[entryI], this)); } break; } else { uDictPtr = uDictPtr->GetUpperDictPtr(); if (uDictPtr == nullptr) { throw vtkFoamError() << "Substituting entry " << identifier << " not found"; } } } } // skip empty entry only with ';' } while (io.Read(currToken) && (currToken.IsStringType() || currToken == ';')); if (currToken.GetType() == vtkFoamToken::TOKEN_ERROR || currToken == '}' || currToken == ')') { return true; } throw vtkFoamError() << "Expected keyword, closing brace, ';' or EOF, found " << currToken; } throw vtkFoamError() << "Expected keyword or identifier, found " << currToken; } catch (const vtkFoamError& err) { if (isSubDict) { throw; } else { io.SetError(err); return false; } } } //------------------------------------------------------------------------------ // Code: vtkFoamIOobject void vtkFoamIOobject::ReadHeader() { this->Superclass::ReadExpecting("FoamFile"); this->Superclass::ReadExpecting('{'); vtkFoamDict headerDict; headerDict.SetStreamOption(this->GetStreamOption()); headerDict.Read(*this, true); // Read as sub-dict. Throw exception in case of error const vtkFoamEntry* eptr; // Essentials if ((eptr = headerDict.Lookup("class")) == nullptr) { throw vtkFoamError() << "No 'class' in FoamFile header"; } this->headerClassName_ = eptr->ToString(); if ((eptr = headerDict.Lookup("object")) == nullptr) { throw vtkFoamError() << "No 'object' in FoamFile header"; } this->objectName_ = eptr->ToString(); if ((eptr = headerDict.Lookup("format")) == nullptr) { // Note (2021-03-19): may make this optional in the future, defaulting to ascii throw vtkFoamError() << "No 'format' (ascii|binary) in FoamFile header"; } this->SetBinaryFormat("binary" == eptr->ToString()); // case sensitive // The arch entry has "label=(32|64) scalar=(32|64)" // If missing/incomplete, use fallback values from reader (defined in constructor and Close) if ((eptr = headerDict.Lookup("arch")) != nullptr) { const std::string archValue(eptr->ToString()); // Match "label=(32|64)" { auto pos = archValue.find("label="); if (pos != std::string::npos) { pos += 6; // Skip past "label=" if (archValue.compare(pos, 2, "32") == 0) { this->SetLabel64(false); } else if (archValue.compare(pos, 2, "64") == 0) { this->SetLabel64(true); } } } // Match "scalar=(32|64)" { auto pos = archValue.find("scalar="); if (pos != std::string::npos) { pos += 7; // Skip past "scalar=" if (archValue.compare(pos, 2, "32") == 0) { this->SetFloat64(false); } else if (archValue.compare(pos, 2, "64") == 0) { this->SetFloat64(true); } } } } } //------------------------------------------------------------------------------ // Code: vtkFoamEntryValue vtkFoamEntryValue::vtkFoamEntryValue(vtkFoamEntryValue& value, const vtkFoamEntry* parent) : vtkFoamToken(value) , IsUniformEntry(value.IsUniform()) , Managed(true) , UpperEntryPtr(parent) { switch (this->Superclass::Type) { case BOOLLIST: case LABELLIST: case SCALARLIST: case VECTORLIST: case STRINGLIST: { this->Superclass::VtkObjectPtr = value.ToVTKObject(); this->Superclass::VtkObjectPtr->Register(nullptr); break; } case LABELLISTLIST: { if (value.LabelListListPtr->IsLabel64()) { this->LabelListListPtr = new vtkFoamLabelListList64(*value.LabelListListPtr); } else { this->LabelListListPtr = new vtkFoamLabelListList32(*value.LabelListListPtr); } break; } case ENTRYVALUELIST: { const size_t nValues = value.EntryValuePtrs->size(); this->EntryValuePtrs = new vtkFoamPtrList(nValues); for (size_t valueI = 0; valueI < nValues; valueI++) { this->EntryValuePtrs->operator[](valueI) = new vtkFoamEntryValue(*value.EntryValuePtrs->operator[](valueI), this->UpperEntryPtr); } break; } case DICTIONARY: { // UpperEntryPtr is null when called from vtkFoamDict constructor if (this->UpperEntryPtr == nullptr) { this->DictPtr = nullptr; } else { this->DictPtr = new vtkFoamDict(*value.DictPtr, this->UpperEntryPtr->GetUpperDictPtr()); this->DictPtr->SetStreamOption(value.GetStreamOption()); } break; } default: break; } } void vtkFoamEntryValue::Clear() { if (this->Managed) { switch (this->Superclass::Type) { case BOOLLIST: case LABELLIST: case SCALARLIST: case VECTORLIST: case STRINGLIST: this->VtkObjectPtr->Delete(); break; case LABELLISTLIST: delete this->LabelListListPtr; break; case ENTRYVALUELIST: delete this->EntryValuePtrs; break; case DICTIONARY: delete this->DictPtr; break; default: break; } } } // general-purpose list reader - guess the type of the list and read // it. only supports ascii format and assumes the preceding '(' has // already been thrown away. the reader supports nested list with // variable lengths (e. g. `((token token) (token token token)).' // also supports compound of tokens and lists (e. g. `((token token) // token)') only if a list comes as the first value. void vtkFoamEntryValue::ReadList(vtkFoamIOobject& io) { this->SetStreamOption(io); vtkFoamToken currToken; currToken.SetStreamOption(io); io.Read(currToken); if (currToken.IsLabel()) { // Use lookahead to guess the list type. // if the first token is of type LABEL it might be either an element of // a labelList or the size of a sublist so proceed to the next token vtkFoamToken nextToken; nextToken.SetStreamOption(io); if (!io.Read(nextToken)) { throw vtkFoamError() << "Unexpected EOF"; } if (nextToken.IsPunctuation('(')) { // A List of List: read recursively this->Superclass::EntryValuePtrs = new vtkFoamPtrList; this->Superclass::EntryValuePtrs->push_back(new vtkFoamEntryValue(this->UpperEntryPtr)); this->Superclass::EntryValuePtrs->back()->SetStreamOption(*this); this->Superclass::EntryValuePtrs->back()->ReadList(io); this->Superclass::Type = vtkFoamToken::ENTRYVALUELIST; } else if (nextToken.IsPunctuation(')')) { // List with only one label element. Eg "(5)" this->MakeLabelList(1, currToken.To()); return; // DONE } else if (nextToken.IsLabel()) { // Start of a list of labels if (this->IsLabel64()) { auto* array = vtkTypeInt64Array::New(); array->InsertNextValue(currToken.To()); array->InsertNextValue(nextToken.To()); this->Superclass::LabelListPtr = array; } else { auto* array = vtkTypeInt32Array::New(); array->InsertNextValue(currToken.To()); array->InsertNextValue(nextToken.To()); this->Superclass::LabelListPtr = array; } this->Superclass::Type = vtkFoamToken::LABELLIST; } else if (nextToken.IsScalar()) { // Start of a list of scalars this->Superclass::ScalarListPtr = vtkFloatArray::New(); this->Superclass::ScalarListPtr->InsertNextValue(currToken.To()); this->Superclass::ScalarListPtr->InsertNextValue(nextToken.To()); this->Superclass::Type = vtkFoamToken::SCALARLIST; } else { throw vtkFoamError() << "Expected number, '(' or ')', found " << nextToken; } } else if (currToken.IsScalar()) { // The first element of a scalar list this->Superclass::ScalarListPtr = vtkFloatArray::New(); this->Superclass::ScalarListPtr->InsertNextValue(currToken.To()); this->Superclass::Type = vtkFoamToken::SCALARLIST; } else if (currToken.GetType() == vtkFoamToken::STRING) { // if the first word is a string we have to read another token to determine // if it is a keyword for the following dictionary vtkFoamToken nextToken; nextToken.SetStreamOption(io); if (!io.Read(nextToken)) { throw vtkFoamError() << "Unexpected EOF"; } if (nextToken.IsPunctuation('{')) { // Dictionary. Use previously read stringToken as the first keyword if (currToken.ToString().empty()) { throw "Empty string is invalid as a keyword for dictionary entry"; } this->ReadDictionary(io, currToken); // the dictionary read as list has the entry terminator ';' so // we have to skip it return; // DONE } else if (nextToken.IsPunctuation(')')) { // List with only one string element. Eg "(wall)" this->Superclass::StringListPtr = vtkStringArray::New(); this->Superclass::StringListPtr->SetNumberOfValues(1); this->Superclass::StringListPtr->SetValue(0, currToken.ToString()); this->Superclass::Type = vtkFoamToken::STRINGLIST; return; // DONE } else if (nextToken.GetType() == vtkFoamToken::STRING) // list of strings { this->Superclass::StringListPtr = vtkStringArray::New(); this->Superclass::StringListPtr->InsertNextValue(currToken.ToString()); this->Superclass::StringListPtr->InsertNextValue(nextToken.ToString()); this->Superclass::Type = vtkFoamToken::STRINGLIST; } else { throw vtkFoamError() << "Expected string, '{' or ')', found " << nextToken; } } else if (currToken == '(' || currToken == '{') { // List of lists or dictionaries: read recursively this->Superclass::EntryValuePtrs = new vtkFoamPtrList; this->Superclass::EntryValuePtrs->push_back(new vtkFoamEntryValue(this->UpperEntryPtr)); this->Superclass::EntryValuePtrs->back()->SetStreamOption(io); if (currToken == '(') { this->Superclass::EntryValuePtrs->back()->ReadList(io); } else // currToken == '{' { this->Superclass::EntryValuePtrs->back()->ReadDictionary(io, vtkFoamToken()); } // read all the following values as arbitrary entryValues // the alphaContactAngle b.c. in multiphaseInterFoam/damBreak4phase // reaquires this treatment (reading by readList() is not enough) do { this->Superclass::EntryValuePtrs->push_back(new vtkFoamEntryValue(this->UpperEntryPtr)); this->Superclass::EntryValuePtrs->back()->SetStreamOption(io); this->Superclass::EntryValuePtrs->back()->Read(io); } while (*this->Superclass::EntryValuePtrs->back() != ')' && *this->Superclass::EntryValuePtrs->back() != '}' && *this->Superclass::EntryValuePtrs->back() != ';'); if (*this->Superclass::EntryValuePtrs->back() != ')') { throw vtkFoamError() << "Expected ')' before " << *this->Superclass::EntryValuePtrs->back(); } // Drop ')' this->Superclass::EntryValuePtrs->remove_back(); this->Superclass::Type = vtkFoamToken::ENTRYVALUELIST; return; // DONE } else if (currToken.IsPunctuation(')')) { // Empty list this->Superclass::Type = vtkFoamToken::EMPTYLIST; return; // DONE } // FIXME: may (or may not) need identifier handling while (io.Read(currToken) && !currToken.IsPunctuation(')')) { if (this->Superclass::Type == vtkFoamToken::LABELLIST) { if (currToken.GetType() == vtkFoamToken::SCALAR) { // Encountered a scalar while reading a labelList - switch representation // Need intermediate pointer since LabelListPtr and ScalarListPtr are in a union vtkDataArray* labels = this->LabelListPtr; const vtkIdType currLen = labels->GetNumberOfTuples(); const bool use64BitLabels = ::Is64BitArray(labels); // <- Same as io.IsLabel64() // Copy, with append auto* scalars = vtkFloatArray::New(); scalars->SetNumberOfValues(currLen + 1); for (vtkIdType i = 0; i < currLen; ++i) { scalars->SetValue(i, static_cast(GetLabelValue(labels, i, use64BitLabels))); } scalars->SetValue(currLen, currToken.To()); // Append value // Replace labels->Delete(); this->Superclass::ScalarListPtr = scalars; this->Superclass::Type = vtkFoamToken::SCALARLIST; } else if (currToken.IsLabel()) { if (currToken.IsLabel64()) { assert(vtkTypeInt64Array::FastDownCast(this->LabelListPtr) != nullptr); static_cast(this->LabelListPtr) ->InsertNextValue(currToken.To()); } else { assert(vtkTypeInt32Array::FastDownCast(this->LabelListPtr) != nullptr); static_cast(this->LabelListPtr) ->InsertNextValue(currToken.To()); } } else { throw vtkFoamError() << "Expected a number, found " << currToken; } } else if (this->Superclass::Type == vtkFoamToken::SCALARLIST) { if (currToken.IsNumeric()) { this->Superclass::ScalarListPtr->InsertNextValue(currToken.To()); } else if (currToken == '(') { vtkDebugWithObjectMacro(nullptr, "Found a list containing scalar data followed " "by a nested list, but this reader only " "supports nested lists that precede all " "scalars. Discarding nested list data."); vtkFoamEntryValue tmp(this->UpperEntryPtr); tmp.SetStreamOption(io); tmp.ReadList(io); } else { throw vtkFoamError() << "Expected a number, found " << currToken; } } else if (this->Superclass::Type == vtkFoamToken::STRINGLIST) { if (currToken.GetType() == vtkFoamToken::STRING) { this->Superclass::StringListPtr->InsertNextValue(currToken.ToString()); } else { throw vtkFoamError() << "Expected a string, found " << currToken; } } else if (this->Superclass::Type == vtkFoamToken::ENTRYVALUELIST) { if (currToken.IsLabel()) { // skip the number of elements to make things simple if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } } if (currToken != '(') { throw vtkFoamError() << "Expected '(', found " << currToken; } this->Superclass::EntryValuePtrs->push_back(new vtkFoamEntryValue(this->UpperEntryPtr)); this->Superclass::EntryValuePtrs->back()->ReadList(io); } else { throw vtkFoamError() << "Unexpected token " << currToken; } } if (this->Superclass::Type == vtkFoamToken::BOOLLIST) { this->Superclass::BoolListPtr->Squeeze(); } else if (this->Superclass::Type == vtkFoamToken::LABELLIST) { this->Superclass::LabelListPtr->Squeeze(); } else if (this->Superclass::Type == vtkFoamToken::SCALARLIST) { this->Superclass::ScalarListPtr->Squeeze(); } else if (this->Superclass::Type == vtkFoamToken::STRINGLIST) { this->Superclass::StringListPtr->Squeeze(); } } // a list of dictionaries is actually read as a dictionary void vtkFoamEntryValue::ReadDictionary(vtkFoamIOobject& io, const vtkFoamToken& firstKeyword) { this->Superclass::DictPtr = new vtkFoamDict(this->UpperEntryPtr->GetUpperDictPtr()); this->DictPtr->SetStreamOption(io); this->Superclass::Type = vtkFoamToken::DICTIONARY; this->Superclass::DictPtr->Read(io, true, firstKeyword); } // Guess the type of the given entry value and read it // Return: // - true on success // - false if end of entry (';') encountered during parsing composite entry value bool vtkFoamEntryValue::Read(vtkFoamIOobject& io) { this->SetStreamOption(io); vtkFoamToken currToken; currToken.SetStreamOption(io); if (!io.Read(currToken)) { throw vtkFoamError() << "Unexpected EOF"; } // List types if (currToken == '{') { this->ReadDictionary(io, vtkFoamToken()); return true; } else if (currToken == '(') { this->ReadList(io); return true; } else if (currToken == '[') { this->ReadDimensionSet(io); return true; } vtkFoamTypes::dataType listDataType(vtkFoamTypes::NO_TYPE); if (currToken == "uniform") { if (!io.Read(currToken)) { throw vtkFoamError() << "Expected a uniform value or a list, found unexpected EOF"; } if (currToken == '(') { this->ReadList(io); } else if (currToken == ';') { this->Superclass::operator=("uniform"); return false; } else if (currToken.IsNumeric() || currToken.GetType() == vtkFoamToken::STRING) { this->Superclass::operator=(currToken); } else // unexpected punctuation token { throw vtkFoamError() << "Expected number, string or ( for uniform entry, found " << currToken; } this->IsUniformEntry = true; } else if (currToken == "nonuniform") { if (!io.Read(currToken)) { throw vtkFoamError() << "Expected list type specifier for nonuniform entry, found EOF"; } this->IsUniformEntry = false; if (currToken.GetType() == vtkFoamToken::STRING) { // List types: "List