Improve makie support

docs/src/geo_and_meshes.md

@@ -11,8 +11,10 @@ using Gmsh
 gmsh.initialize()
 gmsh.option.setNumber("General.Verbosity", 2)
 gmsh.model.add("Sphere")
+# set mesh size
 gmsh.model.occ.addSphere(0,0,0,1)
 gmsh.model.occ.synchronize()
+gmsh.option.setNumber("Mesh.MeshSizeMax", 0.2)
 gmsh.model.mesh.generate(3)
 ents = gmsh.model.getEntities()
 Ω   = Inti.gmsh_import_domain(;dim=3)
@@ -24,16 +26,17 @@ gmsh.finalize()
 - You can plot a `msh` if you have a `Makie` backend (e.g. `GLMakie`)
 
 ```@example gmsh-sphere
-using CairoMakie # or GLMakie
+using CairoMakie # or GLMakie for interactivity if supported
 Γ = Inti.external_boundary(Ω)
-poly(view(msh,Γ);strokewidth=1,color=:lightgray, transparency=true)
+color = [cos(10*x[1]) for x in Inti.nodes(msh)]
+poly(view(msh,Γ);strokewidth=0.2,color, transparency=true)
 ```
 
 You can also plot the volume mesh (see the Documentation of `Makie.poly` for
 more details on possible arguments):
 
 ```@example gmsh-sphere
-poly(view(msh,Ω);strokewidth=1,color=:lightgray, transparency=true)
+poly(view(msh,Ω);color,transparency=true,strokecolor=:lightgray, alpha = 0.1)
 ```
 
 Two-dimensional meshes are very similar:
@@ -44,12 +47,17 @@ Two-dimensional meshes are very similar:
     using CairoMakie
     gmsh.initialize()
     gmsh.option.setNumber("General.Verbosity", 2)
+    gmsh.option.setNumber("Mesh.MeshSizeMax", 0.1)
     gmsh.model.add("Disk")
-    gmsh.model.occ.addDisk(0,0,0,1,1)
+    gmsh.model.occ.addDisk(0,0,0,3,1)
     gmsh.model.occ.synchronize()
     gmsh.model.mesh.generate(2)
     Ω   = Inti.gmsh_import_domain(;dim=2)
     msh = Inti.gmsh_import_mesh(Ω;dim=2)
     gmsh.finalize()
