trixi-framework / Trixi.jl

Trixi.jl: Adaptive high-order numerical simulations of conservation laws in Julia
https://trixi-framework.github.io/Trixi.jl
MIT License
540 stars 110 forks source link

GMSH curved meshes support #1847

Open Arpit-Babbar opened 9 months ago

Arpit-Babbar commented 9 months ago

Thanks to Daniel and everyone involved for bringing the wonderful GMSH mesh support to Trixi in https://github.com/trixi-framework/Trixi.jl/pull/1799!!

I am unable to get the current parser to handle a GMSH generated mesh with curved boundaries. Here is how I generate a curved mesh with GMSH (I only know how to do it with GUI at this moment). I tried this with NACA6412.geo available here and another one here

1) Open the .geo file in GMSH 2) Click Modules -> Mesh -> 2D (Generates the straight sided mesh) 3) Click Modules -> Mesh -> Set order 2 (curves the mesh boundary to be order 2. Higher is also possible) 4) File -> Export -> mesh_NACA6412.inp, making sure to tick "Save all elements" and "Save groups of nodes" (Exporting with these options ticked is the step I cannot do via command line)

If I skip step 3, the parser can handle it fine and I am able to run the available elixir https://github.com/trixi-framework/Trixi.jl/blob/main/examples/p4est_2d_dgsem/elixir_euler_NACA6412airfoil_mach2.jl. However, if I try to run same the elixir when I did the above step 3 to curve the mesh, I get a segmentation fault in the line that reads the mesh, i.e., mesh = P4estMesh{2}(mesh_file, polydeg = polydeg, boundary_symbols = boundary_symbols)

It starts with

[p4est] No elements or nodes found in mesh file.
[p4est] Failed to read /home2/arpit/repositories/cloned/Trixi.jl/examples/p4est_2d_dgsem/NACA6412_curved.inp: pass 1
Full error message P4est segmentation fault ```shell [p4est] No elements or nodes found in mesh file. [p4est] Failed to read /home2/arpit/repositories/cloned/Trixi.jl/examples/p4est_2d_dgsem/NACA6412_curved.inp: pass 1 [595759] signal (11.1): Segmentation fault in expression starting at /home2/arpit/repositories/cloned/Trixi.jl/examples/p4est_2d_dgsem/elixir_euler_NACA6412airfoil_mach2.jl:78 unsafe_load at ./pointer.jl:119 [inlined] unsafe_load at /home2/arpit/.julia/packages/P4est/0rTeV/src/pointerwrappers.jl:86 [inlined] getindex at /home2/arpit/.julia/packages/P4est/0rTeV/src/pointerwrappers.jl:82 [inlined] getindex at /home2/arpit/.julia/packages/P4est/0rTeV/src/pointerwrappers.jl:82 [inlined] p4est_mesh_from_standard_abaqus at /home2/arpit/repositories/cloned/Trixi.jl/src/meshes/p4est_mesh.jl:460 #_#279 at /home2/arpit/repositories/cloned/Trixi.jl/src/meshes/p4est_mesh.jl:376 P4estMesh at /home2/arpit/repositories/cloned/Trixi.jl/src/meshes/p4est_mesh.jl:353 unknown function (ip: 0x7f420e7c66a0) _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 jl_apply at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/julia.h:1982 [inlined] do_call at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:126 eval_value at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:223 eval_stmt_value at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:174 [inlined] eval_body at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:617 jl_interpret_toplevel_thunk at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:775 jl_toplevel_eval_flex at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:934 jl_toplevel_eval_flex at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:877 ijl_toplevel_eval_in at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:985 eval at ./boot.jl:385 [inlined] include_string at ./loading.jl:2070 _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 _include at ./loading.jl:2130 include at ./client.jl:489 unknown function (ip: 0x7f420e7b5b55) _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 jl_apply at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/julia.h:1982 [inlined] do_call at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:126 eval_value at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:223 eval_stmt_value at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:174 [inlined] eval_body at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:617 jl_interpret_toplevel_thunk at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/interpreter.c:775 jl_toplevel_eval_flex at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:934 jl_toplevel_eval_flex at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:877 jl_toplevel_eval_flex at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:877 jl_toplevel_eval_flex at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:877 ijl_toplevel_eval_in at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/toplevel.c:985 eval at ./boot.jl:385 [inlined] eval_user_input at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/REPL/src/REPL.jl:150 repl_backend_loop at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/REPL/src/REPL.jl:246 #start_repl_backend#46 at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/REPL/src/REPL.jl:231 start_repl_backend at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/REPL/src/REPL.jl:228 _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 #run_repl#59 at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/REPL/src/REPL.jl:389 run_repl at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/usr/share/julia/stdlib/v1.10/REPL/src/REPL.jl:375 jfptr_run_repl_91689.1 at /home2/arpit/.local/julia-1.10.0/lib/julia/sys.so (unknown line) _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 #1013 at ./client.jl:432 jfptr_YY.1013_82677.1 at /home2/arpit/.local/julia-1.10.0/lib/julia/sys.so (unknown line) _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 jl_apply at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/julia.h:1982 [inlined] jl_f__call_latest at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/builtins.c:812 #invokelatest#2 at ./essentials.jl:887 [inlined] invokelatest at ./essentials.jl:884 [inlined] run_main_repl at ./client.jl:416 exec_options at ./client.jl:333 _start at ./client.jl:552 jfptr__start_82703.1 at /home2/arpit/.local/julia-1.10.0/lib/julia/sys.so (unknown line) _jl_invoke at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:2894 [inlined] ijl_apply_generic at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/gf.c:3076 jl_apply at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/julia.h:1982 [inlined] true_main at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/jlapi.c:582 jl_repl_entrypoint at /cache/build/builder-amdci4-6/julialang/julia-release-1-dot-10/src/jlapi.c:731 main at julia (unknown line) __libc_start_main at /lib/x86_64-linux-gnu/libc.so.6 (unknown line) unknown function (ip: 0x4010b8) Allocations: 25176877 (Pool: 25169453; Big: 7424); GC: 90 Segmentation fault ```

