Like 0Likes Dislike Area and Volume Calculations

equation int #8729 volume area set flux part closed
Presents a three dimensional volume calculator based on the pressure volume model.

Motivation

In recent months, the indie game programming scene has witnessed the rising popularity of a pressure volume model described in a paper written by Maciej Matyka and Mark Ollila. The paper can be
found here: http://panoramix.ift.uni.wroc.pl/~maq/eng/index.php. In the model's first presentation the authors used an axis aligned
bounding box to approximate the volume of a body. Such a rough approximation is easily computed but highly inaccurate. A second paper followed sometime later describing a two dimensional
implementation that employed Gauss's Theorem to calculate the exact area of a simulated body. The author mentioned the possibility of a three dimensional version using the same technique but
restricted his presentation to the two dimensional case. It would seem that a three dimensional volume calculator should have surfaced shortly afterward; however, no such algorithm has been detailed

Problem Space

The two dimensional case: Given a set of edges that form a closed loop calculate the enclosed area.

The three dimensional case: Given a set of triangles that form a closed surface calculate the enclosed volume.

Derivation (2D)

Gauss's Theorem relates the divergence of a vector field within an area to the flux of a vector field through a closed path by the following:

∫∫a dell ∙F da = ∫l F ∙ dn
equation 1.1

where the path l encloses the area a. The dell vector is defined by:

dell = (∂/∂x, ∂/∂y)
equation 1.2

When F is set to:

F(x, y) = (x, 0)
equation 1.3

the innards of the double integral in equation 1.1 can be reduced to one. This reduction lets us write:

∫∫a 1 da = ∫l F ∙ dn
equation 1.4

which equates to the area of a because:

dell ∙ F = x ∂/∂x + 0 ∂/∂y = ∂x/∂x = 1
equation 1.5

Equation 1.4 tells us that the flux of the vector field F through the closed path l is equal to the area enclosed by l. When l is defined by a set of edges the enclosed
area is equal to the sum of the flux of F through every edge. A straight forward parameterization of the edge {v1, v2} is given by:

e = v2v1
equation 1.6

l(t) = v1 + te
equation 1.7

where the evaluation of l(t) yields any point on the edge for t on [0, 1]. Under this parameterization dn from equation 2.4 can be written as:

dn = cross dl = cross e dt = (-ey, ex)dt
equation 1.8

The flux of F through the edge {v1, v2} now follows:

Φ = ∫ F(l(t)) (-ey, ex)dt
equation 1.9

Φ = ∫ -ey(v1x + tex)dt
equation 1.10

where the integral is evaluated through t on [0, 1]. The analytical solution to equation 1.10 can be written as:

Φ = ½ ∙ (v1yv2y) ∙ (v1x +v2x)
equation 1.11

The right hand side of equation 1.11 yields the flux of F through the edge {v1, v2}. Given a closed and clockwise wound edge set {vi1,
vi2} for i = {0, 1, 2,…,n} the enclosed area is computed by evaluating the following sum:

area = ½ ∑ (vi1yvi2y) ∙ (vi1x +vi2x)
equation 1.12

C++ Implementation (2D)

Listing 1: Area Calculation:


<span class="codecomment">// v: A pointer to the array of vertices

// i: A pointer to the array of indices

// n: The number of indices (multiple of 2)

// This function uses Gauss's Theorem to calculate the area (or field)

// of a body enclosed by a set of edges. The edge set must form a

// closed path in R2 and be outward facing. An outward facing edge set

// wraps clockwise around the body where the x-axis points left and the

// y-axis points up.</span>

<span class="codekeyword">float</span> CalculateArea(<span class="codekeyword">const</span> VECTOR2      *v,

<span class="codekeyword">const unsigned int</span> *i,

<span class="codekeyword">const unsigned int</span>  n)

{

<span class="codekeyword">unsigned int</span> j;

VECTOR2      v1;

VECTOR2      v2;

<span class="codekeyword">float</span>        area = 0.0f;

<span class="codekeyword">for</span>(j = 0; j < n; j+=2)

{

v1 = v[i[j ]];

v2 = v[i[j+1]];

area += (v1.y-v2.y) * (v1.x+v2.x);

}

<span class="codekeyword">return</span> area / 2.0f;

}



