In 3D you can triangulate the cloud of points - here you can start with any tetrahedron (take one where the points are far apart), and then continue to add more points. If the new point is within a tetrahedron you have already constructed (can be checked by barycentric coordinates), the point is ignored. Otherwise, a new tetrahedron is added to the triangulation (you have to make sure that it does not intersect an already added tetrahedron - if it completely covers a generated tetrahedron, the covered one is removed).

After this procedure which is linear in the number of points and roughly O(nT) where T is the number of tetrahedrons needed (quick guess) you are left with a triangulation of the polyhedron whose 3D-volume can be easily computed as sum of the volumes of the tetrahedra (simple formula!). This is the O(T).

Nth root of the volume of the convex enveloppe?!

Can we compute such a volume for an arbitrary irregular polyhedron?

Are you working in 3-D or n-D?

I found this

http://www.jstor.org/stable/1967277?seq=3

a bit old probably

Thanks Patrick. Yes, I am working on 3D. I think that the desired formula is inside the pdf you mentioned, thanks a lot. However, it seems that is it computationally very cost since we have to compute a determinant of order N. I can't understand how do some people think when they down-vote the questions / answer of others!!

In 3D you can triangulate the cloud of points - here you can start with any tetrahedron (take one where the points are far apart), and then continue to add more points. If the new point is within a tetrahedron you have already constructed (can be checked by barycentric coordinates), the point is ignored. Otherwise, a new tetrahedron is added to the triangulation (you have to make sure that it does not intersect an already added tetrahedron - if it completely covers a generated tetrahedron, the covered one is removed).

After this procedure which is linear in the number of points and roughly O(nT) where T is the number of tetrahedrons needed (quick guess) you are left with a triangulation of the polyhedron whose 3D-volume can be easily computed as sum of the volumes of the tetrahedra (simple formula!). This is the O(T).

@Hermann:

After generating a triangulation, the volume can be computed in O(m), where m

Thanks Christopher for your powerful answer. I am not seeing perfectly what you said since it is too technical but the purpose of my study is to characterize a bunch of fluid particles (passive tracers) emitted from point sources in thousands (N here is equal to 2000). The particles are initially very close to each other (at the smallest scale of turbulence) and then under the effect of the fluid motion, they are advected and the distance increases dramatically between each 2 particles. I would like to find a way to characterize this bunch as it evolves in time. I was first interested in the eigenvalues of the moment of inertia tensor and now I would like to see other geometrical properties of the emitted bunch. I think that I am more interested in the overall volume.

Thanks Christopher. Yes, I think that CGAL would be a good option since I am working with C. Best, a

Dear Abdallah,

There are N vertices of your polyhedron. In other words, you are dealing in a situation in (N - 1) dimensions. Have I understood your situation correctly?

Hemanta. The N points are in 3D. I wanted to characterize this bunch in order to do a range of things with it. Thanks!

Dear Abdallah,

In that case, it is a polygon, and not a polyhedron. A polyhedron (known also as a simplex) in an (N -1) dimensional space, has N vertices.

Dear Hemanta,

since the points are in 3D it is a polyhedron. A polygon is in 2D. A simplex in 3D is a tetrahedron. What you probably mean is a polytope, which is the convex hull of finitely many points in nD.

Abdallah, considering your application, as you have explained it in your last mail, and thinking about it I am roughly guessing that you would like to use these N particles to characterize how much a fluid (e.g. a contaminant) is being dispersed within another fluid. The tracers are then some sort of Monte-Carlo-like method to get an idea of what happens.

If this is true, characterizing the polyhedral volume of the particle distribution generated by these tracers is probably not really what you want. Given that the overall behaviour of the turbulent fluid is chaotic, making several runs with different initial distribution of the tracers might lead to significantly different results. Also, it is likely that the size of the polyhedron significantly underestimates the possible spread of the particles emanating from the local source that you are considering. This is since it is rather unlikely that by pure sampling in one run you get more than one of the "extremal particles", i.e., the ones drifting farthest from the barycenter of all particles. Due to this, the statistical properties of the volume of the polyhedron over many runs will likely not resemble the statistical properties of the dispersion of particles in the fluid.

Just a rough fantasizing on my part: You could subdivide your 2000 initial tracers in packets of about 50, putting those in the same set that are closest in the initial configuration at t=0. At t=T those 50 points will have new positions in R³. Then you estimate a covariance model using these 50 points and compute its potential 1/2 (x−μ)C^(-1)(x-μ) + 1/2 log det C, where C is the covariance matrix of the model, and μ is the mean (in your case the barycenter). You sum up all the 40 potentials to get a combined potential V. Then the level sets of this potential provide some model for the probability of a fluid particle ending up within this level set. The volumes of these level sets might be interesting numbers to consider.

Dear Hermann. Thanks a lot for your powerful suggestion. Yes you are right with concern to the extreme events that are present in the dispersion process of these particles: particles which separate too much faster than the mean squared separation and particles that experience very slow motion. This is indeed linked to regions where the fluid velocity gradient is too large or too small. I really appreciate the idea of computing the potential over these 40 samples. Just to be sure, you mean by the covariance matrix C the moment of inertia tensor or it is an other quantity? I will then implement it to get new insights. Thanks!

