6
$\begingroup$

Bug introduced in version 12.1.1 or earlier, and persisting through 13.2.0 or later.

Why would that not work:

L = 20;
h = 1;
Needs["NDSolve`FEM`"]
gr = Graphics3D[Cuboid[{-h/2, 0, -h/2}, {h/2, L, h/2}]]
mesh = ToElementMesh[DiscretizeGraphics[gr], "MeshOrder" -> 1]
(* TetGenTetrahedrize: returned error, 2 *)

when with a simple translation of $(h/2, 0, h/2)$ I get no error?

gr = Graphics3D[Cuboid[{0, 0, 0}, {h, L, h}]]
mesh = ToElementMesh[DiscretizeGraphics[gr], "MeshOrder" -> 1]
(* works *)

Edit Filled bug report [CASE:4609730]

Improve this question
$\endgroup$
4
  • 1
    $\begingroup$ You don't need Graphics3D[...] and DiscretizeGraphics[...]. Remove these and mesh["Wireframe"] works. $\endgroup$
    – flinty
    Commented Jul 15, 2020 at 19:31
  • 1
    $\begingroup$ @flinty Ok but it is not a strange behaviour still, that only one works? $\endgroup$
    – anderstood
    Commented Jul 15, 2020 at 19:35
  • $\begingroup$ Maybe - there are all sorts of problems like this in the region functionality. $\endgroup$
    – flinty
    Commented Jul 15, 2020 at 19:59
  • $\begingroup$ Seems a bug to me. $\endgroup$
    – Michael E2
    Commented Jul 15, 2020 at 20:27

1 Answer 1

Reset to default
7
$\begingroup$

It seems something with the automatic MaxCellMeasure computation is not quite right. Here is a workaround with a manually specified MaxCellMeasure

mr = DiscretizeGraphics[gr];
(*ToElementMesh[mr]*)
ToElementMesh[mr, MaxCellMeasure -> .1]

Update:

This is a bug in TetGen, the mesh generator used.

(bmesh = ToBoundaryMesh[mr]);
pts = bmesh["Coordinates"];
facets = Partition[
   Join @@ ElementIncidents[bmesh["BoundaryElements"]], 1];
Graphics3D[GraphicsComplex[pts, Map[Polygon, facets]]]

enter image description here

In the below code, if you change the 0.06 to 0.05, TetGen will fail.

Needs["TetGenLink`"]
inInst = TetGenCreate[];
TetGenSetPoints[inInst, pts];
TetGenSetFacets[inInst, facets];

(* changing "pqa0.006" to "pqa0.005" *)
outInst = TetGenTetrahedralize[inInst, "pqa0.006"];
coords = TetGenGetPoints[outInst];
surface = TetGenGetFaces[outInst];
Graphics3D[GraphicsComplex[coords, Polygon[surface]]]

enter image description here

Since TetGen is software that Wolfram Research purchased for inclusion there is little hope that this bug in TetGen will be fixed by myself or other people at Wolfram Research. Since the workaround is easy:

ToElementMesh[Cuboid[{-h/2, 0, -h/2}, {h/2, L, h/2}], 
 "MeshElementType" -> TetrahedronElement]

I do not think this is a critical bug. If you deeply care about this, you could try to make use of the TetGenLink that ships with Wolfram products and compile a newer version of TetGen and link and see if the bug has been fixed. Wolfram Products use TetGen version 1.4.3 and the current version is 1.5.1.

Improve this answer
$\endgroup$
5
  • $\begingroup$ Nice suggestion to the user who didn't purchase TetGen but did purchase or receive his copy of Mathematica :( sorry for the sarcasm :( $\endgroup$
    – Pavel Perikov
    Commented Jul 16, 2020 at 13:13
  • $\begingroup$ @PavelPerikov, and what do you suggest I do? $\endgroup$
    – user21
    Commented Jul 16, 2020 at 13:16
  • $\begingroup$ you did great. Sorry if my comment sounded offensive in any way. I think your explanation was delightful and useful. $\endgroup$
    – Pavel Perikov
    Commented Jul 16, 2020 at 14:05
  • 1
    $\begingroup$ @PavelPerikov, thanks for the clarification - I just wanted to be honest here, this is a bug and I should take care of it. At the same time I also want to be honest about my assessment of the difficulty/importance to fix this issue by saying that it is unlikely that I will fix this seen it is external code and not that critical. Now, I could be wrong. If enough people feel that this a vital thing to fix then the best way to proceed is to put it on this list and I will try to make it happen. $\endgroup$
    – user21
    Commented Jul 16, 2020 at 14:14
  • $\begingroup$ Wow. Thanks for the link!!! $\endgroup$
    – Pavel Perikov
    Commented Jul 16, 2020 at 14:28

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.