# Finding points in a Delaunay triangulation
**URL:** https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154
**Category:** Modelling & Simulations
**Tags:** mesh
**Created:** [November 8, 2021, 3:24pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154 "2021-11-08T15:24:02Z")
**Posts on this page:** 15
**Page:** 1
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [November 8, 2021, 3:24pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/1 "2021-11-08T15:24:02Z")
Hello everyone,
Is there a specific package/function that computes the id’s of the simplices that contain a given (probably large) set of points in a n-d Delaunay triangulation? Something equivalent to Matlab’s `tsearchn` function.
I am using the Delaunay.jl package but found no useful function.
Thank you all!
---
### Author:  [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [November 8, 2021, 3:38pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/2 "2021-11-08T15:38:49Z")
We have some search algorithms in Meshes.jl that you can take a look: [Neighbor search · Meshes.jl](https://juliageometry.github.io/Meshes.jl/stable/algorithms/neighborsearch.html)
Currently refactoring some of these, but the stable version may already have what you need?
Check the test suite for examples of usage.
---
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [November 16, 2021, 12:07pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/3 "2021-11-16T12:07:15Z")
Thanks for your answer. I went through the package but didn’t find what I was looking for. Thank you anyway!!
---
### Author:  [@jd-foster](https://discourse.julialang.org/u/jd-foster)
#### Post date: [November 22, 2021, 12:06pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/4 "2021-11-22T12:06:49Z")
This thread [[ANN] Announcing CGAL.jl](https://discourse.julialang.org/t/ann-announcing-cgal-jl/43260/1) might provide some more leads.
---
### Author:  [@mschauer](https://discourse.julialang.org/u/mschauer)
#### Post date: [November 22, 2021, 12:33pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/5 "2021-11-22T12:33:04Z")
Check `locate` in [https://github.com/JuliaGeometry/VoronoiDelaunay.jl](https://github.com/JuliaGeometry/VoronoiDelaunay.jl)
---
### Author:  [@rgcv](https://discourse.julialang.org/u/rgcv)
#### Post date: [November 24, 2021, 7:18pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/6 "2021-11-24T19:18:29Z")
Howdy!
I’m afraid CGAL.jl can’t provide a solution for OP’s general issue. The mappings I’ve exposed only pertain to 2D Delaunay Triangulations. More work could be done to make 3D ones available with some effort, but even then it seems too specific regarding what they’re looking for.
AFAIK CGAL does not have generic n-D Delaunay Triangulation algorithms, but it might just be a matter of exploring the plethora of packages deeper. I’d love to be proven wrong!
---
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [December 7, 2021, 8:06pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/7 "2021-12-07T20:06:39Z")
Hello there!
After searching, it seems that there are no available packages for this specific topic yet. I shall wait for upcoming updates!
Thank you all.
---
### Author:  [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [December 7, 2021, 9:05pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/8 "2021-12-07T21:05:04Z")
@jmaffi please provide a MWE with what you expect as the input and output. The definition of the problem is unclear still.
---
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [February 16, 2022, 11:57am UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/9 "2022-02-16T11:57:42Z")
Hello all! I’m sorry for the delayed answer… quite busy months.
I apologise if the original question was not clear enough. Hope this example will help.
Say that I have a [1 3]x[1 3] grid in 2D space:
```julia
x = [1 2 3;1 2 3;1 2 3]
3×3 Matrix{Int64}:
1 2 3
1 2 3
1 2 3
y = [1 1 1;2 2 2; 3 3 3]
3×3 Matrix{Int64}:
1 1 1
2 2 2
3 3 3
```
The complete list of points is:
```julia
9×2 Matrix{Int64}:
1 1
1 2
1 3
2 1
2 2
2 3
3 1
3 2
3 3
```
A Delaunay triangulation would have the following connectivity list:
```julia
8×3 Matrix{Int64}:
1 4 2
4 5 2
2 5 3
5 8 6
3 5 6
4 7 5
6 8 9
5 7 8
```
So, a total of 8 triangles connecting all 9 points of the grid. I plotted it in Matlab:
If I want to know which triangle holds the point say `[2.4, 1.2]`, I expect a function that would take this point as an input and provide the number `6` as a result, since the 6th entry of the connectivity list is the triangle defined by points `[4 7 5]`, which are the `[2 1]`, `[3 1]` and `[2 2]` respectively.
I am looking for such a function, not only for 2D cases, but for a general n-D triangulation. The Matlab equivalent is `tsearchn`.
I hope it’s clearer this time.
Thank you all once again.
---
### Author:  [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [February 16, 2022, 12:24pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/10 "2022-02-16T12:24:39Z")
Yes, it is much clearer now @jmaffi.
This is a solution with Meshes.jl:
```julia
using Meshes, MeshViz
import GLMakie as Mke
grid = CartesianGrid(Point(1.,1.), Point(3.,3.), dims=(2,2))
mesh = triangulate(grid) # triangles of your example
viz(mesh, showfacets = true)
```
```julia
# K-nearest search with K=1
searcher = KNearestSearch(mesh, 1)
search(Point(2.4,1.2), searcher)
1-element view(::Vector{Int64}, 1:1) with eltype Int64:
3
```
This means that the triangle 3 is the closest to the point of interest:
```julia
mesh[3]
Triangle{2,Float64}
└─Point(2.0, 1.0)
└─Point(3.0, 1.0)
└─Point(3.0, 2.0)
```
---
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [February 16, 2022, 7:27pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/11 "2022-02-16T19:27:17Z")
@juliohm thanks a lot for your answer.
Is it possible that `triangulate` works only for 2D grids? That is what I understand from the documentation. I tried creating the following 3D grid:
```julia
using Meshes, MeshViz
grid = CartesianGrid(Point(1,0,0), Point(10,5,5), dims = (10,6,6))
10×6×6 CartesianGrid{3,Float64}
minimum: Point(1.0, 0.0, 0.0)
maximum: Point(10.0, 5.0, 5.0)
spacing: (0.9, 0.8333333333333334, 0.8333333333333334)
import GLMakie as Mke
viz(grid, showfacets = true)
```
This is what I get, as expected:
However, if I try to triangulate it:
```julia
mesh = triangulate(grid)
ERROR: MethodError: no method matching triangulate(::Hexahedron{3, Float64, Vector{Point3}})
Closest candidates are:
triangulate(::Box{2}) at C:\Users\JM\.julia\packages\Meshes\mqr1d\src\discretization.jl:104
triangulate(::Triangle{Dim, T} where {Dim, T}) at C:\Users\JM\.julia\packages\Meshes\mqr1d\src\discretization.jl:106
triangulate(::Quadrangle{Dim, T} where {Dim, T}) at C:\Users\JM\.julia\packages\Meshes\mqr1d\src\discretization.jl:108
...
Stacktrace:
[1] triangulate(mesh::CartesianGrid{3, Float64})
@ Meshes C:\Users\JM\.julia\packages\Meshes\mqr1d\src\discretization.jl:128
[2] top-level scope
@ none:1
```
I don’t know if the package does not compute Delaunay triangulations for higher dimensions, or if I’m using the function in a wrong way. By the error message, it seems as if the code interprets the CartesianGrid as an Hexahedron (which is not the case).
So I would need to create a 3-D or even 4-D Delaunay triangulations and use the `KNearestSearch` function you suggested (or other).
I appreciate your help.
Regards,
---
### Author:  [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [February 16, 2022, 7:45pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/12 "2022-02-16T19:45:46Z")
@jmaffi you cannot split an hexahedron into triangles. Hexahedron are 3D geometries and triangles are 2D geometries. This is what the error message is saying.
You could split the hexahedrons into tetrahedrons but that is not implemented yet. Contributions are welcome.
---
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [February 16, 2022, 8:46pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/13 "2022-02-16T20:46:26Z")
Yes, but the grid is not a hexahedron at all, it’s a rectangular prism, as you can see in the previous image.
This 3D grid can be successfully triangulated into tetrahedra.
The result with Matlab is the following:
---
### Author:  [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [February 16, 2022, 9:13pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/14 "2022-02-16T21:13:51Z")
> [@jmaffi](#):
>
> Yes, but the grid is not a hexahedron at all, it’s a rectangular prism, as you can see in the previous image.
And who is saying that the grid is a hexahedron? The error message is clear: it is saying that in order to triangulate a grid you need to triangulate its elements, which are hexahedra. Hence, the error.
> [@jmaffi](#):
>
> This 3D grid can be successfully triangulated into tetrahedra.
That is what I said, except that it is not called a triangulation anymore, it is a “tetrahedrization”. More generally, you are after simplex decompositions of grids.
---
### Author:  [@jmaffi](https://discourse.julialang.org/u/jmaffi)
#### Post date: [February 16, 2022, 9:26pm UTC](https://discourse.julialang.org/t/finding-points-in-a-delaunay-triangulation/71154/15 "2022-02-16T21:26:45Z")
Ah, I did not get you the first time! Thank you for clarifying.
Guess I’ll be waiting for updates.
Regards,