Dear Abdallah.

What the best choice for C is, is something that has to be tested, and it depends on the amount of nonlinearity that is present in your application. There are several possibilities depending on the geometry of the particle distribution at the end of the simulation. If you take the C to be constant, then the sum of the potentials will be a quadratic function, and the level sets will be ellipsiods, which is likely not what you want. So I thought that C should depend on x. However, something else came to my mind, which might make implementation easier. You could use as potentials a function f((x−μ)^TC^(-1)(x-μ) + 1/2 log det C) where f:R->R with f(0)=0 and f'(0)=1, e.g., f(x)=log(1+x). In that case C would be the constant 3x3 matrix given by the formula $1/N \sum_{i=1}^N x_i x_i^T - μμ^T$ where N=50, x_i are the positions of the 50 tracers of one group in R^3, and μ is the barycenter of the x_i, i.e., $1/N \sum_{i=1}^N x_i$.

But I told you, this is just guessing... You might need to play around with the function f.

Then there is an alternative, which is to use a bigger matrix, e.g. the 6x6 matrix C which is given by $1/N \sum_{i=1}^N y_i y_i^T - μμ^T$ and here y_i is the 6-dimensional vector consisting of the position x_i of the tracers and their squares x_i.^2. The μ would then be the average of the y_i. The potential would then depend on the y_i, but since the y_i depend nonlinearly on the x_i this potential captures more of the non-linearities of the flow of the fluid (level sets of one term of the potential look then more like bananas as opposed to ellipses for the above mentioned potential).

You can put even more nonlinear functions of x into y but then you probably need to increase the number of tracers per group (about 5/2*s(s+1), where s is the number of components of y) in order to get a useful matrix C.

Thanks Hermann too much for your useful clarifications. I have followed the method you described in order to compute the potential for each of the groups composed by the 50 particles. What I found is that all the curves nicely collapse to a single one (the plot in attach there is the mean potential for the 2,000 particles versus time).

The first term of the potential (which contain the inverse of C) is found to be constant in time. The potential is then uniquely controlled by the logarithmic term. I have just used the identity function for f (as defined in one of your previous answers). I am afraid to use the other definition for f because the potential could be too negative.

For times in the dissipation range, i.e. for the smallest scales, the potential is negative. The transition to the positive potentials occurs in the intermediate scales up to the integral scale. I think that it is time now to use the other formulation you suggested in order to include the nonlinearities. Thanks again for your help! best

a

Dear Abdallah,

the first term of the potential should be constant if you plug the x values of the tracer particles used for constructing C back into the potential function. Then the log det C is giving information on the size of the ellipsoids describing the level sets. If this term increases with time that tells you that the tracer particles cover an ever increasing area in the fluid, and if all the curves look similar, then the particles tend to spread somewhat similarly, no matter where you start.

However, I made a small mistake, meaning f((x−μ)^TC^(-1)(x-μ)) + 1/2 log det C (the closing parenthesis was in the wrong place). Then the f term cannot become negative since C is a positive (semi)definite matrix, making the argument of f nonnegative.

The point with the potential is that if you look at the level sets of the function V(x) = \sum_{i=1}^50 f((x-μ_i)^TC_i^(-1)(x-μ_i)) + 1/2 log det C_i (with a free variable x) that might provide some interesting information about the probable distribution of the particles (you can get it with contour3 in Matlab).

Thanks Hermann too much for your clarifications and ideas. Indeed, you helped me to make a good progress in my project. Thanks!

a

Uhm, maybe you need (from time to time) to estimate the volume directly, and the rest of the time ticks, you can estimate its hmm, dispersion and volume change indirectly.

But about the direct measurement: My suggestion is to use an octree with a cheap predicate, which tells which of the points p are "boundary points". Let's call this predicate B(p). This does not mean the same as convex closure but simply that there is an "open view" outwards along an axis from the point. Open view is a spatial cone defined with one parameter, the cone angle. (cone angle == pi ---> approximates the bounding box).

Predicate B(p) can be made rather cheap, O(n^(1/3) log(n)) easily and O(log(n)) with some effort.

The construction of the octree can be tweaked in that way that the an "inner-point-to-be" (permanently not-B(p)) is deleted from the tree. That way the tree stays rather small. (I assume you will not have internal cavities which should be taken into account.)

In the end you simply sum up the volume of the tree. Boundary points can have a local adaptation to make the mapping smoother. The required neighborhood search is rather cheap, and can be made cheaper than presented in classical papers by adding some extra linking layer. But if you do some sort of numerical simulation, your point set evolves and you could smooth the results afterwards?

If cavities are possibility, then there could be a volume summation using e.g. three-stage process, where you:

- approximate the local point density and perhaps check the local dimensionality,

- estimate an effective radius for volume summation. Properties of this estimation could be learned by Machine Learning style if you can provide a set of typical samples

- sum up V=a*\sum_i effectiveRadius_i^3

This approach stands and falls upon how similar your cases are from run to run.

This is just a set of wild areas. Most of this can be done, given some money and time, but dunno how relevant it is for your particular problem. Anyways, this is very interesting topic.