This is the second part of a series of blog articles about vtkHyperTreeGrid usage and implementation in VTK. The first part, an introduction about HTG, can be found here, the third part, HTG: Using Masks, can be found here, the fourth part, HTG: Specific Filters, can be found here.
After quickly focusing on the TB-AMR (tree-based Adaptive Mesh Refinement) technique, we are now exploring the construction of a HyperTreeGrid (the TB-AMR representation in VTK).
Let’s illustrate the construction of such a representation by following an example written in Python:
First step : Declare and construct the rectilinear grid
The definition of a vtkHyperTreeGrid in Python is done as follows:
import vtk htg = vtk.vtkHyperTreeGrid()
We continue with the description of the global rectilinear grid:
Then the coordinates of the points along each axis (here, TB-AMR 2D), starting with:
Where xValues is a vtkDoubleArray of 4 values. Respectively, this work needs to be done along Y (3 values) by SetYCoordinates and Z (one value) by SetZCoordinates.
The branch factor still needs to be defined at the global level by the following call:
Second step : Declare and set a cursor on a cell of the rectilinear grid
One of the construction methods is to go only through each root cell of this grid once in order to describe the refined tree of this cell, if necessary.
To do so, we create:
- A non-oriented cursor. Non-oriented because we want to be able to move from a cell to a child cell:
ToChild(iChild), or to the parent cell:
- The global offset of the root cell tree being processed.
cursor = htg.vtkHyperTreeGridNonOrientedCursor() globalOffset = 0
To identify a root cell in the grid by a global index tree, we can use the method:
htg.GetIndexFromLevelZeroCoordinates(global_index_tree, index_tree_i, index_tree_j, index_tree_k)
index_tree_i (resp. j and k) starting at 0 as in C language.
We then need to position the cursor following the first cell of this decomposition tree as follows:
htg.InitializeNonOrientedCursor(cursor, global_index_tree, True)
create parameter is set here to
True to indicate that we are requesting the creation of a decomposition tree. It does not need to be specified (
False by default) when creating a cursor to explore this same tree.
Still in the case of creation, we must give the global offset of the root cell tree.
The cursor then points to the root cell tree. Initially, this cell is created and considered as leaf cell, unrefined.
Third step : Use this cursor to build the tree with successive refinements
From there, we can decide to subdivide this cell:
Then we need to move to one of these daughter cells :
Where iChild is the index of the child cell, which must be between 0 and b^d, b being the branch factor and d being the number of spatial dimension..
It is also possible to move to the parent cell:
At any time, information is accessible from this cursor:
- The global offset, which identifies the index in an array where the value of the current cell must be stored, or how to access its value,
idx = cursor.GetGlobalNodeIndex() scalar.InsertTuple1(idx, val) scalar.GetTuple1(idx, val)
- The current level (with level 0 as the root cell tree):
iLevel = cursor.GetLevel()
Without having to move the cursor to the root cell tree, you can decide to finish building this refined tree. All you have to do is to go to the next tree that has not yet been defined.
Warning: This construction method does not allow you to complete an already described TB-AMR, as it requires a global index definition which will be covered in the next article.
But first, you must update the global offset of the root cell tree from the construction of the current tree :
globalOffset += cursor.GetTree().GetNumberOfVertices()
Fourth step : Move the cursor to another rectilinear grid cell, an unrefined cell.
Before positioning on the next tree to be subdivided down by resetting the cursor:
htg.GetIndexFromLevelZeroCoordinates(new_global_index_tree, new_index_tree_i, new_index_tree_j, new_index_tree_k) htg.InitializeNonOrientedCursor(cursor, new_global_index_tree, True) cursor.SetGlobalIndexStart(globalOffset)
Then we continue as in step 3.
Build view: an unstructured representation
To view this representation, you can convert it to an unstructured representation with the following call:
geometry = vtk.vtkHyperTreeGridGeometry() geometry.SetInputData(htg)
- The topology of the unstructured mesh produced is not correct but is more than enough to perform the rendering
- In 3D, only the mesh surface is rendered.
Currently, in 2D, this geometry filter can take a camera as a parameter in order to control an adapted building of the unstructured mesh. This filter roughly builds a pixel level or offscreen cells.
As we arrive to the end of this second article, we need to send a warning: The vtkHyperTreeGrid representation is compact; hence a filter application producing another representation as output (unstructured for example) will generate a memory blast. This is why it is preferable to use, or create, specific filters for this representation.
Simple yet complete example
Here is a small complete Python code whose result will give:
#!/usr/bin/env python import vtk # Declare htg htg = vtk.vtkHyperTreeGrid() htg.Initialize() # Declare scalar scalarArray = vtk.vtkDoubleArray() scalarArray.SetName('MyScalar') scalarArray.SetNumberOfValues(0) htg.GetCellData().AddArray(scalarArray) htg.GetCellData().SetActiveScalars('MyScalar') # Construct the rectilinear grid htg.SetDimensions([4, 3, 1]) htg.SetBranchFactor(2) ## X coordinates xValues = vtk.vtkDoubleArray() xValues.SetNumberOfValues(4) xValues.SetValue(0, -1) xValues.SetValue(1, 0) xValues.SetValue(2, 1) xValues.SetValue(3, 2) htg.SetXCoordinates(xValues) ## Y coordinates yValues = vtk.vtkDoubleArray() yValues.SetNumberOfValues(3) yValues.SetValue(0, -1) yValues.SetValue(1, 0) yValues.SetValue(2, 1) htg.SetYCoordinates(yValues) ## Z coordinates zValues = vtk.vtkDoubleArray() zValues.SetNumberOfValues(1) zValues.SetValue(0, 0) htg.SetZCoordinates(zValues) # Declare cursor cursor = vtk.vtkHyperTreeGridNonOrientedCursor() # Initialise values offset offsetIndex = 0 # Set cursor on cell #0 and set the start index htg.InitializeNonOrientedCursor(cursor, 0, True) cursor.SetGlobalIndexStart(offsetIndex) # cell 0 # Assigns an index for the value of the cell pointed to by the current cursor state idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 1) # Decides to subdivide cell pointed by current cursor state cursor.SubdivideLeaf() # cell 0 # In this example, 4 child cells were constructed # Move the cursor to the first child cursor.ToChild(0) # cell 0/0 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 7) cursor.ToParent() # cell 0 # Move the cursor to the second child cursor.ToChild(1) # cell 0/1 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 8) cursor.ToParent() # cell 0 # And next cursor.ToChild(2) # cell 0/2 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 9) cursor.ToParent() # cell 0 # And next cursor.ToChild(3) # cell 0/3 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 10) cursor.ToParent() # cell 0 # Refine the cell 0 offsetIndex += cursor.GetTree().GetNumberOfVertices() # Set cursor on cell #3 and set the start index htg.InitializeNonOrientedCursor(cursor, 1, True) cursor.SetGlobalIndexStart(offsetIndex) # cell 1 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 2) # Refined the cell 1 offsetIndex += cursor.GetTree().GetNumberOfVertices() # Set cursor on cell #3 and set the start index htg.InitializeNonOrientedCursor(cursor, 3, True) cursor.SetGlobalIndexStart(offsetIndex) # cell 3 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 4) # Refined the cell 3 offsetIndex += cursor.GetTree().GetNumberOfVertices() # Set cursor on cell #5 and set the start index htg.InitializeNonOrientedCursor(cursor, 5, True) cursor.SetGlobalIndexStart(offsetIndex) # cell 5 idx = cursor.GetGlobalNodeIndex() scalarArray.InsertTuple1(idx, 6) # Refined the cell 5 offsetIndex += cursor.GetTree().GetNumberOfVertices() # Which ends the construction of my mesh # Use a Geometry filter geometry = vtk.vtkHyperTreeGridGeometry() geometry.SetInputData(htg) # Add, for example, a Shrink filter shrink = vtk.vtkShrinkFilter() shrink.SetInputConnection(geometry.GetOutputPort()) shrink.SetShrinkFactor(.8) # LookupTable lut = vtk.vtkLookupTable() lut.SetHueRange(0.66, 0) lut.UsingLogScale() lut.Build() # Mappers mapper = vtk.vtkDataSetMapper() mapper.SetInputConnection(shrink.GetOutputPort()) mapper.SetLookupTable(lut) mapper.SetColorModeToMapScalars() mapper.SetScalarModeToUseCellFieldData() mapper.SelectColorArray('MyScalar') # Actors actor = vtk.vtkActor() actor.SetMapper(mapper) # Renderer renderer = vtk.vtkRenderer() renderer.SetActiveCamera(camera) renderer.AddActor(actor) # Render window renWin = vtk.vtkRenderWindow() renWin.AddRenderer(renderer) renWin.SetSize(600, 400) # Render window interactor iren = vtk.vtkRenderWindowInteractor() iren.SetRenderWindow(renWin) # render the image renWin.Render() iren.Start()
Compatible with VTK 9.0.1 and ParaView 5.9.1
A complete example in Python was made by Sébastien Jourdain during a “Hackathon” in January 2019.
This work was supported by CEA
CEA, DAM, DIF, F-91297 Arpajon, France