Note that the .geo file does have enough information to curve the mesh, see images in spoilers below. Here is gmsh showing a straight sided mesh

Straight sided GMSH mesh ![image](https://github.com/trixi-framework/Trixi.jl/assets/7436733/bde7fec6-79ea-4e6a-ae40-a26fa186f891)

The curved airfoil already indicates that it is possible to curve the mesh. And, by setting a high order with Set order 3, I get

GMSH mesh of order 3 ![image](https://github.com/trixi-framework/Trixi.jl/assets/7436733/3a71ae0f-bf7d-4d13-a414-666c5e985167)

HOHQ generated meshes also give curved boundaries already.

Here is a gist with .inp file generated with curved boundaries - https://[gist.githubusercontent.com/Arpit-Babbar/d9238c194b3e87574fdcbbd39fa59180/raw/6c985a65b58e5eda38630b6a8d3e54652293650a/NACA6412_curved.inp](https://gist.githubusercontent.com/Arpit-Babbar/d9238c194b3e87574fdcbbd39fa59180/raw/6c985a65b58e5eda38630b6a8d3e54652293650a/NACA6412_curved.inp)

Here is a gist with .inp file generated with straight sided boundaries - https://[gist.githubusercontent.com/Arpit-Babbar/f4271dc9c57724bf62681c725aa16cd6/raw/e83e50d4f67ee1853f189d6b7269e8fb1cde5ab0/NACA6412.inp](https://gist.githubusercontent.com/Arpit-Babbar/f4271dc9c57724bf62681c725aa16cd6/raw/e83e50d4f67ee1853f189d6b7269e8fb1cde5ab0/NACA6412.inp)

andrewwinters5000 commented 9 months ago

Can you post the .inp file that gmsh generates are you curve it? It would be helpful to see what the curved ABAQUS file format looks like.

The issue is with the creation of the connectivity versus the curvature. The connectivity is handled by p4est which can only handle straight-sided elements, as discussed briefly here. The curving of the mesh is handled internally within Trixi. We just need to work out how to read this boundary point data into appropriate arrays in order to create the boundary interpolants and compute the metric terms.

Arpit-Babbar commented 9 months ago

Can you post the .inp file that gmsh generates are you curve it? It would be helpful to see what the curved ABAQUS file format looks like.

Here a gist with .inp file generated with curved boundaries - https://[gist.githubusercontent.com/Arpit-Babbar/d9238c194b3e87574fdcbbd39fa59180/raw/6c985a65b58e5eda38630b6a8d3e54652293650a/NACA6412_curved.inp](https://gist.githubusercontent.com/Arpit-Babbar/d9238c194b3e87574fdcbbd39fa59180/raw/6c985a65b58e5eda38630b6a8d3e54652293650a/NACA6412_curved.inp)

Here is a gist with .inp file generated with straight sided boundaries - https://[gist.githubusercontent.com/Arpit-Babbar/f4271dc9c57724bf62681c725aa16cd6/raw/e83e50d4f67ee1853f189d6b7269e8fb1cde5ab0/NACA6412.inp](https://gist.githubusercontent.com/Arpit-Babbar/f4271dc9c57724bf62681c725aa16cd6/raw/e83e50d4f67ee1853f189d6b7269e8fb1cde5ab0/NACA6412.inp)

The issue is with the creation of the connectivity versus the curvature. The connectivity is handled by p4est which can only handle straight-sided elements, as discussed briefly here. The curving of the mesh is handled internally within Trixi. We just need to work out how to read this boundary point data into appropriate arrays in order to create the boundary interpolants and compute the metric terms.

I see, thanks for sharing that!

andrewwinters5000 commented 9 months ago

This strategy of having p4est create and manage the connectivity and Trixi handle the curvature was the reason for the special .inp mesh file format for HOHQMesh, see here. Basically, all the curved boundary information is encoded as comments in the ABAQUS format such that the p4est parser ignores it and we can let Trixi handle everything curved boundary related.

DanielDoehring commented 9 months ago

We just need to work out how to read this boundary point data into appropriate arrays in order to create the boundary interpolants and compute the metric terms.

That would actually be a nice student project for which relatively little knowledge on math is required!