This page is an advanced supplement to the visual demo of Monte Carlo integration. This page requires knowledge of calculus; the other does not.
Imagine a cube with the length of each side being one unit. The volume of the cube is one unit cubed. Define x, y, and z perpendicular axes centered at the center of the cube, perpendicular to the cube's faces, parallel to various edges. Use a drill of diameter one unit to drill out the cube along the x axis. The drill is centered in two parallel cube faces and tangent to four faces. After drilling, there would be four identical pieces, just barely split from one another (but imagine we keep the pieces from separating). How much of the original volume did the drill remove? How much volume remains? Well, this is equivalent to a circle inscribed in a square. The fraction removed is pi/4 (~0.7854); the fraction remaining is 1-(pi/4), or about 0.2146.
Now we drill out the four remaining pieces along the y axis, centered in that face, tangent to four faces. Lastly we drill out the eight remaining pieces along the z axis, centered in that face, tangent to four faces. We end up with eight identical corner pieces: eight octants, eight corners. How much volume remains now? Let's call the volume remaining Vr.
We have the cube centered on the origin and aligned to the axes. This problem is mirror symmetric in each of the x, y, and z directions. The volumes of the eight corner pieces are equal to one another. Vr is equal to eight times the volume of any one corner. Let's use this to simplify the math. A cube with sides of length two has eight times the volume of the unit cube. So Vr for the unit cube drilled out by drills of unit diameter equals the volume remaining in a single corner of a cube with sides of length two drilled out by drills of unit radius. We will use the octant in which x, y, and z are non-negative: 0<=x<=1; 0<=y<=1; 0<=z<=1. Let's look at that corner.
Let's look down the line x=y=z toward the origin.
We see that the corner has three flat faces: one in the x=1 plane; one in the y=1 plane; and one in the z=1 plane. Less easy to visualize, this corner can be divided into six equivalent parts by cutting it with the following three planes: x=y; y=z; x=z. Those three planes intersect at the line x=y=z. So Vr equals six times the volume of one of those six equivalent sliver pieces of the corner (of the double size cube). Let's take the piece that has a face in the z=1 plane and includes the edge (x, 1, 1), i.e., the intersection of the y=1 plane and the z=1 plane. Let's look at the z=1 face of this sliver.
This hole was drilled perpendicular to this z=1 plane, and the the plane x=y sliced the corner perpendicular to the z=1 plane. The unseen boundary in the z direction is the plane y=z. Therefore, slicing this corner piece perpendicular to the x axis at some given x value results in an isoceles right triangle with z-depth equal to y-thickness. Let's look up the positive x axis at slices of this sliver at six different x values: x = 0.2, 0.3, 0.4, 0.5, 0.6, 0.7: six isoceles right triangles each with one corner at y=z=1.
So, for a slice of this sliver at any given x value, the area of the isoceles right triangular slice perpendicular to the x axis is always half the square of the y-thickness for that x value. The volume integral is now simply an integral in one dimension (x). Let's rename the y-thickness at an x value as h(x). Then Vr can be written as follows, where h is a function of x.
Function h(x) changes form at x = y = sqrt(2)/2, so let's separate the integral into two integrals as follows.
For 0 <= x <= sqrt(2)/2, h(x) = 1 - y(x) = 1 - sqrt(1 - x2), since the circle is x2 + y2 = 1. Substituting, we get the following.
For sqrt(2)/2 <= x <= 1, h(x) = 1 - x. Substituting, we get the following.
Expanding the square terms in the integrands we get the following.
Integrating, we get the following.
Substituting the limits of integration, we get the following.
Simplifying, we get the following.
Multiplying, we get the following, final analytic answer.
The volume remaining after drilling out the unit cube in three directions is exactly 1+sqrt(2)-(3*pi/4). Evaluating this formula to nine decimal places we get 0.058019072. Rounding to six decimal places we get 0.058019.
Go here to see a non-calculus, visual demo of the use of computer-based Monte Carlo integration for this problem. Note that this Monte Carlo simulation does not lock onto six decimal place precision even after a billion random points generated.
Friend Links
Recommended Links