Error exporting expression as a Silo file.
Greg Burton reported this bug. Exporting the data from an expression defining q criterion in a Silo file (HDF5) results in a file that can't be plotted by VisIt. If you export the data as a VTK file it works (ASCII XML).
The data is in /usr/workspace/brugger/greg_burton_data on poodle.
To reproduce the bug do the following.
- Open the file.
- Run in the debug queue with 16 processors (the data has 16 domains).
- Define the expression
q_critto beq_criterion(gradient(U[0]), gradient(U[1]), gradient(U[2])). - Plot a Pseudocolor plot of q_crit. This will take several minutes.
- Export the data as a Silo file named
q_crit1, writing it as HDF5. - Exit VisIt (don't save the expression).
- Restart VisIt.
- Open
q_crit1.silo. - Run in serial.
- Create a Pseudocolor plot of q_crit.
You will get the following error:
Pseudocolor: ()
viewer: Zone type -1 is not supported by VTK.
I then looked at the Silo file corresponding to the first domain. h5dump q_crit1.0.silo > q_crit1.0.txt.
Looking at the contents of q_crit1.0.txt file I found the following:
The portion showing the fields for the zonelist.
SOFTLINK "zonelist_nodelist" {
LINKTARGET "/.silo/#000001"
}
SOFTLINK "zonelist_shapecnt" {
LINKTARGET "/.silo/#000002"
}
SOFTLINK "zonelist_shapesize" {
LINKTARGET "/.silo/#000003"
}
SOFTLINK "zonelist_shapetype" {
LINKTARGET "/.silo/#000004"
}
If you then look at #000002 you see.
DATASET "#000002" {
DATATYPE H5T_STD_I32LE
DATASPACE SIMPLE { ( 25208 ) / ( 25208 ) }
DATA {
(0): 39, 1, 8, 1, 4, 1, 17, 1, 9, 1, 13, 1, 15, 1, 24, 1, 43, 1, 62,
(19): 1, 10, 1, 87, 1, 1527, 4, 4, 1, 1, 8, 1, 9, 1, 2, 1, 18, 1, 8,
(38): 1, 46, 1, 1, 1, 4, 20, 1, 11, 1, 25, 1, 29, 1, 1, 1, 1, 1, 11,
(57): 1, 16, 1, 2, 1, 1, 1, 5, 3, 6, 2, 2, 1, 3, 1, 2, 3, 13, 1, 6,
which looks good or at least not obviously incorrect. Now if you look at #000003 you see.
DATASET "#000003" {
DATATYPE H5T_STD_I32LE
DATASPACE SIMPLE { ( 25208 ) / ( 25208 ) }
DATA {
(0): 8, 13, 8, 9, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13,
(18): 8, 13, 8, 13, 8, 13, 8, 13, 8, 9, 13, 8, 13, 8, 13, 8, 13, 8,
(36): 13, 8, 13, 8, 13, 9, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 9,
(54): 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13, 8, 13,
which is definitely wrong. This is an all hex mesh, yet we have zone counts of 13 and 9. Now if you look at #000004 you see.
DATASET "#000004" {
DATATYPE H5T_STD_I32LE
DATASPACE SIMPLE { ( 25208 ) / ( 25208 ) }
DATA {
(0): 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1,
(16): 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, -1, 38, -1,
(31): 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, -1, 38, -1, 38, -1,
(46): 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38,
(61): -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1, 38, -1,
which is also definitely wrong. Wherever there was a zone count of 13 or 9, the zone type is -1. This is the cause of the error message.
Looking at the VTU file I noticed that the cell types that were of type -1 in the Silo file were of type 41 in the VTU file. Type 41 is a VTK_CONVEX_POINT_SET, which the header describes as
// Special class of cells formed by convex group of points
Mystery solved. It's a crazy unsupported cell type.
Great diagnosis, @brugger1. I wonder where that cell type is coming from. Its either something that gets generated in VTK filter we use or there is code in VisIt that creates it directly. So, I searched our sources...
mark@mypi5:~/visit/repo/src$ find . -type f -exec grep VTK_CONVEX_POINT_SET {} \; -print
interfaceCellType = VTK_CONVEX_POINT_SET;
//Mats[m].cellTypes.push_back( VTK_CONVEX_POINT_SET );
nextCell.type = VTK_CONVEX_POINT_SET;
./avt/MIR/Youngs/vtkYoungsMaterialInterface.C
case VTK_CONVEX_POINT_SET: SET_VAL(CPS); break;
./avt/Expressions/General/avtZoneTypeLabelExpression.C
case VTK_CONVEX_POINT_SET: rank = 3; break;
./avt/Expressions/General/avtZoneTypeRankExpression.C
if (ugrid != NULL && ugrid->GetCellType(i) == VTK_CONVEX_POINT_SET)
./visit_vtk/full/vtkVertexFilter.C
cerr << " N=0 (default): any shape, handled as general VTK_CONVEX_POINT_SET" << endl;
./databases/lata/LataFilter.C
type_cell=VTK_CONVEX_POINT_SET;
} else if (type_cell==VTK_CONVEX_POINT_SET) {
if (type_cell!=VTK_CONVEX_POINT_SET)
./databases/lata/avtlataFileFormat.C
int cellType = VTK_CONVEX_POINT_SET;
if (cellType == VTK_CONVEX_POINT_SET)
./databases/OpenFOAM/visit_vtkOpenFOAMReader.cxx
I don't think user is using Lata or OpenFOAM. Does q-crit use YoungMaterialInterface of ZoneTypeLabel or ZoneTypeRank?
The data files are VTU files, so in theory they could be in the file. If I open the VTU file corresponding to domain 0, do a mesh plot and export to VTK in ASCII XML I get a file with zones with zone type 41.
It is possible to read this datatype into VisIt using VTK but when we try to export to Silo the results aren't readable by VisIt. We should figure out why the Silo exporter doesn't error for this cell type.
@JustinPrivitera do you have an example VTK file with this cell type that VisIt reads and plots correctly? If so, can it be attached here?
@markcmiller86
The data is in
/usr/workspace/brugger/greg_burton_dataon poodle.
Ok, it looks like a VTK_CONVEX_POINT_SET is a special case of an arbitrary polyhedron. Apparently, with just the points alone, and nothing else, a VTK_CONVEX_POINT_SET specifies the faces (and edges) of the resulting points using an algorithm to extract the faces of the convex hull.
Now, Silo supports arbitrary polyhedra but that support involves specifying both a list of points and a list of faces. A list of points alone, even if producer promises they define a convex hull, are not sufficient to define a Silo arbitrary polyhedron.
I suspect Silo's behavior is the common case and the special case represented by VTK's VTK_CONVEX_POINT_SET is very uncommon across file formats that handle polyhedra.
Do we have policy and code enforcing that policy for when an exported database uses stuff the export format does not support? Other examples are mixing materials, adaptive mesh, high order elements, 2D topological/rectilinear surfaces in 3D space, etc.
My first thought is the general situation would be handled by an avtExportTransformManager that would do the necessary conversions (similar to existing avtTransformManager used when reading exotic stuff into VisIt). This would be better than having each and every writer that does not support, for example, VTK_CONVEX_POINT_SET having to do it themselves.
Setting doing any conversions aside, can we detect and warn user about the issue? To do that, we'd have to, just as for the case where we would actually do conversions, we have to have some relatively detailed knowledge about what formats do and do not support. This in itself is a complex problem because just coming up with a general description of format features alone involves high level data modeling know-how.
Best I think we can do in short term...
- Silo: Formalize Silo's
unknownzone type by addingDB_ZONETYPE_UNKNOWNtosilo.h - VisIt: Enhance Silo plugin to use missing data or ghost functionality to mask out any unknown zones.
- VisIt: Enhance Silo plugin to issue warning message about encountering unknown zones.
In this case, I think he was exporting to cache some results.
We tend to default to Silo, but using VTK which would already support this type of data is a good solution for this case.
I agree it would be good to handle transforms in one place, but this is a very niche case.
I think the best step is to add a warning to Silo, and maybe other commonly suggested outputs like Blueprint to let folks know their data can not be fully exported to what they requested.