tetrahedral meshing of very thin box

Hi,

I have an application where there can be simulatenously distances in the mm range, together with features and thin layers in the range of nm. This causes several issues in netgen:

1) In the attached example (thin.geo), the meshing will fail completely if the boundingbox statement is commented out. The reason is that netgen uses a default bounding box [-1000, 1000]^3 if the boundingbox statement is missing.
2) Meshing the attached example (thin.geo) with the default "moderate" mesh granilarity will fail with the error message "ERROR: Problem in Surface mesh generation".
3) Meshing the attached example (thin.geo) with "very fine" mesh granilarity will cause a crash during the first SwapImprove step.

The attached patch (improve3.cpp) fixes that crash. However, it probably makes sense to apply that patch in two separate steps, first only the part which detects whether to topology is already broken, and then the part which tries to avoid that SwapImprove itself breaks the topology.

It can also happen (not reproducible with the attached file) that the topology is already broken when SwapImprove is called the first time. I have no patch for this yet, hopefully this is related to the issue 2) above.

Another issue (not reproducible with the attached file) was that GeomSearch3d :: Create exceeded the 32 GB memory of my computer and started to swap during creation of the cell based search structure. The attached patch (geomsearch.cpp) is how I fixed this issue for me.


I will probably continue to debug netgen to find out whether the remaining issues can be fixed too, but it made sense to me to first report what I have found out so far. (I have also created a github account, but just attaching the patches seemed easier to me at the moment, since they are not "polished" anyway.)


Hope reporting those issues here is fine, and maybe I could get some help at some later point, if I should get totally stuck while debugging.

> It can also happen (not reproducible with the attached file) that the topology is already broken when SwapImprove is called the first time.

This is actually not true, I just didn't notice the SwapImprove call which broke it. I attached a simple test case (high.zip, contains high.geo and high.msz, i.e. a mesh-size-file) which triggers the problem quite early. (It is a normal box, no extreme geometry this time, just a fine local resolution).

I have a better fix for SwapImprove now (patch.zip, contains meshtool.hpp.patch, meshtool.cpp.patch, and improve3.cpp.patch), which directly fixes what actually goes wrong in my case. I replaced 1e24 in CalcTetBadness by a macro INFINIBAD, which I defined to be 1e48.

Even with this patch, I can still run into a broken topology in SwapImprove sometimes, especially if I use "very fine" mesh granularity to work around the "ERROR: Problem in Surface mesh generation" reported under 2) in the previous message. Probably will try to find a solution for that issue next, before trying to figure out why the topology can still be boken in SwapImprove, despite my cool patch.

Hi,
the error "ERROR: Problem in Surface mesh generation" can be avoided if I specify maxh properties for the planes on the side. It occurs because Netgen fails in reducing the meshsize enough there and then stops with that error. When the surface mesh is fine, then meshing works for me even without your patches. Do you have a reason why you don't want prismatic elements for the thin layer? They would have nicer properties... I have attached a file generating a prismatic mesh using the CloseSurfaces function. Both are Python files, but they should translate easily to geo files as well.

Best
Christopher

Thanks for the help.

Specifying -maxh for the planes on the side basically works, at least with respect to avoiding the "ERROR: Problem in Surface mesh generation". However, automatic surface identification seems to interact with this feature in hard to control ways. On top of the thin box to which -maxh is applied, there is a huge box. But depending on which box is defined first, the final mesh gets a totally different density. I attached two larger examples (maxh_upper.geo, maxh_lower.geo), where the effect is easy to see.

I tried to fix this with the attached patch (maxh.patch.txt), but the density on some edges still depends on the orders, but I don't know exactly why yet, i.e. where exactly it goes wrong. (Maybe segment.surfnr1 and segment.surfnr2 are not as reliable as they seem.)


The functionality with the prismatic elements is nice, but it probably won't help me. My main structure is normally embedded into the thin layer. (I tested whether I could use prisms surrounding my main structure, but netgen then stops because the surface mesh gets inconsistent.) I attached a real non-simplified example from my application (structure.geo, it still fails in SwapImprove in the end, even with my patches). All that huge radial symmetric part is just a poor man's boundary condition, in a certain sense. (I do have a solver to that radial symmetric part, just the coupling to the FEM solver has room for improvement.)

> The functionality with the prismatic elements is nice, but it probably won't help me. ... (I tested whether I could use prisms surrounding my main structure, but netgen then stops because the surface mesh gets inconsistent.)

I did some more tests (attached as prism_tests.zip). The trick is to create a separate top level object for each prism layer region. Otherwise it just fails, even for trivial geometries.

I new tried to generate a border with prism layers for my application (attached an example structure.geo), and it seems to work well so far.