-    poly(view(msh,Ω);strokewidth=2)
+    color = [cos(20*x[1]) for x in Inti.nodes(msh)]
+    fig,ax,p = poly(view(msh,Ω);strokewidth=1,color)
+    colsize!(fig.layout, 1, Aspect(1, 3))
+    resize_to_layout!(fig)
+    fig
 ```

ext/IntiMakieExt.jl

@@ -27,96 +27,15 @@
     return coords
 end
 
-function tomakie_dim2(msh::Inti.AbstractMesh{N,T}) where {N,T}
-    coords = Makie.Point{N,T}[]
-    connec = Int[]
-    for E in Inti.element_types(msh)
-        iter = Inti.elements(msh, E)
-        D = Inti.domain(E)
-        if D isa Inti.ReferenceTriangle
-            vtxs = Inti.vertices(D)
-            for el in iter
-                for vtx in vtxs
-                    push!(coords, el(vtx))
-                    push!(connec, length(coords))
-                end
-            end
-        elseif D isa ReferenceSquare
-            # split square in two triangles for visualization
-            vtxs_down = Inti.vertices(Inti.ReferenceTriangle())
-            vtxs_up   = map(v -> -v .+ 1, vtx_down)
-            for el in iter
-                for vtxs in (vtxs_down, vtxs_up) # the two triangles
-                    for vtx in vtxs
-                        push!(coords, el(vtx))
-                        push!(connec, length(coords))
-                    end
-                end
-            end
-        end
-    end
-    return coords, connec
-end
-
-function tomakie_dim3(msh::Inti.AbstractMesh{N,T}) where {N,T}
-    coords = Makie.Point{N,T}[]
-    connec = Int[]
-    for E in Inti.element_types(msh)
-        iter = Inti.elements(msh, E)
-        D = Inti.domain(E)
-        @assert D isa Inti.ReferenceTetrahedron
-        vtxs = Inti.vertices(D)
-        for el in iter
-            for nf in 1:4 # four faces
-                for (i, vtx) in enumerate(vtxs)
-                    i == nf && continue # i-th face exclude the i-th vertex
-                    push!(coords, el(vtx))
-                    push!(connec, length(coords))
-                end
-            end
-        end
-    end
-    return coords, connec
-end
-
-function Makie.convert_arguments(P::Type{<:Makie.Lines}, msh::Inti.AbstractMesh)
+function Makie.convert_arguments(::Type{<:Makie.Lines}, msh::Inti.AbstractMesh)
     @assert Inti.geometric_dimension(msh) == 1 "Lines only supported for meshes of geometric dimension 1"
     coords = tomakie_dim1(msh)
     return (coords,)
 end
 
 function Makie.convert_arguments(P::Type{<:Makie.Poly}, msh::Inti.AbstractMesh)
-    gdim = Inti.geometric_dimension(msh)
-    if gdim == 2
-        coords, connec = tomakie_dim2(msh)
-    elseif gdim == 3
-        coords, connec = tomakie_dim3(msh)
-    else
-        error("Poly only supported for meshes of geometric dimension 2 or 3")
-    end
-    return Makie.convert_arguments(P, coords, connec)
-end
-
-function Makie.convert_arguments(
-    P::Type{<:Makie.Arrows},
-    msh::Inti.AbstractMesh{N,T},
-) where {N,T}
-    gdim = Inti.geometric_dimension(msh)
-    adim = Inti.ambient_dimension(msh)
-    codim = adim - gdim
-    @assert codim == 1 "Arrows only supported for meshes of codimension 1"
-    coords  = Makie.Point{N,T}[]
-    normals = Makie.Point{N,T}[]
-    for E in Inti.element_types(msh)
-        iter = Inti.elements(msh, E)
-        dom = Inti.domain(E)
-        xc = Inti.center(dom)
-        for el in iter
-            push!(coords, el(xc))
-            push!(normals, normal(el, xc))
-        end
-    end
-    return Makie.convert_arguments(P, coords, normals)
+    connec = Inti.triangle_connectivity(msh)
+    return Makie.convert_arguments(P, Inti.nodes(msh), connec)
 end
 
 end # module

src/domain.jl

@@ -93,3 +93,5 @@
 Base.iterate(Ω::Domain, state = 1) = iterate(entities(Ω), state)
 
 Base.isempty(Ω::Domain) = isempty(entities(Ω))
+
+Base.in(ent::AbstractEntity, Ω::Domain) = in(ent, entities(Ω))

src/mesh.jl

@@ -84,8 +84,17 @@
 elements(msh::LagrangeMesh) = msh.etype2mat
 ent2tags(msh::LagrangeMesh) = msh.ent2tags
 
+entities(msh::LagrangeMesh) = keys(msh.ent2tags)
+
+"""
+    domain(msh::LagrangeMesh)
+
+Set of all entities covered by the mesh.
 """
-    dom2elt(m::LagrangeMesh,Ω,E)
+domain(msh::LagrangeMesh) = Domain(entities(msh))
+
+"""
+    dom2elt(m::LagrangeMesh,Ω,E)::Vector{Int}
 
 Compute the element indices `idxs` of the elements of type `E` composing `Ω`, so
 that `elements(m)[idxs]` gives all the elements of type `E` meshing `Ω`.
@@ -181,6 +190,9 @@
 
 geometric_dimension(msh::SubMesh) = geometric_dimension(msh.domain)
 
+nodes(msh::SubMesh) = nodes(msh.parent)
+domain(msh::SubMesh) = msh.domain
+
 element_types(msh::SubMesh) = keys(msh.etype2etags)
 
 # ElementIterator for submesh
@@ -206,3 +218,67 @@
     state > length(iter) && (return nothing)
     return iter[state], state + 1
 end
+
+"""
+    triangle_connectivity(msh::AbstractMesh)
+
+Return a vector of integers `connec` providing a triangulation of the elements
+in the `mesh`. The `connec` vector should be read in groups of three to get the
+triangles.
+
+For elements of `LagrangeTriangle` type, the first three nodes are returned. For
+elements of `LagrangeSquare` type, the triangulation is done by splitting the
+element in two triangles. For elements of `LagrangeTetrahedron` type, the four
+faces are returned. This method will error if the element type is not supported.
+"""
+triangle_connectivity(msh::SubMesh) = _triangle_connectivity(msh.parent, msh.domain)
+triangle_connectivity(msh::LagrangeMesh) = _triangle_connectivity(msh, domain(msh))
+function _triangle_connectivity(msh::Inti.LagrangeMesh{N,T}, Ω::Inti.Domain) where {N,T}
+    connec = Int[]
+    for E in Inti.element_types(msh)
+        el_idxs = Inti.dom2elt(msh, Ω, E)::Vector{Int}
+        isempty(el_idxs) && continue
+        tags = msh.etype2mat[E]::Matrix{Int}
+        if E <: Inti.LagrangeTriangle
+            # extract the first three tags
+            for n in el_idxs
+                push!(connec, tags[1, n])
+                push!(connec, tags[2, n])
+                push!(connec, tags[3, n])
+            end
+        elseif E <: Inti.LagrangeSquare
+            for n in el_idxs
+                # lower triangle
+                push!(connec, tags[1, n])
+                push!(connec, tags[2, n])
+                push!(connec, tags[3, n])
+                # upper triangle
+                push!(connec, tags[3, n])
+                push!(connec, tags[4, n])
+                push!(connec, tags[1, n])
+            end
+        elseif E <: Inti.LagrangeTetrahedron
+            for n in el_idxs
+                # four faces
+                push!(connec, tags[1, n])
+                push!(connec, tags[2, n])
+                push!(connec, tags[3, n])
+                #
+                push!(connec, tags[1, n])
+                push!(connec, tags[2, n])
+                push!(connec, tags[4, n])
+                #
+                push!(connec, tags[1, n])
+                push!(connec, tags[3, n])
+                push!(connec, tags[4, n])
+                #
+                push!(connec, tags[2, n])
+                push!(connec, tags[3, n])
+                push!(connec, tags[4, n])
+            end
+        else
+            error("element type $E not supported")
+        end
+    end
+    return connec
+end