Theory

In many engineering problems one needs to integrate a function in a domain. MonteCarlo integration may be an intuitive and easy solution to do so.

The idea of Monte Carlo integration of \(\int_a^b f(x)dx\) is to use the mean-value theorem from calculus, which states that the integral \(\int_a^b f(x)dx\) equals the length of the integration domain, here \(b-a\) , times the average value of \(f\), \(\bar{f}\) , in \([a,b]\). The average value can be computed by sampling f at a set of random points inside the domain and take the mean of the function values. In higher dimensions, an integral is estimated as the area/volume of the domain times the average value, and again one can evaluate the integrand at a set of random points in the domain and compute the mean value of those evaluations.

Let us introduce some quantities to help us make the specification of the integration algorithm more precise. Suppose we have some two-dimensional integral \(\int_\Omega f(x,y)dxdy\), where \(\Omega\) is a two-dimensional domain defined via a help function \(g(x,y)\): \(\Omega = {(x,y)|g(x,y) > 0}\).

That is, points \((x,y)\) for which \(g(x,y)>0\) lie inside \(\Omega\), and points for which \(g(x,y)<\Omega\) are outside \(\Omega\). The boundary of the domain ∂Ω is given by the implicit curve \(g(x,y)=0\). Such formulations of geometries have been very common during the last couple of decades, and one refers to g as a level-set function and the boundary \(g=0\) as the zero-level contour of the level-set function. For simple geometries one can easily construct \(g(x,y)\) by hand, while in more complicated industrial applications one must resort to mathematical models.

Let \(A(\Omega)\) be the area of a domain \(\Omega\). We can estimate the integral by this Monte Carlo integration method:

  1. embed the geometry \(\Omega\) in a rectangular area \(R\)

  2. draw a large number of random points \((x,y)\) in \(R\)

  3. count the fraction \(q\) of points that are inside \(\Omega\)

  4. approximate \(A(\Omega)/A(R)\) by q, i.e., set \(A(\Omega)=qA(R)\)

  5. evaluate the mean of \(f\), \(\bar{f}\) , at the points inside \(\Omega\)

  6. estimate the integral as \(A(\Omega)\bar{f}\)

Note that \(A(R)\) is trivial to compute since \(R\) is a rectangle, while \(A(\Omega)\) is unknown. However, if we assume that the fraction of \(A(R)\) occupied by \(A(\Omega)\) is the same as the fraction of random points inside \(\Omega\), we get a simple estimate for \(A(\Omega)\).