mesh_top1 and large numbers of geometric model edges/points
Opened this issue · 10 comments
Hello,
I'm trying to generate a mesh using jigsaw with a geometric model that has roughly 16k edges and hitting problems with memory consumption and runtime when enabling mesh_top1.
The test repo here:
https://github.com/cwsmith/jigsaw_gis/tree/519e7e1437fbabbbe3037cee4212721f079bddb4
contains gis.py that loads the geometric model (.msh) and size field (hfun) and runs jigsaw in a few seconds to generate a mesh of 358k triangles that doesn't exactly 'recover' the geometric model edges/vertices in the resulting mesh as mesh_top1 is disabled.
When I add opts.mesh_top1 = True to gis.py jigsaw runs out of memory on a gnu/linux machine with 64GB of RAM after about 10mins. The log from the failed run is pasted below.
I'm a novice jigsaw user and am concerned I have a problematic combination of settings and mesh size field that is causing the large memory usage and long runtime. Do you see anything obviously wrong? Any suggestions are appreciated.
Thank you.
verbose output from failed run
Call libJIGSAW: gis
#------------------------------------------------------------
#
# ,o, ,o, /
# ` ` e88~88e d88~\ /~~~8e Y88b e /
# 888 888 88 88 C888 88b Y88b d8b /
# 888 888 "8b_d8" Y88b e88~-888 Y888/Y88b/
# 888 888 / 888D C88 888 Y8/ Y8/
# 88P 888 Cb \_88P "8b_-888 Y Y
# \_8" Y8""8D
#
#------------------------------------------------------------
# JIGSAW: an unstructured mesh generation library.
#------------------------------------------------------------
JIGSAW VERSION 1.0.0
Reading CFG. data...
CFG. data summary...
GEOM-FILE =
MESH-FILE =
HFUN-FILE =
INIT-FILE =
TRIA-FILE =
BNDS-FILE =
NUMTHREAD = 20
GEOM-SEED = 8
GEOM-PHI1 = 6.00e+01
GEOM-PHI2 = 6.00e+01
GEOM-ETA1 = 5.00e+00
GEOM-ETA2 = 5.00e+00
GEOM-FEAT = TRUE
INIT-NEAR = 1.00e-08
HFUN-SCAL = ABSOLUTE
HFUN-HMAX = inf
HFUN-HMIN = 0.00e+00
BNDS-KERN = BND-TRIA
MESH-KERN = DELFRONT
MESH-TOP1 = TRUE
MESH-TOP2 = FALSE
MESH-ITER = MAXINT
MESH-DIMS = 2
MESH-SIZ1 = 1.33e+00
MESH-SIZ2 = 1.31e+00
MESH-SIZ2 = 1.30e+00
MESH-EPS1 = 3.33e-01
MESH-EPS2 = 3.33e-01
MESH-RAD2 = 1.05e+00
MESH-RAD3 = 2.05e+00
MESH-OFF2 = 9.00e-01
MESH-OFF3 = 1.10e+00
MESH-SNK2 = 2.00e-01
MESH-SNK3 = 3.33e-01
MESH-VOL3 = 0.00e+00
OPTM-KERN = ODT+DQDX
OPTM-ITER = 16
OPTM-COST = AREA-LEN
OPTM-BETA = 4.95e-01
OPTM-ZETA = 8.25e-01
OPTM-QTOL = 1.00e-04
OPTM-QLIM = 9.38e-01
OPTM-ZIP_ = TRUE
OPTM-DIV_ = TRUE
OPTM-TRIA = TRUE
OPTM-DUAL = FALSE
Done. (6.90e-05sec)
#------------------------------------------------------------
Reading GEOM data...
Done. (1.36e-03sec)
#------------------------------------------------------------
Forming GEOM data...
GEOM data summary...
EUCLIDEAN-MESH
|NDIMS.| = 2
|COORD.| = 16675
|EDGE-2| = 16674
|SEEDS.| = 0
|BOUND.| = 0 (0)
Done. (3.66e-03sec)
#------------------------------------------------------------
Reading HFUN data...
Done. (1.54e-03sec)
#------------------------------------------------------------
Forming HFUN data...
HFUN data summary...
EUCLIDEAN-GRID
|NDIMS.| = 2
.MIN(H). = 3.00e+03
.MAX(H). = 3.00e+04
|MASKED| = 0
|XGRID.| = 832
|YGRID.| = 1408
Done. (1.15e-03sec)
#------------------------------------------------------------
Generate rDT MESH...
#------------------------------------------------------------
# |ITER.| |DEL-1| |DEL-2|
#------------------------------------------------------------
10000 16485 28345
25000 16485 58345
50000 16485 108345
75000 16485 158345
100000 16807 207745
250000 25094 472549
500000 25236 930619
750000 25321 1389276
1000000 25496 1847653
2500000 25957 4583110
5000000 26593 9142783
7500000 27115 13692413
10000000 27355 18256201
25000000 28722 45592494
50000000 30873 91114451
75000000 32284 136651641
100000000 33472 182152897
250000000 38587 455269920
Terminated
I think I see a problem with the input data... closing this for now while I dig deeper.
@cwsmith let me know if you're still having issues on this --- when you ask for mesh_top1 = TRUE though, jigsaw really will try to recover the topology of the geometry, even in cases that might not make sense (self intersections, duplicates, etc, etc). This is hopefully the root of the non-convergence you're seeing.
Will do. Thanks for checking in.
The input geometric model was not in the same coordinate frame (different bounding box and scale) as the size field (hfun) so I was effectively getting two meshes and that was definitely not what was wanted. I'm working on a fix now and will let you know if there are problems.
@dengwirda I fixed the coordinate frame issue and ran what I think is a valid configuration (1, below) to use a size field (hfun) and recover the geometric model vertices (mesh_top1). Unfortunately, I'm seeing long run times and high memory usage. I was unsure if feature detection (geom_feat) and/or the size field was somehow causing problems with recovery so I attempted (without much success) to disable them. Details on those additional configurations (2,3) are also listed below.
Any suggestions are appreciated. Thanks in advance.
1 - feature detection enabled, top1 enabled, hfun enabled
result: killed after 45 mins, memory usage was about 50GB
code: https://github.com/cwsmith/jigsaw_gis/blob/72921f669aee1fea4b2eb89e7d240b422506e7ce/gis.py
log: https://github.com/cwsmith/jigsaw_gis/blob/72921f669aee1fea4b2eb89e7d240b422506e7ce/log
2 - feature detection enabled, top1 enabled, hfun disabled
result: killed after 35 mins, memory usage was about 60GB
code: https://github.com/cwsmith/jigsaw_gis/blob/21a49d1960cfd4aed8ec23baba9fdbf182d55f75/gis.py
log: https://github.com/cwsmith/jigsaw_gis/blob/21a49d1960cfd4aed8ec23baba9fdbf182d55f75/feat5deg_top1.log
3 - feature detection disabled, top1 enabled, hfun disabled
result: mesh generated but does not recover the geometric model vertices/edges, there is also a warning at the end of execution
code: https://github.com/cwsmith/jigsaw_gis/blob/8d8f0a3ae6f0744aae5ca36483eaff9006e786de/gis.py
log: https://github.com/cwsmith/jigsaw_gis/blob/8d8f0a3ae6f0744aae5ca36483eaff9006e786de/top1.log
screenshot: https://github.com/cwsmith/jigsaw_gis/blob/8d8f0a3ae6f0744aae5ca36483eaff9006e786de/top1.png
- created by loading input model entities and 'gis.vtk' (written by
gis.pyafter mesh generation) - black points are vertices from the input geometric model
- all the input geometric model vertices were not recovered in the mesh - it seems like feature detection is required
@dengwirda
I found two more issues with the .msh in the test repo:
- it has two points that are nearly coincident, and
- there are two model edge loops (sequences of adjacent edges that begins and ends at the same vertex) and the second one references a duplicate first vertex instead of closing back on the correct first vertex.
Closing this again while I clean up the model. Sorry for the noise.
Success! Cleaning up the model issues above, and, for now, removing the BOUND data section, gives me the result I was looking for:
the red dots/points are vertices of the input geometric model. This mesh was generated with cwsmith/jigsaw_gis@454dea1.
Looks good to me @cwsmith --- with default settings an approximation to the geometry will be made based on the sizing function alone (as per your figure on the right) but you can force additional refinement (to recover topology, geometric accuracy, etc, etc) with other settings, as per the mesh on the left.
@dengwirda I'm hitting problems (long time and high memory use in the Generate rDT MESH...) when I attempted to add a bounding rectangle to the above and run with mesh_top1. Do you see anything obviously wrong with my modifications (diff is pasted below) to the .msh file?
original geometry, meshes with mesh_top1 enabled, see #54 (comment): mesh_12_14.msh
geometry with bounding rectangle: mesh_12_14_withBdry.msh
If not, is there a way to extract some debug info from jigsaw about where the geometry may be giving problems to (I think it is stuck in mesh::rdel_mesh_2d(...))?
Thanks in advance.
(ins)(pyEnv) cwsmith@checkers: /space/cwsmith/compassLandice/jigsaw_gis (main)$ diff mesh_12_14.msh mesh_12_14_withBdry.msh
4c4
< POINT=16665
---
> POINT=16669
16670c16670,16674
< EDGE2=16665
---
> -713000;-3396000;0
> 949000;-3396000;0
> 949000;-582000;0
> -713000;-582000;0
> EDGE2=16669
33335a33340,33348
> 16665;16666;1
> 16666;16667;1
> 16667;16668;1
> 16668;16665;1
> BOUND=4
> 1;16665;20
> 1;16666;20
> 1;16667;20
> 1;16668;20
@cwsmith: interesting --- I was able to reproduce the non-convergence and (at a first look at least) it actually appears to be related to setting geom_eta1 = 5.0, which controls what's considered to be a "sharp" feature in the geometry. With geom_eta1 left at its default I didn't run into any convergence issues. Your geometry appears to be correctly defined to mesh both the land and ocean part of the domain.
I'm not 100% on whether this behaviour is a bug or not at this stage, so will keep this open for now.
@dengwirda Thank you.
Is there a way to 'restrict' the mesh to conform to the geometry and 'disable' sharp feature detection?
The screenshot below shows three meshes with mesh_top1 enabled for all but with different settings of geom_feat and geom_eta1. With geom_feat=False some features are missed/coarsened out and, when its enabled, a relatively low geom_eta1 setting is needed to avoid some coarsening.
