Python-Created HyperTreeGrid File Performance on Paraview

Hi,

I am converting TreeBased-AMR CFD simulation data into VTK HyperTreeGrid using a python script similar the one of this tutorial (simplified script version at the end).

The created htg file from relatively large data snapshot (>2GB) is handled well in Paraview when the AMR is refined in one region (see the slice in fig. 1 below).

However, when the refined regions are more spread, even for small initial data snapshot (200MB, see fig. 2 below), Paraview struggles a lot to handle the file (opening/ rendering/ or any action).

I am just starting to use htg and VTK for visualisation. Is there a way to improve the performance in this case? I am a bit surprised that in the second case Paraview struggles with only 200MB data.

Best wishes,
Nicolas

Here the simplified script:

import numpy as np
import vtk
from vtkmodules.vtkCommonCore import (
    vtkDoubleArray,
    vtkLookupTable,)
from vtkmodules.vtkCommonDataModel import (
    vtkHyperTreeGrid,
    vtkHyperTreeGridNonOrientedCursor,)



def do_refine_root(xf,i,j,k,level,levelmax=3,data=[]):
    if (i==j)&(i==k)&(level<levelmax):
        return True
    else:
        return False

def refine_cell(cursor,x,y,z,max_level,data=[]):
    ii=0
    cursor.SubdivideLeaf()
    ii1 = 0
    for i1 in range(2):
        for j1 in range(2):
            for k1 in range(2):       
                cursor.ToChild(ii1)
                iLevel = cursor.GetLevel()
                iBounds = cursor.GetGrid()
                print(iLevel)
                print(iBounds)    
                idx = cursor.GetGlobalNodeIndex()
                print(iLevel)
                refine = do_refine_root(0,i1,j1,k1,iLevel,levelmax=max_level)
                if refine:
                    scalarArray.InsertTuple1(idx, ii+ii1)
                    vecArray.InsertComponent(idx, 0, ii+ii1)
                    vecArray.InsertComponent(idx, 1, ii+ii1)
                    vecArray.InsertComponent(idx, 2, ii+ii1)
                    refine_cell(cursor,x,y,z,max_level, data=[])
                else:
                    scalarArray.InsertTuple1(idx, ii+ii1)
                    vecArray.InsertComponent(idx, 0, ii+ii1)
                    vecArray.InsertComponent(idx, 1, ii+ii1)
                    vecArray.InsertComponent(idx, 2, ii+ii1)
                cursor.ToParent()  
                ii1 += 1
    return ii


# =========================================== #
# ===          IC PARAMETERS              === #
# =========================================== #
N_base = 4
Lbox = 10.
dim=3
max_level=3
# =========================================== #
x_c = np.linspace(0.5*Lbox/N_base,Lbox-0.5*Lbox/N_base,N_base)
x_f = np.linspace(0.,Lbox,N_base+1)
Ntot_root = N_base**dim

# =========================================== #
# ===        INITIALIZE THE GRID          === #
# =========================================== #
htg = vtk.vtkHyperTreeGrid()
htg.Initialize()

scalarArray = vtkDoubleArray()
scalarArray.SetName('scalar')
scalarArray.SetNumberOfValues(0)
htg.GetCellData().AddArray(scalarArray)
htg.GetCellData().SetActiveScalars('scalar')

vecArray = vtkDoubleArray()
vecArray.SetName('vector')
vecArray.SetNumberOfValues(0)
vecArray.SetNumberOfComponents(3)
htg.GetCellData().AddArray(vecArray)
htg.GetCellData().SetActiveVectors('vector')

if (dim==2):
    htg.SetDimensions(N_base+1, N_base+1, 1)
elif(dim==3):
    htg.SetDimensions(N_base+1, N_base+1, N_base+1)
htg.SetBranchFactor(2)

# Rectilinear root grid coordinates
xValues = vtkDoubleArray()
xValues.SetNumberOfValues(N_base+1)
yValues = vtkDoubleArray()
yValues.SetNumberOfValues(N_base+1)
for i in range(N_base+1):
    xValues.SetValue(i,x_f[i])
    yValues.SetValue(i,x_f[i])
htg.SetXCoordinates(xValues)
htg.SetYCoordinates(yValues)
if dim==3:
    zValues = vtkDoubleArray()
    zValues.SetNumberOfValues(N_base+1)
    for i in range(N_base+1):
        zValues.SetValue(i,x_f[i])
    htg.SetZCoordinates(zValues)
else:
    zValues = vtkDoubleArray()
    zValues.SetNumberOfValues(1)
    zValues.SetValue(0,x_c[0])
    htg.SetZCoordinates(zValues)

# =========================================== #
# ===             FILL THE GRID           === #
# =========================================== #
# Let's split the various trees
cursor = vtkHyperTreeGridNonOrientedCursor()
offsetIndex  =0

# ROOT CELLS
ii=0
for i in range(N_base):
    for j in range(N_base):
        for k in range(N_base):
            #htg.GetIndexFromLevelZeroCoordinates(global_index_tree, i, j, k)
            htg.InitializeNonOrientedCursor(cursor, ii, True)
            cursor.SetGlobalIndexStart(offsetIndex)
            idx = cursor.GetGlobalNodeIndex()     
            iLevel = cursor.GetLevel()
            iBounds = cursor.GetGrid()
            print(iLevel)
            print(iBounds)
            refine = do_refine_root(0,i,j,k,iLevel,levelmax=max_level)
            if refine:
                refine_cell(cursor,0,0,0,max_level=max_level)
            else:
                scalarArray.InsertTuple1(idx, ii)
                vecArray.InsertComponent(idx, 0, ii)
                vecArray.InsertComponent(idx, 1, ii)
                vecArray.InsertComponent(idx, 2, ii)
            offsetIndex += cursor.GetTree().GetNumberOfVertices()
            ii += 1

# =========================================== #
# ===             SAVE THE GRID           === #
# =========================================== #
writer = vtk.vtkXMLHyperTreeGridWriter()
writer.SetFileName('test_'+str(dim)+'D.htg')
writer.SetInputData(htg)
writer.Write()

@Louis_Gombert @Jean_Fechter