Commit 104bcfaa authored by Joachim Pouderoux's avatar Joachim Pouderoux Committed by Kitware Robot
Browse files

Merge topic 'Xdmf2ReaderWriterOptimizations'

90bd6dc1 Add test and fix reader
a5e481dc fix header test
7bb9fd44 Fix error with java wrapping
5d3c7d68 code cleaning
54f52ce1

 Enhance and optimize Xdmf2 reader and writer to support static meshes.
Acked-by: Kitware Robot's avatarKitware Robot <kwrobot@kitware.com>
Reviewed-by: default avatarDavid E DeMarle <dave.demarle@kitware.com>
Merge-request: !347
parents 87dc9cac 90bd6dc1
ExternalData_Expand_Arguments(VTKData _
"DATA{${VTK_TEST_INPUT_DIR}/XDMF/,REGEX:.*}"
)
vtk_add_test_cxx(${vtk-module}CxxTests tests
NO_DATA NO_VALID
XdmfTestVTKIO.cxx
TestTemporalXdmfReaderWriter.cxx,NO_VALID
XdmfTestVTKIO.cxx,NO_VALID
)
vtk_test_cxx_executable(${vtk-module}CxxTests tests)
/*=========================================================================
Program: Visualization Toolkit
Module: TestTemporalXdmfReaderWriter.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.
=========================================================================*/
// Description:
// This tests temporal reading and writing of static meshes using
// vtkXdmfReader and vtkXdmfWriter.
#include "vtkInformation.h"
#include "vtkMultiBlockDataSet.h"
#include "vtkNew.h"
#include "vtkPointData.h"
#include "vtkStreamingDemandDrivenPipeline.h"
#include "vtkTestUtilities.h"
#include "vtkThreshold.h"
#include "vtkUnstructuredGrid.h"
#include "vtkXdmfReader.h"
#include "vtkXdmfWriter.h"
#include "vtksys/SystemTools.hxx"
#define ASSERT_TEST(_cond_, _msg_) \
if (!(_cond_)) { std::cerr << _msg_ << std::endl; return VTK_ERROR; }
int TestStaticMesh(vtkXdmfReader* reader)
{
reader->UpdateInformation();
vtkInformation* outInfo = reader->GetExecutive()->GetOutputInformation(0);
int steps = (outInfo->Has(vtkStreamingDemandDrivenPipeline::TIME_STEPS())) ?
outInfo->Length(vtkStreamingDemandDrivenPipeline::TIME_STEPS()) : 0;
ASSERT_TEST(steps == 3, "Read data does not have 3 time steps as expected!");
double* timeSteps =
outInfo->Get(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
vtkPoints* geometryAtT0 = 0;
vtkCellArray* topologyAtT0 = 0;
for (int i = 0; i < steps; i++)
{
double updateTime = timeSteps[i];
outInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP(),
updateTime);
reader->Update();
vtkMultiBlockDataSet* mb =
vtkMultiBlockDataSet::SafeDownCast(reader->GetOutputDataObject(0));
ASSERT_TEST(mb, "Root data is not a multiblock data set as expected!");
ASSERT_TEST(mb->GetNumberOfBlocks() == 2, "Root multiblock data is supposed to have 2 blocks!");
vtkUnstructuredGrid* grid =
vtkUnstructuredGrid::SafeDownCast(mb->GetBlock(0));
ASSERT_TEST(grid, "Block 0 is not an unstructured grid as expected!");
if (i == 0)
{
geometryAtT0 = grid->GetPoints();
topologyAtT0 = grid->GetCells();
}
ASSERT_TEST(grid->GetPoints() == geometryAtT0, "Geometry is not static over time as expected!");
ASSERT_TEST(grid->GetCells() == topologyAtT0, "Topology is not static over time as expected!");
}
return 0;
}
int TestTemporalXdmfReaderWriter(int argc, char *argv[])
{
// Read the input data file
char *filePath = vtkTestUtilities::ExpandDataFileName(argc, argv,
"Data/XDMF/temporalStaticMeshes.xmf");
vtkNew<vtkXdmfReader> reader;
reader->SetFileName(filePath);
if (TestStaticMesh(reader.Get()) == VTK_ERROR)
{
std::cerr << "Error while reading " << reader->GetFileName() << std::endl;
return VTK_ERROR;
}
// Write the input data to a new Xdmf file
std::string outFilePath = "temporalStaticMeshesTest.xmf";
vtkNew<vtkXdmfWriter> writer;
writer->SetFileName(outFilePath.c_str());
writer->WriteAllTimeStepsOn();
writer->MeshStaticOverTimeOn();
writer->SetInputConnection(reader->GetOutputPort());
writer->Write();
// Test written file
vtkNew<vtkXdmfReader> reader2;
reader2->SetFileName(outFilePath.c_str());
if (TestStaticMesh(reader2.Get()) == VTK_ERROR)
{
std::cerr << "Error while reading " << reader2->GetFileName() << std::endl;
return VTK_ERROR;
}
return 0;
}
......@@ -43,6 +43,8 @@
#include <deque>
#include <cassert>
#include <libxml/tree.h>
using namespace xdmf2;
static void vtkScaleExtents(int in_exts[6], int out_exts[6], int stride[3])
......@@ -137,7 +139,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadData()
}
//----------------------------------------------------------------------------
vtkDataObject* vtkXdmfHeavyData::ReadData(XdmfGrid* xmfGrid)
vtkDataObject* vtkXdmfHeavyData::ReadData(XdmfGrid* xmfGrid, int blockId)
{
if (!xmfGrid || xmfGrid->GetGridType() == XDMF_GRID_UNSET)
{
......@@ -151,7 +153,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadData(XdmfGrid* xmfGrid)
{
// grid is a temporal collection, pick the sub-grid with matching time and
// process that.
return this->ReadTemporalCollection(xmfGrid);
return this->ReadTemporalCollection(xmfGrid, blockId);
}
else if (gridType == XDMF_GRID_COLLECTION ||
gridType == XDMF_GRID_TREE)
......@@ -160,7 +162,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadData(XdmfGrid* xmfGrid)
}
// grid is a primitive grid, so read the data.
return this->ReadUniformData(xmfGrid);
return this->ReadUniformData(xmfGrid, blockId);
}
//----------------------------------------------------------------------------
......@@ -189,7 +191,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadComposite(XdmfGrid* xmfComposite)
if (!child_is_leaf || !distribute_leaf_nodes ||
(number_of_leaf_nodes % this->NumberOfPieces) == this->Piece)
{
vtkDataObject* childDO = this->ReadData(xmfChild);
vtkDataObject* childDO = this->ReadData(xmfChild, cc);
if (childDO)
{
multiBlock->SetBlock(cc, childDO);
......@@ -204,7 +206,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadComposite(XdmfGrid* xmfComposite)
//----------------------------------------------------------------------------
vtkDataObject* vtkXdmfHeavyData::ReadTemporalCollection(
XdmfGrid* xmfTemporalCollection)
XdmfGrid* xmfTemporalCollection, int blockId)
{
assert(xmfTemporalCollection->GetGridType() & XDMF_GRID_COLLECTION &&
xmfTemporalCollection->GetCollectionType() == XDMF_GRID_COLLECTION_TEMPORAL
......@@ -254,7 +256,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadTemporalCollection(
std::deque<XdmfGrid*>::iterator iter;
for (iter = valid_children.begin(); iter != valid_children.end(); ++iter)
{
vtkDataObject* childDO = this->ReadData(*iter);
vtkDataObject* childDO = this->ReadData(*iter, blockId);
if (childDO)
{
child_data_objects.push_back(childDO);
......@@ -286,7 +288,7 @@ vtkDataObject* vtkXdmfHeavyData::ReadTemporalCollection(
//----------------------------------------------------------------------------
// Read a non-composite grid. Note here uniform has nothing to do with
// vtkUniformGrid but to what Xdmf's GridType="Uniform".
vtkDataObject* vtkXdmfHeavyData::ReadUniformData(XdmfGrid* xmfGrid)
vtkDataObject* vtkXdmfHeavyData::ReadUniformData(XdmfGrid* xmfGrid, int blockId)
{
assert(xmfGrid->IsUniform() && "Input must be a uniform xdmf grid.");
......@@ -300,6 +302,79 @@ vtkDataObject* vtkXdmfHeavyData::ReadUniformData(XdmfGrid* xmfGrid)
// Read heavy data for grid geometry/topology. This does not read any
// data-arrays. They are read explicitly.
XdmfTopology* topo = xmfGrid->GetTopology();
XdmfGeometry* geom = xmfGrid->GetGeometry();
xmlChar* filePtr;
bool caching = true;
XdmfDOM* topoDom = topo->GetDOM();
XdmfXmlNode topoNode = topo->GetElement();
XdmfXmlNode topoNodeDataItem = topoDom->FindElement("DataItem", 0, topoNode);
std::string topoFilename = "NULL";
if (topoNodeDataItem && caching)
{
filePtr = topoNodeDataItem->children->content;
if (filePtr != NULL)
{
topoFilename = reinterpret_cast<char*>(filePtr);
}
else
{
vtkErrorWithObjectMacro(this->Reader, << "Cannot find DataItem element in topology xml, no caching possible");
caching = false;
}
}
else
{
caching = false;
}
XdmfDOM* geomDom = geom->GetDOM();
XdmfXmlNode geomNode = geom->GetElement();
XdmfXmlNode geomNodeDataItem = geomDom->FindElement("DataItem", 0, geomNode);
std::string geomFilename = "NULL";
if (geomNodeDataItem && caching)
{
filePtr = geomNodeDataItem->children->content;
if (filePtr != NULL)
{
geomFilename = reinterpret_cast<char*>(filePtr);
}
else
{
vtkErrorWithObjectMacro(this->Reader, << "Cannot find DataItem element in geometry xml, no caching possible");
caching = false;
}
}
else
{
caching = false;
}
vtkXdmfReader::XdmfReaderCachedData& cache =
vtkXdmfReader::SafeDownCast(this->Reader)->GetDataSetCache();
vtkXdmfReader::XdmfDataSetTopoGeoPath& cachedData = cache[blockId];
if (caching &&
(cachedData.topologyPath == topoFilename) && (cachedData.geometryPath == geomFilename))
{
vtkDataSet* ds = vtkDataSet::SafeDownCast(
vtkDataObjectTypes::NewDataObject(cachedData.dataset->GetDataObjectType()));
ds->ShallowCopy(cachedData.dataset);
this->ReadAttributes(ds, xmfGrid);
return ds;
}
if (caching)
{
cachedData.topologyPath = topoFilename;
cachedData.geometryPath = geomFilename;
if (cache[blockId].dataset != NULL)
{
cache[blockId].dataset->Delete();
cache[blockId].dataset = NULL;
}
}
XdmfInt32 status = xmfGrid->Update();
if (status == XDMF_FAIL)
{
......@@ -335,6 +410,11 @@ vtkDataObject* vtkXdmfHeavyData::ReadUniformData(XdmfGrid* xmfGrid)
return 0;
}
if (caching)
{
cache[blockId].dataset = vtkDataSet::SafeDownCast(dataObject);
dataObject->Register(0);
}
return dataObject;
}
......
......@@ -55,7 +55,7 @@ public:
~vtkXdmfHeavyData();
// Description:
vtkDataObject* ReadData(xdmf2::XdmfGrid* xmfGrid);
vtkDataObject* ReadData(xdmf2::XdmfGrid* xmfGrid, int blockId = -1);
// Description:
vtkDataObject* ReadData();
......@@ -70,10 +70,11 @@ public:
// of points possible.
static int GetNumberOfPointsPerCell(int vtk_cell_type);
private:
// Description:
// Read a temporal collection.
vtkDataObject* ReadTemporalCollection(xdmf2::XdmfGrid* xmfTemporalCollection);
vtkDataObject* ReadTemporalCollection(xdmf2::XdmfGrid* xmfTemporalCollection, int blockId);
// Description:
// Read a spatial-collection or a tree.
......@@ -82,7 +83,7 @@ private:
// Description:
// Read a non-composite grid. Note here uniform has nothing to do with
// vtkUniformGrid but to what Xdmf's GridType="Uniform".
vtkDataObject* ReadUniformData(xdmf2::XdmfGrid* xmfGrid);
vtkDataObject* ReadUniformData(xdmf2::XdmfGrid* xmfGrid, int blockId);
// Description:
// Reads the topology and geometry for an unstructured grid. Does not read any
......@@ -156,7 +157,6 @@ private:
// Creates a new dataset with egdes selected by the set, extracting them from
// the input dataset.
vtkDataSet* ExtractEdges(xdmf2::XdmfSet* xmfSet, vtkDataSet* dataSet);
};
#endif
......@@ -19,6 +19,7 @@
#include "vtkCharArray.h"
#include "vtkCompositeDataPipeline.h"
#include "vtkDataObjectTypes.h"
#include "vtkDataSet.h"
#include "vtkExtentTranslator.h"
#include "vtkInformation.h"
#include "vtkInformationVector.h"
......@@ -116,6 +117,8 @@ vtkXdmfReader::~vtkXdmfReader()
delete this->CellArraysCache;
delete this->GridsCache;
delete this->SetsCache;
this->ClearDataSetCache();
}
//----------------------------------------------------------------------------
......@@ -319,8 +322,13 @@ int vtkXdmfReader::RequestInformation(vtkInformation *, vtkInformationVector **,
outInfo->Set(vtkDataObject::SIL(), domain->GetSIL());
// * Publish time information.
std::vector<double> time_steps(domain->GetTimeSteps().begin(),
domain->GetTimeSteps().end());
const std::map<int, XdmfFloat64>& ts = domain->GetTimeStepsRev();
std::vector<double> time_steps(ts.size());
std::map<int, XdmfFloat64>::const_iterator it = ts.begin();
for (int i = 0; it != ts.end(); i++, ++it)
{
time_steps[i] = it->second;
}
if (time_steps.size() > 0)
{
......@@ -386,6 +394,10 @@ int vtkXdmfReader::RequestData(vtkInformation *, vtkInformationVector **,
}
this->LastTimeIndex = this->ChooseTimeStep(outInfo);
if (this->LastTimeIndex == 0)
{
this->ClearDataSetCache();
}
vtkXdmfHeavyData dataReader(this->XdmfDocument->GetActiveDomain(), this);
dataReader.Piece = updatePiece;
......@@ -615,3 +627,24 @@ vtkGraph* vtkXdmfReader::GetSIL()
}
return 0;
}
//----------------------------------------------------------------------------
void vtkXdmfReader::ClearDataSetCache()
{
XdmfReaderCachedData::iterator it = this->DataSetCache.begin();
while (it != this->DataSetCache.end())
{
if (it->second.dataset != NULL)
{
it->second.dataset->Delete();
}
++it;
}
this->DataSetCache.clear();
}
//----------------------------------------------------------------------------
vtkXdmfReader::XdmfReaderCachedData& vtkXdmfReader::GetDataSetCache()
{
return this->DataSetCache;
}
......@@ -35,6 +35,7 @@
#include "vtkIOXdmf2Module.h" // For export macro
#include "vtkDataReader.h"
#include <map> // for caching
class vtkXdmfArraySelection;
class vtkXdmfDocument;
......@@ -140,6 +141,21 @@ public:
virtual vtkGraph* GetSIL();
//BTX
class XdmfDataSetTopoGeoPath
{
public:
XdmfDataSetTopoGeoPath() : dataset(0), topologyPath(), geometryPath() {}
vtkDataSet* dataset;
std::string topologyPath;
std::string geometryPath;
};
typedef std::map<int, XdmfDataSetTopoGeoPath> XdmfReaderCachedData;
// Description
// Get the data set cache
XdmfReaderCachedData& GetDataSetCache();
protected:
vtkXdmfReader();
~vtkXdmfReader();
......@@ -182,11 +198,16 @@ protected:
vtkXdmfArraySelection* SetsCache;
int SILUpdateStamp;
XdmfReaderCachedData DataSetCache;
private:
// Description:
// Prepares the XdmfDocument.
bool PrepareDocument();
void ClearDataSetCache();
// Description:
// Returns the time-step index requested using the UPDATE_TIME_STEPS from the
// information.
......
......@@ -218,7 +218,7 @@ vtkXdmfDomain::vtkXdmfDomain(XdmfDOM* xmlDom, int domain_index)
this->NumberOfGrids = this->XMLDOM->FindNumberOfElements("Grid", this->XMLDomain);
this->XMFGrids = new XdmfGrid[this->NumberOfGrids+1];
XdmfXmlNode xmlGrid = this->XMLDOM->FindElement("Grid", 0, this->XMLDomain);
XdmfXmlNode xmlGrid = this->XMLDOM->FindElement("Grid", 0, this->XMLDomain);
XdmfInt64 cc=0;
while (xmlGrid)
{
......@@ -297,14 +297,14 @@ int vtkXdmfDomain::GetVTKDataType(XdmfGrid* xmfGrid)
{
return VTK_MULTIBLOCK_DATA_SET;
}
if (xmfGrid->GetTopology()->GetClass() == XDMF_UNSTRUCTURED )
if (xmfGrid->GetTopology()->GetClass() == XDMF_UNSTRUCTURED )
{
return VTK_UNSTRUCTURED_GRID;
}
}
XdmfInt32 topologyType = xmfGrid->GetTopology()->GetTopologyType();
if (topologyType == XDMF_2DSMESH || topologyType == XDMF_3DSMESH )
{
return VTK_STRUCTURED_GRID;
{
return VTK_STRUCTURED_GRID;
}
else if (topologyType == XDMF_2DCORECTMESH ||
topologyType == XDMF_3DCORECTMESH)
......@@ -325,7 +325,13 @@ int vtkXdmfDomain::GetVTKDataType(XdmfGrid* xmfGrid)
//----------------------------------------------------------------------------
int vtkXdmfDomain::GetIndexForTime(double time)
{
std::set<XdmfFloat64>::iterator iter = this->TimeSteps.upper_bound(time);
std::map<XdmfFloat64, int>::const_iterator iter = this->TimeSteps.find(time);
if (iter != this->TimeSteps.end())
{
return iter->second;
}
iter = this->TimeSteps.upper_bound(time);
if (iter == this->TimeSteps.begin())
{
// The requested time step is before any available time. We will use it by
......@@ -337,7 +343,7 @@ int vtkXdmfDomain::GetIndexForTime(double time)
iter--;
}
std::set<XdmfFloat64>::iterator iter2 = this->TimeSteps.begin();
std::map<XdmfFloat64, int>::iterator iter2 = this->TimeSteps.begin();
int counter = 0;
while (iter2 != iter)
{
......@@ -437,7 +443,7 @@ bool vtkXdmfDomain::GetOriginAndSpacing(XdmfGrid* xmfGrid,
XdmfGeometry *xmfGeometry = xmfGrid->GetGeometry();
if (xmfGeometry->GetGeometryType() == XDMF_GEOMETRY_ORIGIN_DXDYDZ )
{
{
// Update geometry so that origin and spacing are read
xmfGeometry->Update(); // read heavy-data for the geometry.
XdmfFloat64 *xmfOrigin = xmfGeometry->GetOrigin();
......@@ -542,7 +548,7 @@ void vtkXdmfDomain::CollectMetaData()
this->Grids->clear();
// We have aborted collecting grids information since it was too numerous to
// be of any use to the user.
// be of any use to the user.
this->SILBuilder->Initialize();
blocksRoot = this->SILBuilder->AddVertex("Blocks");
hierarchyRoot = this->SILBuilder->AddVertex("Hierarchy");
......@@ -623,7 +629,7 @@ void vtkXdmfDomain::CollectNonLeafMetaData(XdmfGrid* xmfGrid,
XdmfInt32 numChildren = xmfGrid->GetNumberOfChildren();
for (XdmfInt32 cc=0; cc < numChildren; cc++)
{
XdmfGrid* xmfChild = xmfGrid->GetChild(cc);
XdmfGrid* xmfChild = xmfGrid->GetChild(cc);
this->CollectMetaData(xmfChild, silVertex);
}
......@@ -641,7 +647,12 @@ void vtkXdmfDomain::CollectNonLeafMetaData(XdmfGrid* xmfGrid,
XdmfTime* xmfTime = xmfGrid->GetTime();
if (xmfTime && xmfTime->GetTimeType() != XDMF_TIME_UNSET)
{
this->TimeSteps.insert(xmfTime->GetValue());
int step = static_cast<int>(this->TimeSteps.size());
if (this->TimeSteps.find(xmfTime->GetValue()) == this->TimeSteps.end())
{
this->TimeSteps[xmfTime->GetValue()] = step;
this->TimeStepsRev[step] = xmfTime->GetValue();
}
}
}
}
......@@ -719,7 +730,7 @@ void vtkXdmfDomain::CollectLeafMetaData(XdmfGrid* xmfGrid, vtkIdType silParent)
// XdmfInt32 setCenter = xmfSet->GetSetType();
// Not sure if we want to create separate lists for different types of sets
// or just treat all the sets as same. For now, we are treating them as
// or just treat all the sets as same. For now, we are treating them as
// the same.
this->Sets->AddArray(name);
}
......@@ -728,7 +739,12 @@ void vtkXdmfDomain::CollectLeafMetaData(XdmfGrid* xmfGrid, vtkIdType silParent)
XdmfTime* xmfTime = xmfGrid->GetTime();
if (xmfTime && xmfTime->GetTimeType() != XDMF_TIME_UNSET)
{
this->TimeSteps.insert(xmfTime->GetValue());
int step = static_cast<int>(this->TimeSteps.size());
if (this->TimeSteps.find(xmfTime->GetValue()) == this->TimeSteps.end())
{
this->TimeSteps[xmfTime->GetValue()] = step;
this->TimeStepsRev[step] = xmfTime->GetValue();
}
}
}
......
......@@ -203,8 +203,9 @@ private:
vtkXdmfArraySelection* CellArrays;
vtkXdmfArraySelection* Grids;
vtkXdmfArraySelection* Sets;
std::set<XdmfFloat64> TimeSteps; //< Only discrete timesteps are currently
std::map<XdmfFloat64, int> TimeSteps; //< Only discrete timesteps are currently
// supported.
std::map<int, XdmfFloat64> TimeStepsRev;
public:
//---------------------------------------------------------------------------
......@@ -244,8 +245,10 @@ public:
//---------------------------------------------------------------------------
// Description:
// Returns the timesteps.
const std::set<XdmfFloat64>& GetTimeSteps()
const std::map<XdmfFloat64, int>& GetTimeSteps()
{ return this->TimeSteps; }
const std::map<int,XdmfFloat64>& GetTimeStepsRev()
{ return this->TimeStepsRev; }
//---------------------------------------------------------------------------
// Description:
......@@ -257,18 +260,8 @@ public:
// Returns the time value at the given index.
XdmfFloat64 GetTimeForIndex(int index)
{
std::set<XdmfFloat64>::iterator iter;
int cc=0;
for (iter = this->TimeSteps.begin(); iter != this->TimeSteps.end();
iter++, cc++)
{
if (cc == index)
{
return *iter;
}
}
// invalid index.
return 0.0;
std::map<int, XdmfFloat64>::iterator iter = this->TimeStepsRev.find(index);
return (iter != this->TimeStepsRev.end()) ? iter->second : 0.0;
}
//---------------------------------------------------------------------------
......
......@@ -41,6 +41,7 @@
#include "vtkStructuredGrid.h"
#include "vtkTypeTraits.h"
#include "vtkUnstructuredGrid.h"
#include "vtksys/SystemTools.hxx"
#include "XdmfArray.h"
#include "XdmfAttribute.h"
......@@ -56,6 +57,7 @@
#include <algorithm>
#include <map>
#include <stdio.h>
#include <sstream>
#include <vector>
#include <libxml/tree.h> // always after std::blah stuff
......@@ -186,6 +188,7 @@ vtkXdmfWriter::vtkXdmfWriter()
this->TopTemporalGrid = NULL;
this->DomainMemoryHandler = NULL;
this->SetNumberOfOutputPorts(0);
this->MeshStaticOverTime = false;
}
//----------------------------------------------------------------------------
......@@ -248,6 +251,10 @@ int vtkXdmfWriter::Write()
// always write even if the data hasn't changed
this->Modified();
this->TopologyAtT0.clear();
this->GeometryAtT0.clear();
this->UnlabelledDataArrayId = 0;