Bug with TriangularSurfaceMesh?

Hello,

I am writing an application where I need to make a custom output format from free-surface simulation results. I tried to access vertex ID’s in the following way:

[code=“cpp”]
// Adapted from the writeResults function from the dam break example case
std::vector triangles;
isoSurfaceMarchingCube(triangles, fields->volumeFraction, isoLevels, fields->volumeFraction.getBoundingBox());
TriangleSet triangleSet(triangles);

DEFscaledMesh* defMesh = new DEFscaledMesh(triangleSet);
std::vector<Array<T,3> > vertexList = defMesh->getVertexList();
TriangularSurfaceMesh tempMesh = defMesh->getMesh();
delete defMesh;

tempI = tempMesh.getNumTriangles();
for (plint i = 0; i < tempI; i++) {
for (plint j = 0; j < 3; j++) {
int index = (int)(tempMesh.getVertexId(i, j));
gzwrite(gzf, &index, sizeof(int)); // Binary GZIP output from library
}
}



However, this resulted in wildly unpredictable index numbers (lots of zeros or +/- ~INT_MAX) in the first output files, but an increasing amount of correct numbers as the simulation progressed. I triple-checked the 32- vs. 64-bit sizes in binary output, but that wasn't the issue - rather, it sounded like some kind of memory allocation issue since it became fixed when the meshes were as large as they had ever been. Then, I dissected the TriangularSurfaceMesh source code and found a way around the problem using only DEFscaledMesh:

[code="cpp"]
...

DEFscaledMesh<T>* defMesh = new DEFscaledMesh<T>(triangleSet);
std::vector<Array<T,3> > vertexList = defMesh->getVertexList();
std::vector<Edge> edgeList = defMesh->getEdgeList();
delete defMesh;

...

tempI = triangles.size();
for (plint i = 0; i < tempI; i++) {
    for (plint j = 0; j < 3; j++) {
        int index = edgeList[3*i + (j+2)%3].pv; // Same algorithm as in the TriangularSurfaceMesh source code
        gzwrite(gzf, &index, sizeof(int));
    }
}

… and that worked perfectly from the start. Because I effectively did the same thing as in the original code, I assume that there is some kind of allocation bug in the TriangularSurfaceMesh code, but since I already found a work-around, I’m not digging any deeper. So, just thought I’d let you know!

Cheers,

  • AriH

Dear AriH,

thank you very much for the detailed post and your explanations. Since you are working directly with the “internals” of Palabos, let me say that what you observed is correct, but there is no bug in the implementation of TriangularSurfaceMesh (at least I believe so…). What happens is the following: the data of the TriangularSurfaceMesh are “owned” by the DEFscaledMesh, so in your first code extract, when you “delete defMesh” you corrupt the memory, and therefore you get the strange results. Your work-around is the correct way to go.

Thanks once again for your post!

Best,
Dimitris