Basic mesh examples

[mancellin]

Inti's mesh type looks very powerful, but seems a bit complicated for my basic needs.

Could you provide some basic examples of just defining a 2D boundary from a list of point, or a 3D surface from a list of triangles of shape nb_faces × 3 × 3?

[maltezfaria]

At the moment you can't really create a mesh by simply passing e.g. points and connectivities, but that should be easy to add. One issue I see is that it would not support many of the operations currently supported by our mesh type, but you would be able to create a quadrature and solve an integral equation on it, so it may be enough 👍

[maltezfaria]

But what about using the gmsh API? You can just import meshes from there:

using Inti
using Gmsh
filename = joinpath(Inti.PROJECT_ROOT,"docs", "assets", "piece.msh")
msh = Inti.import_mesh(filename)

[mancellin]

One issue I see is that it would not support many of the operations currently supported by our mesh type

How many methods would I need to implement to have my own DumbMesh type and pass it to Quadrature and single_double_layer()?

But what about using the gmsh API?

That's a possible workaround, but it would be nice to test my hydrodynamical meshes without a detour by Gmsh.

[tanderson92]

Doesn't this issue expand the scope of Inti beyond integral equations? The docs are explicit that we view meshing as outside our scope.

[maltezfaria]

Doesn't this issue expand the scope of Inti beyond integral equations? The docs are explicit that we view meshing as outside our scope.

I think this is well within the scope of stuff we can/should cover.

[maltezfaria]

How many methods would I need to implement to have my own DumbMesh type and pass it to Quadrature and single_double_layer()?

You would need to implement element_types(msh) to return a list of element types of your mesh, and elements(msh,E) to return an iterable collection of elements of a given type. For simple meshes containing only one type, this is an overkill, but needed since the internals support meshes with elements of various types.

But what about using the gmsh API?

That's a possible workaround, but it would be nice to test my hydrodynamical meshes without a detour by Gmsh.

It would be easy to create a constructor like simple_mesh(nodes::Vector{SVector{N,T}}, connectivity::Matrix{Int}) returning a mesh of a single element type (but still using our internal mesh format). You probably would still need to pass an element type, since e.g. nodes in 3d and a 4 x n connectivity matrix could mean a surface mesh of quads or a volume mesh of tetras.

But even simpler: is there a mesh format you have in mind? FWIW, this could be as easy as extending import_mesh to cover your type. If you have a mesh format and example file, I could look into that.

[mancellin]

You probably would still need to pass an element type, since e.g. nodes in 3d and a 4 x n connectivity matrix could mean a surface mesh of quads or a volume mesh of tetras.

Call the function simple_surface_mesh(...) and there is no ambiguity anymore?

But even simpler: is there a mesh format you have in mind? FWIW, this could be as easy as extending import_mesh to cover your type. If you have a mesh format and example file, I could look into that.

I have a couple of STL files, so you can start by providing some interface with https://github.com/JuliaIO/MeshIO.jl
The more niche hydro qmesh format are basically a list of quadrilateral faces under the form Tuple{Vector{SVector{3, Float64}}, Vector{SVector{4, Int}}} or Vector{SVector{4, SVector{3, Float64}}}}. I'll need to write some IO for another project, so we can do some interfacing then.

[maltezfaria]

You probably would still need to pass an element type, since e.g. nodes in 3d and a 4 x n connectivity matrix could mean a surface mesh of quads or a volume mesh of tetras.

Call the function simple_surface_mesh(...) and there is no ambiguity anymore?

Sure, that is possible. But I am not sure how useful it would be. Do people really create a mesh by hand and pass the points? It is more likely the points are already in some file right?

But even simpler: is there a mesh format you have in mind? FWIW, this could be as easy as extending import_mesh to cover your type. If you have a mesh format and example file, I could look into that.

I have a couple of STL files, so you can start by providing some interface with https://github.com/JuliaIO/MeshIO.jl The more niche hydro qmesh format are basically a list of quadrilateral faces under the form Tuple{Vector{SVector{3, Float64}}, Vector{SVector{4, Int}}} or Vector{SVector{4, SVector{3, Float64}}}}. I'll need to write some IO for another project, so we can do some interfacing then.

I think STL files are already supported. I just tried the following on a simple cube.stl mesh file:

using Inti
using Gmsh
msh = Inti.import_mesh("cube.stl")

and it worked. As long as gmsh.open(file) works, Inti should be able to handle the rest.

[mancellin]

Sure, that is possible. But I am not sure how useful it would be. Do people really create a mesh by hand and pass the points? It is more likely the points are already in some file right?

For basic meshes, users (like myself) might want to use their own loaders or mesh-tweaking-function.

I would be nice to have either a simple API or (the initial question of this thread) an example of how to handle Dict{DataType, AbstractVector} with 1 entry: LagrangeElement{ReferenceSimplex{2}, 3, SVector{3, Float64}} => LagrangeElement{ReferenceSimplex{2}, 3, SVector{3, Float64}}[LagrangeElement{ReferenceSimpl….

I think STL files are already supported.

👍 I keep underestimating Gmsh.

[maltezfaria]

Fair enough, I guess I see your use case now...

[maltezfaria]

This should work:

function quadrilateral_surface_mesh(
    nodes::Vector{SVector{3,Float64}},
    connectivity::Matrix{Int},
)
    @assert size(connectivity, 1) == 4 "connectivity matrix must have 4 rows"
    E = Inti.LagrangeSquare{4,SVector{3,Float64}}
    etype2mat = Dict(E => connectivity)
    etype2els = Dict{DataType,AbstractVector}()
    tag = Inti.new_tag(2)
    ent = Inti.GeometricEntity(2, tag, EntityKey[], String[], nothing)
    etags = collect(1:size(connectivity, 2))
    ent2etag = Dict(ent => Dict(E => etags))
    msh = Inti.Mesh{3,Float64}(nodes, etype2mat, etype2els, ent2etag)
    msh.etype2els[E] = Inti.ElementIterator(msh, E)
    return msh
end

[maltezfaria]

With 0ba4726, the internal ordering nodes has been changed to facilitate the implementation of generic $\mathbb{P}_k$ and $\mathbb{Q}_k$ elements. As such, the script below will need to be modified in the next release.

I should add some docs explaining the new ordering for LagrangeElements...

This should work:

function quadrilateral_surface_mesh(
nodes::Vector{SVector{3,Float64}},
connectivity::Matrix{Int},
)
@Assert size(connectivity, 1) == 4 "connectivity matrix must have 4 rows"
E = Inti.LagrangeSquare{4,SVector{3,Float64}}
etype2mat = Dict(E => connectivity)
etype2els = Dict{DataType,AbstractVector}()
tag = Inti.new_tag(2)
ent = Inti.GeometricEntity(2, tag, EntityKey[], String[], nothing)
etags = collect(1:size(connectivity, 2))
ent2etag = Dict(ent => Dict(E => etags))
msh = Inti.Mesh{3,Float64}(nodes, etype2mat, etype2els, ent2etag)
msh.etype2els[E] = Inti.ElementIterator(msh, E)
return msh
end