Generating an unstructured triangle grid

New to Julia
1

Hi, I’m trying to generate an unstructured triangle grid with an area constraint for data defined on irregular polygonal domains. I have coordinate data for the boundaries of two polygons respectively and for measuring points within those polygons. I now need to generate a triangular grid between the measuring points, preferably using the Triangulate.jl package, but I couldn’t find out how to do so.

I’ve found this bit of code and inserted my data, however it doesn’t actually show me a plot, even though the function runs without error.

function example_domain_bcdt_area(; maxarea = 0.05)

    triin = Triangulate.TriangulateIO()

 #I've substituted the values below with my measuring point data 
    triin.pointlist = Matrix{Cdouble}(\[0.0 0.0; 1.0 0.0; 1.0 1.0; 0.6 0.6; 0.0 1.0\]')

#I've subsituted the values below with my polygon boundaries
    triin.segmentlist = Matrix{Cint}(\[1 2; 2 3; 3 4; 4 5; 5 1\]')

    triin.segmentmarkerlist = Vector{Int32}(\[1, 2, 3, 4, 5\])

    area = sprint("%.15f", maxarea)

    (triout, vorout) = triangulate("pa$(area)DQ", triin)

    return plot_in_out(

        triin,

        triout;

        circumcircles = true,

    )

end;

I’m not sure if this is the right way to go about it, or whether I’ve inserted my data in the right places. Is there maybe a different/easier way to generate a triangular grid for my data points that will give me a plotted output? Thank you in advance for the help :slight_smile:

2

For plotting, you can use GridVisualize.jl:

You can use GLMakie, CairoMakie or PyPlot as plotting backend.

Perhaps this should be made more prominent in the package docs.

3
4

Thank you for the replies. I’ve tried to implement the code at the link @j-fu provided, and the function runs without error. However it still does not show me a plot as output. I’m using GLMakie in my script, which I’m assuming it should use as a plotting backend automatically?

I’m currently just confused as to how to get an actual plot out of running the code that looks like the example provided under the link.

I’ve also tried to use the following code with PyPlot:

injupyter()=(isdefined(Main, :IJulia) && Main.IJulia.inited)

if injupyter()
    import PyPlot
end 

function example_convex_hull_with_boundary(;Plotter=nothing, n=10, raster=10)
    triin=Triangulate.TriangulateIO()
    triin.pointlist=hcat(unique([ Cdouble[rand(1:raster)/raster, rand(1:raster)/raster] for i in 1:n])...)
    display(triin)
    (triout, vorout)=triangulate("cQ", triin)
    display(triout)
    @test numberofpoints(triin)==numberofpoints(triout)
    @test numberoftriangles(triout)>0 
    @test numberofsegments(triout)>0
    plot_in_out(Plotter, triin, triout, title="Convex hull with boundary")
end

injupyter()&&  example_convex_hull_with_boundary(Plotter=PyPlot,n=10,raster=10);

But with this I run into the same issue of not getting a plot output. Can anyone help?

5

Please open an issue in triangulate.jl.I can have a look around end of this week. Currently traveling… I didn’t use jupyter for a long time…

6

Let us continue the discussion on Discourse thread "Generating an unstructured triangle grid" · Issue #46 · JuliaGeometry/Triangulate.jl · GitHub

7

Thank you for stucking with me, I have now figured out what my issue concerning the visualization of the triangle grid was.

After running the function I failed to execute it inside a plot. I had just run the function waiting for something to pop up without then operating the function with adequate bounds, maxarea contraints, etc.

This is the code that got me where I wanted to be in the end:

#create function for triangulated mesh 
function MakeMesh(Coordinates::Matrix{Float64}, MaxArea::Float64 = 160.0) 
    triin = Triangulate.TriangulateIO()
    triin.pointlist = Matrix(Coordinates')
    N = size(PolygonCoor, 1)
    triin.segmentlist = Int32[1:N  vcat(2:N, 1)]'
    (triout, _) = Triangulate.triangulate("qa$(MaxArea)", triin)
    return triout
end

#define mesh 
mesh1 = MakeMesh(polygon1, 160.0)
mesh2 = MakeMesh(polygon2, 160.0)

#make function to triangulate mesh
function PlotMesh!(ax, mesh, fill_color, line_color)

    points = [Point2f(mesh.pointlist[1, i], mesh.pointlist[2, i]) for i in 1:size(mesh.pointlist, 2)]
    triangles = [Makie.TriangleFace(mesh.trianglelist[1, i], mesh.trianglelist[2, i], mesh.trianglelist[3, i]) for i in 1:size(mesh.trianglelist, 2)]
    mesh!(ax, points, triangles, color = fill_color)
    edges = Point2f[]
    for i in 1:size(mesh.trianglelist, 2)
        pt1 = points[mesh.trianglelist[1, i]]
        pt2 = points[mesh.trianglelist[2, i]]
        pt3 = points[mesh.trianglelist[3, i]]
        push!(edges, pt1, pt2, pt2, pt3, pt3, pt1)
    end
    linesegments!(ax, edges, color = line_color, linewidth = 0.4)
end

#make mesh plot
f = Figure(); 
ax = Axis(f[1, 1], xlabel="X[m]", ylabel="Y[m]", title="Triangle Mesh")
PlotMesh!(ax, mesh1, (:steelblue3, 0.4), :steelblue3)
PlotMesh!(ax, mesh2, (:seagreen, 0.4), :seagreen)
display(GLMakie.Screen(), f)