Using Gmsh transfinite mehses with curved geometries

[FRUrgorri]

trafficstars

I am trying to launch the function "main" in debug mode using transfine meshes designed with Gmsh.

I have prepared two simple examples: A square-section channel and a pipe. In the first case everything seems to work ok. However, in the second case I obtain a bound error when main is defining the test space V_p (line 83 of src/Main.jl V_p = TestFESpace(Ω,reffe_p))

Perhaps it is the way in which I define the pipe mesh but I do not see the difference between the two cases other than the curvature (which naturally introduce non orthogonal cells).

You can reproduce the error with the test.jl file (you need to add GridapGmsh to the project first). I include the .msh files. I also include the .geo files I used to generete them in case they help
Test.zip

Thanks for the help

[fverdugo]

I can reproduce this error:

julia> include("test.jl")
Info    : Reading 'Channel.msh'...
Info    : 275 nodes
Info    : 432 elements
Info    : Done reading 'Channel.msh'
Info    : No current model available: creating one
Info    : Reading 'Pipe.msh'...
Info    : 2684 nodes
Info    : 3292 elements
Info    : Done reading 'Pipe.msh'
[1.949e+00 s in MAIN] fe_spaces
[1.333e+00 s in MAIN] solve
ERROR: LoadError: BoundsError: attempt to access 1-element Vector{Vector{Int64}} at index [2]

[FRUrgorri]

Yes, that is the error.

The Cartesian geometry (Channel.msh) works but the cylindrical geometry (Pipe.msh) gives the error.

[fverdugo]

An even smaller reproducer:

using Gridap
using GridapGmsh
model = GmshDiscreteModel("Pipe.msh")
writevtk(model,"run_model")
k = 2
reffe = ReferenceFE(lagrangian,Float64,k-1;space=:P)
V = TestFESpace(model,reffe) # Error

the error can be circumvented using

V = TestFESpace(model,reffe,conformity=:L2) # Works

I'll fix it in the main function

[fverdugo]

@FRUrgorri fixed in master

[FRUrgorri]

Thanks @fverdugo. I think that a similar solution needs to be applied to V_j as well since the same error is produced when Main.jl tries to define this space.

[fverdugo]

Raviart-Thomas elements seem to be implemented in Gridap for "oriented" meshes at this moment. Cartesian hex meshes and unstructured tet meshes are oriented, but unstructured hex meshes are potentially not.

It is not a minor fix, but the extension to non-oriented meshes should be feasible.

[FRUrgorri]

It seems that the problem can be avoided by not recombining the transfinite meshes in Gmsh, i.e. using prism meshes. I will do more testing but so far it has worked