# Integral over an irregular triangle mesh

**URL:** https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298
**Category:** New to Julia
**Created:** [October 14, 2024, 4:40pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298 "2024-10-14T16:40:15Z")
**Posts on this page:** 11
**Page:** 1

<div class="post-metadata">

### Author: ![sofiane](https://avatars.discourse-cdn.com/v4/letter/s/e19adc/32.png) [@sofiane](https://discourse.julialang.org/u/sofiane)
#### Post date: [October 14, 2024, 4:40pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/1 "2024-10-14T16:40:15Z")

</div>

Hi everyone, lately I’ve been using Julia to solve a stochastic partial differential equation for a 2D problem. I used a refined triangular mesh for my mesh generation, but I encountered some issues while implementing my log-likelihood for the Bayesian estimation of the range correlation parameter and the standard deviation. Specifically, I’m having trouble performing the integral over the triangular mesh. If anyone can help me with that, I would greatly appreciate it. Thank you!

---

<div class="post-metadata">

### Author: ![DanielVandH](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielvandh/32/31134_2.png) [@DanielVandH](https://discourse.julialang.org/u/DanielVandH)
#### Post date: [October 14, 2024, 4:42pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/2 "2024-10-14T16:42:56Z")

</div>

You will probably get faster help if you include the mathematical details / the code you have thus far

---

<div class="post-metadata">

### Author: ![sofiane](https://avatars.discourse-cdn.com/v4/letter/s/e19adc/32.png) [@sofiane](https://discourse.julialang.org/u/sofiane)
#### Post date: [October 14, 2024, 4:52pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/3 "2024-10-14T16:52:51Z")

</div>

