Filter vtkSubdivideTetra results in negative volumes

When using the vtkSubdivideTetra filter, half of the resulting tetrahedral cells have negative volume, resulting in a zero or close to zero volume for the whole mesh.
This probably has to do with with the ordering of vertex references of the newly created cells, resulting in some “inverted” cells.

python example below:

import vtk
import numpy as np
from vtkmodules.vtkCommonDataModel import VTK_TETRA
from vtkmodules.vtkCommonCore import vtkPoints
from vtkmodules.vtkFiltersModeling import vtkSubdivideTetra
from vtkmodules.vtkFiltersVerdict import vtkCellSizeFilter

# construct UnstructuredGrid with single tetrahedral cell
coords = [[0.0, 0.0, 0.0],
          [2.0, 0.0, 0.0],
          [0.0, 2.0, 0.0],
          [0.0, 0.0, 1.5]]
points = vtkPoints()
for i in range(0, len(coords)):
    points.InsertPoint(i, coords[i])
mesh = vtk.vtkUnstructuredGrid()
mesh.Allocate(1)
mesh.InsertNextCell(VTK_TETRA, 4, [0, 1, 2, 3])
mesh.SetPoints(points)

# subdivide into 12 smaller tetrahedral cells
filter_subdivision = vtkSubdivideTetra()
filter_subdivision.SetInputData(mesh)
subdivided = vtk.vtkUnstructuredGrid()
subdivided.Allocate(12)
filter_subdivision.SetOutput(subdivided)
filter_subdivision.Update()

# compute cell size [single tet cell]
filter_cellsize = vtkCellSizeFilter()
filter_cellsize.SetInputData(mesh)
filter_cellsize.Update()
mesh_sized = filter_cellsize.GetOutput()
data = mesh_sized.GetCellData()

data_idx = 3
if not data.GetArrayName(data_idx) == "Volume":
    raise ValueError(f"Expected volume data at array of index {data_idx}")

# compute cell size [subdivided tet cells]
filter_cellsize = vtkCellSizeFilter()
filter_cellsize.SetInputData(subdivided)
filter_cellsize.Update()
subdivided_sized = filter_cellsize.GetOutput()
subdivided_sized.GetNumberOfCells()
subdivided_data = subdivided_sized.GetCellData()
subdivided_volumes = [subdivided_data.GetArray(data_idx).GetValue(i)
                      for i in range(subdivided_data.GetArray(data_idx).GetNumberOfValues())]

print("single tetrahedral cell:")
print(f"\tnumber of cells: {mesh_sized.GetNumberOfCells()}")
print(f"\tcell {data.GetArrayName(data_idx)}(s): "
      f"{[data.GetArray(data_idx).GetValue(i) for i in range(data.GetArray(data_idx).GetNumberOfValues())]}")

print("subdivided tetrahedral cells:")
print(f"\tnumber of cells: {subdivided_sized.GetNumberOfCells()}")
print(f"\tcell {subdivided_data.GetArrayName(data_idx)}(s): {subdivided_volumes}")
print(f"\ttotal volume: {sum(subdivided_volumes)}")
print(f"\tabsolute volume: {np.abs(subdivided_volumes).sum()}")

Output of the example:

single tetrahedral cell:
number of cells: 1
cell Volume(s): [1.0]
subdivided tetrahedral cells:
number of cells: 12
cell Volume(s): [0.125, -0.125, 0.125, -0.125, -0.0625, 0.0625, -0.0625, 0.0625, -0.0625, 0.0625, 0.0625, -0.0625]
total volume: 0.0
absolute volume: 1.0