Derivation (3D)

Gauss's Theorem relates the divergence of a vector field within a volume to the flux of a vector field through a closed surface by the following:

∫∫∫v dell ∙ F dv = ∫∫s F ∙ da
equation 2.1

where the surface s encloses the volume v. The dell vector is defined by:

dell = (∂/∂x, ∂/∂y, ∂/∂z)
equation 2.2

When F is set to:

F(x, y, z) = (x, 0, 0)
equation 2.3

the innards of the triple integral in equation 1.1 can be reduced to one. This reduction lets us write:

∫∫∫v 1 dv = ∫∫s F ∙ da
equation 2.4

which equates to the volume of v because:

dell ∙ F = x ∂/∂x + 0 ∂/∂y + 0 ∂/∂z = ∂x/∂x = 1
equation 2.5

Equation 2.4 tells us that the flux of the vector field F through the closed surface s is equal to the volume enclosed by s. When s is defined by a set of triangles the
enclosed volume is equal to the sum of the flux of F through every triangle. A straightforward parameterization of the triangle {v1, v2,
v3} is given by:

e1 = v2v1
equation 2.6

e2 = v3v1
equation 2.7

s(u, v)= v1 + ue1 + ve2
equation 2.8

where the evaluation of s(u, v) yields any point on the triangle for u on [0, 1] and v on [0, 1-u]. Under this parameterization da from equation 2.4 can be written as:

dn = (e1 cross e2)dvdu
equation 2.9

The flux of F through the triangle {v1, v2, v3} now follows:

Φ = ∫∫F(s(u, v)) (e1 cross e2)dvdu
equation 2.10

Φ = ∫∫ (v1x + ue1x + ve2x)(e1ye2ze1ze2y)dvdu
equation 2.11

where the integral is evaluated through u on [0, 1] and v on [0, 1-u]. The analytical solution to equation 2.11 can be written as:

Φ = 1/6 ((v2yv1y)(v3zv1z) –
(v2zv1z)(v3yv1y))(v1x +
v2x + v3x)
equation 2.12

The right hand side of equation 2.12 yields the flux of F through the triangle {v1, v2, v3}. Given a closed and clockwise wound
triangle set {vi1, vi2, vi3} for i = {0, 1, 2,…,n} the enclosed volume is computed by evaluating the following sum:

volume = 1/6 ∑ ((vi2yvi1y)(vi3zvi1z) –
(vi2zvi1z)(vi3yvi1y))(vi1x +
vi2x + vi3x)

equation 2.13

C++ Implementation (3D)

Listing 2: Volume Calculation


<span class="codecomment">// v: A pointer to the array of vertices

// i: A pointer to the array of indices

// n: The number of indices (multiple of 3)

// This function uses Gauss's Theorem to calculate the volume of a body

// enclosed by a set of triangles. The triangle set must form a closed

// surface in R3 and be outward facing. Outward facing triangles index

// their vertices in a counterclockwise order where the x-axis points

// left, the y-axis point up and the z-axis points toward you (rhs).</span>

<span class="codekeyword">float</span> CalculateVolume(<span class="codekeyword">const</span> VECTOR3      *v,

<span class="codekeyword">const unsigned int</span> *i,

<span class="codekeyword">const unsigned int</span>  n)

{

<span class="codekeyword">unsigned int</span> j;

VECTOR3      v1;

VECTOR3      v2;

VECTOR3      v3;

<span class="codekeyword">float</span>        volume = 0.0f;

<span class="codekeyword">for</span>(j = 0; j < n; j+=3)

{

v1 = v[i[j ]];

v2 = v[i[j+1]];

v3 = v[i[j+2]];

volume += ((v2.y-v1.y)*(v3.z-v1.z)–

(v2.z-v1.z)*(v3.y-v1.y) )*(v1.x+v2.x+v3.x);

}

<span class="codekeyword">return</span> volume / 6.0f;

}



Conclusion

And there you have it, a fast and precise way to find the area or volume of a body given its edge or triangle set. I have tested the code on a variety of basic primitives and am 99% certain it is
bug free. You can contact me with question, comments and bug reports at kwoxrf@umr.edu.