Thank you for your reply,  
So basically i used the function Trianglemesh.jl to generate a refined triangulation over a space, now i just need to estimate the numerical integral of my values at each points (red circle) in this mesh like presented here

 ![mesh_map-1](https://global.discourse-cdn.com/julialang/original/3X/d/a/daa4332376723e0a5771f270bf687cb1adee92ea.jpeg)

---

<div class="post-metadata">

### Author: ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [October 14, 2024, 5:39pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/4 "2024-10-14T17:39:57Z")

</div>

> [@sofiane](#):
>
> now i just need to estimate the numerical integral of my values at each points (red circle) in this mesh like presented here

If you have values only at the red dots, then you’ll need to decide how you want to interpolate them in order to integrate.

One option would be to create a triangulation using _only_ the red dots ([no refinement](https://konsim83.github.io/TriangleMesh.jl/latest/man/examples/#Meshing-a-Point-Cloud-1)), and then do bilinear interpolation on each triangle. This corresponds to multiplying the area of each triangle by the average of the values at its three vertices.

(I don’t see how your “refined” triangulation will be relevant, however, unless you define some other interpolation scheme based on the refined mesh.)

Another option would be to interpolate using some other method, e.g. radial basis functions, and then to integrate the resulting smooth interpolant with a generic quadrature scheme.

It’s also not clear to me over what domain you want to integrate. If you want to integrate over the whole rectangular domain that you showed, you’ll have to decide on some extrapolation scheme for how to define the value of your integrand outside of the convex hull of your red dots.

---

<div class="post-metadata">

### Author: ![DanielVandH](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/danielvandh/32/31134_2.png) [@DanielVandH](https://discourse.julialang.org/u/DanielVandH)
#### Post date: [October 14, 2024, 5:47pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/5 "2024-10-14T17:47:26Z")

</div>

I’ll note that another good smooth interpolant would be most of the interpolants from NaturalNeighbours.jl (e.g. `Farin()`, just not `Nearest()` or probably `Triangle()`) that are defined directly using triangulations (technically it’s actually based on Voronoi tessellations)

Indeed though you might need to be a bit clearer on what you’re interpolating over / your domain if @stevengj’s approach is not clear.

---

<div class="post-metadata">

### Author: ![brainandforce](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/brainandforce/32/211054_2.png) [@brainandforce](https://discourse.julialang.org/u/brainandforce)
#### Post date: [October 14, 2024, 10:32pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/6 "2024-10-14T22:32:01Z")

</div>

Would the package [MeshIntegrals.jl](https://github.com/mikeingold/MeshIntegrals.jl) be of use here?

---

<div class="post-metadata">

### Author: ![mike.ingold](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mike.ingold/32/203749_2.png) [@mike.ingold](https://discourse.julialang.org/u/mike.ingold)
#### Post date: [October 15, 2024, 3:50pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/7 "2024-10-15T15:50:18Z")

</div>

Thanks for the suggestion 😃

It’s not clear to me exactly what @sofiane is trying to accomplish here based on the description and image, i.e. is this:

- calculate the value of a surface integral \int\_0^5\int\_0^5 f(x,y) \text{d}x\,\text{d}y where only the red points have known values?
- at each red point, formulate some integral problem and solve it?

In this first case, it sounds like most of the work will be in figuring out how to calculate/estimate values for other points in the domain, whether by some kind of interpolation scheme, direct calculation, etc. Some kind of post-hoc Monte Carlo integration scheme might also be possible, but results might not be particularly accurate, especially given how non-uniform the sampling distribution appears to be.

In any event, if you have an integral problem where the domain can be defined via a Meshes.jl geometry, and you have a field-like integrand that can be defined as a function of position within that geometry, then MeshIntegrals.jl might be helpful.

---

<div class="post-metadata">

### Author: ![jkbest2](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jkbest2/32/7350_2.png) [@jkbest2](https://discourse.julialang.org/u/jkbest2)
#### Post date: [October 15, 2024, 4:48pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/8 "2024-10-15T16:48:22Z")

</div>

Are you implementing the SPDE method to approximate a Matern process ([Lindgren et al. 2011](https://rss.onlinelibrary.wiley.com/doi/10.1111/j.1467-9868.2011.00777.x))? If so, one strategy I’ve used is to find the [Voronoi polygons](https://en.wikipedia.org/wiki/Voronoi_diagram) of the mesh and use these as weights for the integration. This is also the strategy recommended in [Simpson et al. (2012)](https://www.sciencedirect.com/science/article/pii/S2211675312000048), where the motivation for the integral is to estimate a log-Gaussian Cox process.

On a side note, have you seen [`SPDEPrecisionMatrices.jl`](https://github.com/ElOceanografo/SPDEPrecisionMatrices.jl)?

For others’ context (in case there is a better way to do this integration!), this model uses the mesh with linear elements to sparsely approximate the precision matrix of a Matern Gaussian random field. The value at each vertex of the mesh is a parameter, so the values for the entire domain will be available for the integral during estimation.

---

<div class="post-metadata">

### Author: ![sofiane](https://avatars.discourse-cdn.com/v4/letter/s/e19adc/32.png) [@sofiane](https://discourse.julialang.org/u/sofiane)
#### Post date: [October 16, 2024, 9:51am UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/9 "2024-10-16T09:51:52Z")

</div>

Thank you for your Replay and yes it’s exactly what i’am trying to do. Thank you for the paper but Can you please maybe send a piece of code that Can Do that. If you Can of course

---

<div class="post-metadata">

### Author: ![jkbest2](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jkbest2/32/7350_2.png) [@jkbest2](https://discourse.julialang.org/u/jkbest2)
#### Post date: [October 16, 2024, 11:50pm UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/10 "2024-10-16T23:50:53Z")

</div>

Unfortunately I don’t have an example. It looks TriangleMesh.jl can create a Voronoi Diagram from your mesh. The trick then is to get the area of each of the Voronoi polygons and use these as the integration weights.

---

<div class="post-metadata">

### Author: ![juliohm](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/juliohm/32/215266_2.png) [@juliohm](https://discourse.julialang.org/u/juliohm)
#### Post date: [October 17, 2024, 6:10am UTC](https://discourse.julialang.org/t/integral-over-an-irregular-triangle-mesh/121298/11 "2024-10-17T06:10:53Z")

</div>

GeoStats.jl already provides interpolation/simulation methods on top of Meshes.jl including the LindgrenProcess based on SPDEs.

Unlike the alternatives it works with all kinds of meshes, not just triangle meshes.
