Estimating pi - Iowa State University



Estimating [pic] through Monte Carlo Methods in n-dimensions

Brian Reece

Iowa State University

MSM Creative Component

June 30, 2006

Heather Thompson, Co-Major Professor

Irvin R. Hentzel, Co-Major Professor

Alejandro Andreotti, Committee Member

Introduction

The number [pic] has intrigued mathematicians for more than 4000 years. In the pursuit of[pic], mathematicians have used geometry, trigonometry, algebra, calculus and probability. This interaction between different branches of mathematics is one of things that makes [pic] interesting to study. This paper begins by looking at different ways of calculating[pic] that have been used throughout history. Next, it looks at the Monte Carlo method and its applications. Next, a method of estimating [pic] with Monte Carlo methods in n-dimensional space is presented. The paper ends with a discussion of how the topics can be incorporated into the classroom.

The number [pic]has not always had a name. It wasn’t until 1706, when William Jones wrote[pic], that the world had a symbol for the ratio of a circle’s circumference to its diameter. Euler adopted the symbol in 1737 and it became the standard symbol.

The ratio of a circle’s circumference to its diameter has long been understood as a constant.

Then he made the molten sea; it was round, ten cubits from brim to brim, and five cubits high, and a line of thirty cubits measure its circumference (I Kings 7:23)

This passage from the Bible places that constant at[pic]. This estimate was most likely achieved through measurement, so it is quite reasonable given the likely precision of measurement. An Egyptian document called the Rhind papyrus, which dates to 1650 B.C., gives a value of[pic], or 3.1605. It was written by a scribe named Ahmes, who copied it from a document that is 200 years older. Some even believe that the Rhind papyrus can trace its origin back to 3400 B.C.

Sets of numbers

To understand what makes [pic] special, one should look at where it fits into the realm of all numbers. The number [pic] is a transcendental number. The transcendental numbers are so called because they transcend the algebraic numbers.

Algebraic numbers are roots of a polynomial equation of the form[pic], where [pic] are integers.

For example,[pic], an irrational number, is a solution to the quadratic equation[pic]. Algebraic numbers also include complex numbers. Using the simplest example, [pic]is a solution of the polynomial equation[pic].

Transcendental numbers are not the root of any integer polynomial equation. In 1844, Joseph Liouville proved that transcendental numbers exist by finding examples. Here is an example of one such number, called Liouville’s constant in honor of its discoverer:

[pic]

This number contains 0’s in every place except for the n! decimal places. So the 1st, 2nd, 6th, 24th, … decimal places contain 1’s. The next 1 will occur in the 5! = 120th decimal place. In 1850, Liouville proved that this number cannot be a solution of a polynomial equation with integer coefficients and thus is a transcendental number.

The number [pic] was proven to be transcendental in 1873 by Charles Hermite. Ferdinand Lindemann used Hermite’s results to prove that [pic] was transcendental in 1882. This proof by Lindemann also had a powerful byproduct. Since[pic] is transcendental, it follows that [pic] cannot be the solution to a polynomial equation with integer coefficients; therefore, one of the earliest math mysteries was proven an impossibility.

Since antiquity, mathematicians had tried to construct a square of area equal to the area of a given circle using only a straightedge and compass. Lindemann’s discovery that [pic] is not the root of any polynomial equation proves this exercise impossible. In order to do this, one must construct[pic], which is only a possibility if [pic] is an algebraic number. With Lindemann proving that [pic] is a transcendental number, he showed that the problem of squaring a circle was unsolvable in Euclidean geometry.

Archimedes’ Method of Trapping [pic]

Archimedes was the first to develop a formal approximation of[pic]. He constructed two series of polygons, one that inscribed a given circle and one that superscribed the given circle. Using these polygons he created an interval for[pic].

Archimedes first showed that the area of a regular polygon is[pic]. An n-gon would have exactly n congruent triangles drawn from the circumcenter of the polygon to each vertex. Figure 1 shows the first three regular n-gons divided into n congruent triangles.

[pic][pic][pic]

Figure 1: Regular n-gons divided into n congruent triangles.

The area of a triangle is, of course[pic]. When this formula is applied to one of the congruent triangles, the area of the triangle is[pic], where [pic]is a side length and [pic]is the apothem. The apothem is the perpendicular distance from the center to one side of the polygon. The area of the entire n-gon can be found then by multiplying this area by[pic], resulting in

[pic].

The perimeter of a regular n-gon is just[pic], so substituting results in

[pic].

Figure 2 helps illustrate that as [pic],[pic]and[pic], and therefore

[pic].

This result was Archimedes’ formula for the area of a circle.

[pic] [pic]… [pic]

Figure 2: Regular triangle, square and 30-gon with apothem and radius of the circumscribed circle marked.

From Euclid, Archimedes knew that the ratio of a circle’s area to the square of its diameter was a constant. Let [pic] be that constant. Then

[pic].

Substituting the previous result for the area of a circle, one gets

[pic][pic],

which results in [pic](Linn and Neal, March 2006).

This constant [pic] is of course now known as[pic].

Consider the unit circle. The diameter would be 2 and [pic]. Using construction, one can inscribe a regular polygon of [pic]sides, with semiperimeter [pic], and superscribe a regular polygon of [pic]sides, with semiperimeter [pic]. By constructing these regular polygons, one obtains an increasing sequence [pic]and a decreasing sequence [pic]such that as [pic] and[pic]. This traps [pic] between two numbers,[pic].

Starting with[pic], one has a regular polygon with 6 sides, or a hexagon. Choosing a point A on the circle and copying the radius, one can find two more points of the inscribed hexagon by drawing an arc centered at A. Figure 3 shows the results of this construction, points B and D. Repeating this process with each resulting point as a center for the arc results in the inscribed hexagon.

[pic]

Figure 3: Construction of vertices of inscribed hexagon.

To get the superscribed hexagon, a line is constructed perpendicular to the radius at each of the vertices of the inscribed hexagon. Where the perpendicular lines intersect are the vertices of the superscribed hexagon. Figure 4 illustrates the construction of one such vertex. Point H is a vertex of the superscribed hexagon.

[pic]

Figure 4: Construction of a vertex of the superscribed hexagon.

Figure 5 shows the circle with the inscribed and superscribed hexagons.

[pic]

Figure 5: Archimedes’s first attempt at trapping[pic].

A side length in the inscribed hexagon would be 1. So the perimeter is 6, and[pic]. Figure 6 shows the length of one half-side of the inscribed hexagon expressed in terms of a trigonometric ratio. Figure 7 shows a similar result for the superscribed hexagon.

[pic]

Figure 6: Inscribed regular hexagon with radius and half-side length represented.

[pic]

Figure 7: Superscribed regular hexagon with radius and half-side length represented.

If [pic]is increased to 2, a dodecagon must be constructed. First, one has to draw a ray from the center, O, to a vertex of the superscribed hexagon. The point of intersection between this ray and the circle is a point on the inscribed dodecagon. Drawing 6 of these rays results in 6 new points, or 12 in total. These 12 points are the vertices of the inscribed dodecagon. Figure 8 illustrates the construction of one of these vertices.

[pic]

Figure 8: Construction of a new vertex of the inscribed regular dodecagon.

Figure 9 illustrates the construction of a vertex of the superscribed dodecagon. First the midpoint of [pic] and the midpoint of[pic]are constructed. Then two rays are drawn from the center of the circle, O, through points B and D. The point of intersection of [pic] and [pic]is point A and the intersection of[pic]and [pic]is point E. Next a line perpendicular to [pic] is constructed through point A; the same is done for [pic] through point E. The intersection of these two lines is point I, a vertex on the superscribed dodecagon.

[pic]

Figure 9: Construction of a vertex on the superscribed regular dodecagon.

This process can be repeated from each [pic]-gon to gain vertices of the [pic]-gons.

In general terms [pic] and [pic]can be expressed as

[pic]

[pic]

Using these expressions for [pic] and[pic], the first five terms of each sequence can be found.

[pic]

[pic]

What makes Archimedes’ work incredible is that he did not have the use of trigonometry or algebra to develop his interval for[pic]. These had not been developed before his lifetime. Beginning with a hexagon and increasing the number of sides by a factor of two, Archimedes used a 96-gon, when n=5, to obtain this interval:[pic].

For many years Archimedes’s method remained the only formal method for calculating[pic]. The next section looks at mathematicians’ attempts to calculate [pic] using infinite series.

Gregory’s Infinite Series

James Gregory (1638-1675) and Gottfried Leibniz (1646-1716) both independently developed an infinite series to calculate[pic]. They used the Taylor series expansion of arctangent to find an infinite series that equals[pic]. The Taylor series polynomial for arctangent is[pic]. Choosing[pic], one has

[pic].

Proof That the Series Converges

The alternating series test, or the Leibniz criterion, states that for an alternating series of the form[pic] if the sequence[pic], [pic] is monotone decreasing, and [pic], then the series [pic] converges. Here [pic] ;[pic]; and therefore, the series [pic]converges. In particular, it converges to[pic] .

Table 1 shows the partial sums for the first n terms of the Gregory/Leibniz series multiplied by 4.

|[pic] |[pic] |

|100 |3.1(315929035585537) |

|1000 |3.14(0592653839794) |

|1,000,000 |3.14159(16535897743) |

|1,000,000,000 |3.1415926(2525880504) |

|2,000,000,000 |3.14159265(8505181) |

Table 1

As one can see, the series does converge to[pic], but it does so rather slowly.

More Gregory/Leibniz, John Machin and Machin-Like Formulae

In 1706 John Machin improved upon the Gregory/Leibniz series by finding a series with an increased rate of convergence. The formula that Machin developed was

[pic] .

Machin combined his formula with the Taylor series expansion of arctangent to calculate [pic] (Jovanovic, 2005).

[pic]

Table 2 shows that this series converges much quicker than the Gregory/Leibniz series. By term 11, the limits of the Java primitive data types have been reached. There are two ways of representing decimals in Java. The two data types are called floating point data types, which in Java, are float and double. Data of type float can take values from [[pic],[pic]], resulting in approximately seven digits of accuracy. Variables of type double can assume values [[pic],[pic]], resulting in approximately 17 digits of accuracy. Therefore, the best one could hope to do with Java is 17 digits. However, with rounding error the optimal number of decimal places is actually 16.

|[pic] |[pic] |

|2 |3.14(05970293260603) |

|4 |3.14159(1772184773) |

|10 |3.14159265358979(22) |

|11 |3.14159265358979(4) |

Table 2

The current record for the most digits of[pic], 1,241,100,000,000, was obtained in 2002 by Yasumasa Kanada of Tokyo University. Kanada used a 64-node Hitachi supercomputer with 1 terabyte of main memory to evaluate a Machin-like formula.

[pic]

The previous three methods all dealt with calculating[pic]. While one may use the previous methods to calculate[pic], one may also use simulations to estimate[pic]. In the following section discussion proceeds as to how and why simulations can be used to estimate[pic].

Monte Carlo Method

A Monte Carlo method can be described as a method that solves a problem by generating suitable random numbers and noting that proportion of the numbers obeying some property. Named after the city in Monaco known for its casinos, the method was dubbed Monte Carlo method by Nick Metropolis. Credit for the invention of the method generally goes to Stanislaw Ulam, but there is evidence of its use before the 1940s when it was formalized.

Ulam, inspired by the recent development of the second electronic computer (the ENIAC), suggested to John von Neumann that statistical sampling could be used with more efficiency to solve problems. In particular, von Neumann outlined a statistical approach to solve the problem of neutron diffusion in fissionable material.

The most common application of the Monte Carlo method is Monte Carlo integration. Consider the following integral.

[pic]

[pic]

This integral represents a non-random problem, but suppose there is a random vector U, which is uniformly distributed over the region of integration. Applying the function to U, random variable [pic] is obtained and the expected value of [pic] is[pic], where [pic] is the probability density function of U. Being a pdf, [pic]will equal 1 over the region of integration and, therefore, we have[pic]. Using substitution, [pic] (Monte Carlo Method).

As an example, one can use Monte Carlo integration to find the area under the curve [pic] on the interval[pic]. This is a problem that can be solved with analytical methods and then can be compared to the results when Monte Carlo integration is applied to approximate the area.

The area of the bounded region is given by

[pic][pic]=20.

For Monte Carlo integration, we will determine the maximum value of the function on the interval[pic]. This will correspond to the height of the rectangle. The width of the rectangle will be 2, the width of the interval over which we are integrating. Figure 10 shows the rectangle from which random points will be selected. The number of points that fall under the curve will be counted.

[pic]

Figure 10: Example of Monte Carlo integration.

Monte Carlo integration suggests that

[pic],

where [pic]is the proportion of points found below the curve and [pic]is the area of the rectangular region. Using Monte Carlo integration with 10,000 points, the area is 19.845.

Of course, the true usefulness of Monte Carlo integration comes in integrating a continuous function that is not integrable by analytical means. For example, the normal curve is given by the function,[pic]. Consider the standard normal curve, when [pic] and[pic]. This function cannot be integrated. The following demonstrates the use of Monte Carlo integration to estimate[pic].

The maximum of g(x) occurs at g(0) [pic].3989. Figure 11 illustrates the integrable region bounded by a 1 X 1 square.

[pic]

Figure 11: Standard normal curve. The square represents the space from which random points are chosen for Monte Carlo integration.

Using 10,000 points generated inside the square,

[pic][pic].336.

The empirical rule states that approximately 68% of the area in a normal curve falls within one standard deviation of the mean. With symmetry this means that approximately 34% of the area is between the mean and one standard deviation above the mean. The Monte Carlo method performed very well in both of our examples.

There are two ways to use the Monte Carlo method. One way applies the Monte Carlo method to solve a problem that is deterministic. This was illustrated by the previous two examples of Monte Carlo integration. The second takes advantage of the law of large numbers. Suppose there is an event,[pic], for which one wishes to determine the probability. Let [pic] be the theoretical probability that event occurs. The empirical probability of event A,[pic], can be found by designing a controlled experiment that will reduce or eliminate bias. Suppose that [pic] trials of the experiment are performed. The law of large numbers implies that as[pic],[pic].

One possible source of bias in Monte Carlo methods relates to the random number generator that is used. Numbers used in Monte Carlo methods are not truly random but are pseudorandom. This means that they use an algorithm of some sort to generate the random numbers. A random number generator usually needs a seed value to create a pseudorandom number. The random number generator in Java uses the clock for the seed value.

The following example is presented to illustrate how Monte Carlo methods can be used to estimate[pic].

Buffon’s Needle

In 1733, French naturalist Comte de Georges Buffon posed this question:

Suppose a floor is made of parallel strips of wood, each the same width, and a needle is dropped onto the floor. What is the probability that the needle will lie across a line between two strips?

Buffon presented a solution in 1777. The following is a treatment of the mathematics involved in that solution.

There are two random variables, the angle at which the needle falls ([pic]) and the distance from the midpoint of the needle to the closest line ([pic]). The needle crosses one of the lines if and only if s [pic]. And it follows that as[pic],[pic]. However after graphing[pic], two problems become apparent.

[pic] at [pic] and [pic] when [pic] . Both of these statements are not true. Regardless of the distance between the lines and the length of the needle, [pic] is still a possible outcome of the random variable[pic]. This means that the needle would lie parallel to the two parallel lines, which does not result in a cross.

Also, in order for [pic] to be a valid probability function, the following inequality must be satisfied:[pic]. Therefore, a different rule for[pic]must be found, when[pic].

Define[pic]. When[pic], the needle will cross one of the parallel lines. Figure 13 shows the feasible region bounded by the rectangle of width[pic] and height[pic]. The favorable region, when the needle crosses one of the lines, is shaded. Finding the area of this shaded region will reveal the probability of a cross.

[pic]

Figure 13: The feasible region bounded by the rectangle of width[pic] and height[pic].

The shaded region can be separated into three regions, two congruent parts and a rectangle. Figure 14 shows the rectangle in black surrounded on the left and right by the two congruent regions. The height of the rectangle is obviously[pic]. The width of the rectangle is the difference between the x-coordinates of the intersection of [pic] and[pic]. This width can be found by subtracting the combined width of the congruent regions from the width of the feasible region or[pic].

[pic]

Figure 14: The area below the curve [pic] bounded by the horizontal lines [pic] and [pic] is separated into three regions, the rectangle in black surrounded on the left and right by two congruent regions.

The area of the two congruent regions can be found by finding the area of the first region and doubling it. The area of the first region is found by integrating the curve [pic] from 0 to[pic]. Thus the area of the favorable region is

[pic]

[pic]

[pic]

[pic]

The area of the total feasible region is [pic].

Substituting[pic], one gets

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

Combining this result with the previous expression for [pic] when [pic] yields a piece-wise function. The function is graphed in Figure 15(Weisstein, Buffon’s Needle Problem and Kunkel, September, 4 2003).

[pic]

[pic]

Figure 15: Probability function for Buffon’s needle, where n is the ratio of needle length to distance between lines.

Monte Carlo simulation and this probability function can be utilized to obtain an estimate of[pic]. Of course, this is not a terribly efficient means of finding[pic], especially if one is actually dropping needles. To ensure independence of each trial one would have to remove the dropped needle before dropping the next needle. There are many Java applets online which allow one to simulate the needle dropping with a great deal more efficiency.

An experiment of this nature was conducted by O.C. Fox in 1864 while recovering from wounds sustained during the Civil War. His results are shown in the table below (Monte Carlo Method).

|[pic] drops |[pic]crosses |[pic] inches |[pic] inches |Surface |Estimate of [pic] |

|500 |236 |3 |4 |Stationary |3.1780 |

|530 |253 |3 |4 |Rotated |3.1423 |

|590 |939 |5 |2 |rotated |3.1416 |

Table 3

Using an online applet to simulate the dropping of 100,000,000 needles an estimate of [pic]=3.1415(528912708703) was found, accurate to five digits (Pilana, n.d.).

Results and Mathematics

This paper illustrates a different way to use the Monte Carlo method. Starting by selecting random “points” in a unit n-cube, one counts the number of points that fall inside of a partial n-sphere, centered at the origin. Using this proportion one can obtain an estimate of [pic] by using newly developed functions. In two-dimensional space, one actually chooses points from a unit square with its left vertex at the origin and calculates the proportion of points that fall inside of the ¼ circle with center (0, 0) and radius 1.

[pic]

Figure 16: Quarter circle centered at the origin with radius 1 bounded by a unit square.

In three-dimensional space, one chooses points from a unit 3-cube with a vertex at the origin and three edges along the[pic]and [pic]-axes. Then one must calculate the proportion of points that fall inside the 1/8-sphere. In general, this method involves choosing points from a unit n-cube with a vertex at the origin and with n-edges along the [pic] axes and calculating the proportion of points that fall within the [pic] n-sphere.

First, find the volume of each n-sphere, starting with n=2.

[pic]

Using integration by parts, let [pic] and[pic].

[pic] and [pic]

[pic]

[pic]

Now evaluate the definite integral as follows.

[pic]

=[pic]

The volume of a 2-sphere is[pic]. This is the area of a circle. A 2-sphere is in fact a circle and the “volume” is the area. A sphere is generated when a semicircle is rotated about the x-axis. Consider the function [pic] rotated about the x-axis. The volume of the sphere is given by

[pic]

[pic]

[pic]

=[pic]

[pic]

Figure 17: Cross-sections of a 4-sphere (Fuquay, 2002).

The disk method can be generalized to any n-sphere if one knows a formula for the volume of a cross-section. The volume formula from the (n-1)-sphere is used to represent the cross-sectional volume. Figure 17 shows various cross-sections for the 4-sphere. These cross sections would be 3-spheres. The volume of a 4-sphere can be obtained by using [pic] as the formula for the cross-sectional volume.

[pic]

Substituting [pic] and[pic],

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

Using the power reducing identity[pic], obtains the following results

[pic]

[pic]

[pic][pic][pic][pic][pic]

[pic]

=[pic](Weisstein, Hypersphere)

Unfortunately as the dimensions increase it becomes harder if not impossible to visualize what is taking place. One more example may help clarify this concept. For a 5-sphere, the volume is given by

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

So what does this look like in general? Can it even be generalized? To develop a formula for the volume of an n-sphere, begin by stating that

[pic][pic],

where [pic]is the surface area of an n-sphere with a radius of 1. The unit n-sphere must satisfy

[pic]

[pic](Weisstein, Hypersphere)

A function can now be introduced that will be very helpful in the development of the volume of an n-sphere. The gamma family of distributions is a flexible family of distributions on [0,[pic]).The gamma function is defined as

[pic].

Letting [pic], one has

[pic].

[pic], so

[pic][pic]

[pic],

substituting results in

[pic]

Substituting[pic], results in

[pic]

[pic]

(McDonald, 2003)

Table 4 summarizes the work to the 10th dimension.

|[pic] |[pic] |

|2 |[pic] |

|3 |[pic] |

|4 |[pic] |

|5 |[pic] |

|6 |[pic] |

|7 |[pic] |

|8 |[pic] |

|9 |[pic] |

|10 |[pic] |

Table 4

The development of the volume of an n-cube is a much easier venture. One must first look at the volume of a 2-cube or a square. The volume of a 2-cube is an area, and the area of a square is[pic], where [pic] is the length of an edge. Similarly, the volume of a 3-cube is[pic]. In general, the volume of an n-cube is represented by[pic]. While visualizing n-cubes beyond n=3, two figures have been included to help with that visualization. Figure 18 shows the 4-cube projected into the plane, and Figure 19 shows how Mathematica represents n-cubes from n=1 to n=5.

[pic]

Figure 18: The graph of the 4-cube (tesseract) symmetrically projected into the plane (Weisstein, Hypercube).

[pic]

Figure 19: n-cube, from n=1 to n=5 (Weisstein, Hypercube).

Probabilities

Notation must first be clarified. Let [pic] be the theoretical probability that a random point chosen from U[0, 1) along the unit n-cube will fall inside or on the [pic]n-sphere centered at the origin.

[pic]= The empirical probability that a random point chosen from U[0, 1) along the unit n-cube will fall inside or on the [pic] hyper-sphere centered at the origin.

[pic]

[pic]

[pic]

Probabilities of[pic], [pic], [pic], have been calculated so far.

One will notice that the probabilities decrease as n increases. The volume of the n-cube is then[pic]. It’s easy to see that when[pic], as[pic]. When[pic], the volume of an n-cube remains constant at 1. What happens to the volume of the n-sphere as[pic]? For the ease of calculations, let[pic].

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

After the 5-sphere, the volume of the n-sphere begins to decrease. So while the volume of the unit n-cube is fixed at 1, the volume of the unit n-sphere decreases for dimensions higher than five. Figure 20 shows the results graphically.

[pic]

Figure 20: The volume of an n-sphere as a function of n.

This result is potentially troubling. As the dimensions increase, the number of points that are available to be chosen also increases. A greater number of these points are likely to fall outside of the n-sphere. This means that the number of trials required will increase as the dimensions increase. Is the value of [pic] estimated to be the same in higher dimensions?

After developing probabilities through 10-space one can develop formulas to calculate[pic]. Once the program has generated [pic] random points, the number of points that are within one unit of the origin are counted, represented by C. The proportion of points that fall inside the [pic] n-sphere can now be calculated:[pic].

The law of large numbers implies that as[pic]. Substituting [pic] for [pic] into the probabilities formulas found earlier, one can solve for [pic] to obtain a random variable function that will estimate[pic].

Let [pic] be the estimate for the value of [pic]in the [pic] dimension.

[pic], solving for[pic], we get[pic].

[pic][pic] [pic].

[pic][pic][pic]

[pic]

[pic]

[pic]

[pic]

[pic]

[pic]

Table 5 shows the results of the procedure. The computer was programmed to choose 1 million random points and calculate an estimate for[pic]. This process was repeated for 100 trials to obtain the overall estimate for[pic], second column. The third column represents the closest estimate for [pic] in the 100 trials.

|100,000,000 random points |

|Dimension |[pic] |Closest estimate |

|2 |3.14(33749) |3.1415(04) |

|3 |3.14(247652) |3.141(672) |

|4 |3.14(6113371) |3.1415(76674219491) |

|5 |3.1(38861693) |3.1415(410231286174) |

|6 |3.1(068389807461223) |3.1(1653955588751) |

|7 |3.1(09551964626265) |3.1(19758209812801) |

|8 |3.1(177043676024003) |3.1(321549271198452) |

|9 |3.1(16178352559674) |3.1(368318687296557) |

|10 |3.1(21780390406692) |3.141(414241794067) |

Table 5

As evident from the table, this process does not perform well in dimensions higher than five. This seems to fit with the previous calculations that the volume of the sphere decreases in dimensions higher than five. Though the estimates seem to be getting stronger again, this could easily be attributed to the increase in variation as more variables are added.

So does choosing a larger number of points improve these estimates for higher dimensions? If it does not, then does the value of [pic]remain constant for higher dimensions? It does not seem possible that a constant could change value. But again it seems counterintuitive that the volume of a sphere decreases as one increases dimension beyond[pic]. If one believes that [pic] is in fact still 3.14159265358979.., then how can one improve this process to obtain a closer estimate of[pic]?

A blog was found that discussed n-dimensional cubes and spheres (Lippert, May 11, 2005). One participant mentioned an article that dealt with a similar problem by restricting the size of the n-cube from which points are selected (Goldstein, Platt, and Burges). After contemplating this notion, experiments were constructed in which the section of the Java code that selects the random values that make up the coordinates of the random point was altered. By restricting the values of each coordinate to the intervals represented in the third column of Table 6, instead of the interval [0, 1), one is able to obtain a closer estimate of[pic].

|[pic] |[pic] |interval |

|6 |3.14(3589458142209) |[0, .9946] |

|7 |3.14(5864865018313) |[0, .995] |

|8 |3.141(8563412346634) |[0, .9951] |

|9 |3.154878413601465 |[0, .9952] |

|10 |3.1556855946042437 |[0, .995] |

Table 6

Use in the Classroom

Geometry

For a geometry course, there are two lessons that can be used from this topic. In a chapter on geometric probability, students can be motivated to find the analytical solution to a probability question through numerical means. The use of the Monte Carlo method in two dimensions can be a way of gaining this numerical solution. This provides for an opportunity to teach through problem solving.

Using constructing software (Geometer’s Sketchpad), students can duplicate Archimedes’s method of trapping[pic]. This lesson appears in the appendix.

AP Statistics

In an AP Statistics course, this topic can be used as a way of approaching two concepts that are vital in AP Statistics, the use of Monte Carlo methods to approximate solutions and the idea of confidence intervals. By having students construct confidence intervals for their estimate, students are given practice in these constructions. In addition, as a class some students’ confidence intervals will contain the true value of [pic] and some won’t. This allows me to reemphasize what a C-level of confidence actually means. Since students actually know the parameter that they are trying to estimate, the concept becomes more concrete. Normally throughout the unit on confidence intervals students have no idea what the parameter is. This is good; it’s what they would encounter in practice using statistical methods. However, it is very easy for a student to think that 95% confidence means that there is a 95% chance that once a confidence interval is found, the parameter is inside of that interval. I constantly have to repair that misconception with students. Using the lesson plan that is attached in the appendix allows them to see my point that the parameter is either in the interval or it is not. The confidence comes in the procedure that is used to construct the confidence interval.

When I used the lesson attached this last spring, 13/15 students had 95% confidence intervals that captured[pic]. I was able to talk about sample size and how that affects our ability to gain a good estimate of the parameter. Students could see that by increasing the sample size, we were able to decrease the width of the confidence interval. By this time, most students had a firm grasp of that concept, but this is another way of bringing across that point.

We also discussed the limits of this method for calculating[pic]. On the TI-83 Plus calculator, using as few as 500 points causes an overload in the memory. The random number generator on the TI-83 calculator also contains some flaws that are revealed when conducting the lesson. This spring, three students obtained the exact same results. This allowed me to talk about seeds and why there might be bias in our class wide results. By using the TI-Navigator system, I was able to quickly collect everyone’s data in a few seconds.

AP Computer Science

The major goals of the AP Computer Science program are laid out as:

▪ Design and implement computer-based solutions to problems in a variety of application areas

▪ Use and implement well-known algorithms and data structures

▪ Develop and select appropriate algorithms and data structures to solve problems

▪ Code fluently in an object-oriented paradigm using the programming language Java. Students are expected to be familiar with and be able to use standard Java library classes from the AP Java subset.

▪ Read and understand a large program consisting of several classes and interacting objects. Students should be able to read and understand a description of the design and development process leading to such a program

(College Board)

The assignment that I have included in the appendix addresses several of these goals. Students are asked to design and implement a solution based on a description of the problem. Students need to use “helper” classes to fully complete the assignment. The Point class is from a previous assignment. Students create a variety of classes that deal with geometry. Point is the first class that the students create and they are well versed in its implementation. The Random class is a standard Java library class from the AP Java subset. It is important that students can read the documentation for standard classes and be able to implement them into their programs.

This project asks that they work in the object-oriented paradigm. It is an object-oriented assignment and not procedural.

References

Bohr, J. (2000, August 29). Machin’s method of approximating pi. Retrieved May 20, 2006, from

Buffon’s Needle. (2006, April 29). Retrieved May 20, 2006, from

wiki/Buffon's_needle

Clawson, C. C. (1996). Mathematical mysteries: The beauty and magic of numbers. Cambridge, MA: Perseus Books.

Coexeter, H. S. M. (1973). Regular Polytopes (3rd ed.). New York: Dover.

CollegeBoard AP. (2006). Computer Science Course Description. Accessed at

Fuquay, J. (2002). The Hypersphere. Retrieved June 2, 2006, from

Goldstein, Platt, and Burges. (n.d.). Redundant bit vectors for quickly searching high-dimensional regions. Microsoft Research.

Groleau, R. (2003, September). Approximating pi. Retrieved May 24, 2006, from

Gusmorino, P.A., III. (1998). Paul’s page of pi: Gregory, Leibniz, and Machin. Retrieved June 1, 2006, from

John Machin. (2006, March 12). Retrieved May 16, 2006, from

wiki/John_Machin

Jovanovic, R. (2005, January). Machin’s Formula. Retrieved June 10, 2006, from

Khamis, H. (n.d.). Buffon’s Needle problem. Retrieved June 2, 2006, from Wright State University’s Web site:

buffons_needle_problem/index.htm

Kunkel, P. (2003, September 4). Buffon’s Needle. Retrieved May 20, 2006, from

Larson, Hostetler, and Edwards. (2006). Calculus of a Single Variable Eighth Edition. Boston, MA: Houghton Mifflin Company.

Linn, S., and Neal, D. (2006). Approximating π with the Golden Ratio. Mathematics Teacher, 99(7), 472–477.

Lippert, E. (2005, May 13). High-dimensional spaces are counterintuitive, part two. Message posted to

2005/05/13/417250.aspx

McDonald, K. T. (2003, February). Volume and surface area of an N-Sphere. Lecture given at Princeton University, Princeton, NJ.

Metropolis, N. (1987). Los Alamos Science. [Special issue].

Monte Carlo Method. (n.d.). Retrieved May 13, 2006, from

O’Connor, J. J, and Robertson, E. F. (2001, August). A history of pi. Retrieved April 30, 2006, from

Ortiz, A. L. (1998, February 23). Formula for the surface area of a sphere in Euclidean N-Space. Retrieved May 17, 2006, from University of Waterloo David R. Cheriton School of Computer Science’s Web site:

Peterson, I. (2003). The bias of random-number generators. MAA Online. Retrieved May 3, 2006, from

Pilana, S. (n.d.) Buffon’s Needle experiment. Retrieved May 22, 2006, from

Reese, G. (n.d.). Buffon’s Needle: An analysis and simulation. Retrieved May 22, 2006, from the Office of Mathematics, Science, and Technology Education at the University of Illinois at Urbana-Champaign’s Web site:

Smoller, L. (2001, February). The amazing history of pi. Retrieved June 3, 2006, from University of Arkansas at Little Rock Department of History’s Web site:

Weissten, E. (1999). Hyperspace. MathWorld—A Wolfram Web Resource. Retrieved May 1, 2006 from

Weisstein, E. (1999). Machin’s Formula. MathWorld—A Wolfram Web Resource. Retrieved May 1, 2006 from

Weisstein, E. (1999). Monte Carlo Method. MathWorld—A Wolfram Web Resource. Retrieved May 1, 2006 from

MonteCarloMethod.html

Weisstein, E. (2004, January 25). Machin-like formulas. MathWorld—A Wolfram Web Resource. Retrieved May 3, 2006 from

Weisstein, E. (2005, December 12). Buffon’s Needle Problem. MathWorld—A Wolfram Web Resource. Retrieved May 7, 2006, from

Weisstein, E. (2006, March 15). Hypercube. MathWorld—A Wolfram Web Resource. Retrieved May 3, 2006 from

Weisstein, E. (2006, April 17). Transcendental numbers. MathWorld—A Wolfram Web Resource. Retrieved May 1, 2006 from

TranscendentalNumbers.html

Weisstein, Eric. (2006, April 28). Algebraic numbers. MathWorld—A Wolfram Web Resource. Retrieved May 1, 2006 from

APPENDIX

Archimedes’s Estimation of [pic]

Geometric Constructions Project

Brian Reece

Johnston High School

brian.reece@johnston.k12.ia.us

NCTM Standards

Geometry

▪ Analyze characteristics and properties of two-dimensional geometric shapes and develop mathematical arguments about geometric relationships

▪ Use visualization, spatial reasoning, and geometric modeling to solve problems

Archimedes was the first to formalize an approximation of[pic]. You are going to duplicate the first few steps of his method.

1. Begin by opening Geometer’s Sketchpad

2. Construct a circle with an arbitrary radius

3. Construct a hexagon that is inscribed inside of the circle that you constructed. (You may want to hide some objects as you go)

[pic]

4. Construct a hexagon that superscribes the circle

[pic]

5. Using the measurement tool of Sketchpad, calculate the semi-perimeter of each of the two hexagons that you constructed.

6. Now using what you’ve already done construct a dodecagon that is inscribed in the circle and one that superscribes the circle.

7. Calculate the semiperimeters of the two dodecagons.

Archimedes gave the area of a circle as[pic]. Using this formula for the area of a circle and the more familiar formula for the area of a circle ([pic]) find an equation for[pic]. Then use the semiperimeters you found and any other information you need to find an estimate for[pic].

Estimating Pi through Simulation

Brian Reece

Johnston High School

brian.reece@johnston.k12.ia.us

NCTM Standards:

Data Analysis and Probability

▪ use simulations to explore the variability of sample statistics from a known population and to construct sampling distributions

▪ use simulations to construct empirical probability distributions

Geometry

▪ use geometric models to gain insights into, and answer questions in, other areas of mathematics

AP Statistics Content covered:

“Law of large numbers” concept

Simulation of probability distributions

Large sample confidence interval for a proportion

Materials needed:

TI-83 or TI-84

Pencil

TI-Navigator System

1. The program SIM is distributed to students from Instructor using the TI-Navigator System. This program is included in this handout.

2. Talk about what the calculator is doing. I share the program with the class so that they realize how the process is working. The calculator is selecting random x-values between 0 and 1 and random y-values between 0 and 1. The calculator is then going to count the number of times a point is in the region bounded by the first quadrant of the x-y plane and a circle centered at the origin with a radius of 1.

3. Students are asked to run 10 trials of size 50 on their calculator. Students are then asked to record the number of points in this bounded region for each simulation in list L4.

4. Next students run 10 simulations of size 100 on their calculator. Record the number of points in the bounded region in list L5.

[pic]

[pic]

Depending on the objectives of your class there are different directions you can take at this point. The handout I have included is for my AP Statistics course and uses the concept of confidence intervals. This could also be used in geometry by removing the construction of confidence intervals and instead focusing on the probability of selecting a point that is inside of the bounded region (Ä/4). I also focus on the drawbacks of this process (memory, time, efficiency, etc.)

PROGRAM: SIM

0(L

0(I

Disp “HOW MANY POINTS?”

Prompt N

rand(N)(L1

rand(N)(L2

ClrDraw

FnOff

0(Xmin

1(Xmax

.1(Xscl

0(Ymin

1(Ymax

.1(Yscl

Plot1(Scatter,L1,L2)

Circle(0,0,1)

Pause

(((L1 - 0)2 + (L2 - 0)2)(L3

For(I,1,N,1)

If (L3(I) < 1 or L3(I) = 1)

Then

L+1 (L

End

End

ClrHome

Disp “THE NUMBER

Disp “OF POINTS IN

Disp “THE REGION

Disp “IS

Disp L

Estimating through Simulation

Student Handout

1. Get program SIM from Instructor

2. The calculator is selecting random x values between 0 and 1 and random y-values between 0 and 1. The calculator is then going to count the number of times a point is in the region bounded by the first quadrant of the x-y plane and a circle centered at the origin with a radius of 1. You are to record the number of points in this bounded region for each of ten simulations in list L4.

3. You are going to run 10 trials of size 50 on your calculator.

4. Now run 10 simulations of size 100 on your calculator. Record the number of points in the bounded region in list L5.

[pic]

[pic]

5. Using your 10 simulations of size 50 as one large sample of size 500, construct a 95% confidence interval for the proportion of points that can be found inside the bounded region.

6. Using your 10 simulations of size 100 as one large sample of size 1000, construct a 95% confidence interval for the proportion of points that can be found inside the bounded region.

7. Multiply the bounds of each confidence interval by 4. What have we just estimated? Why?

Estimating Pi through Simulation

AP Computer Science

Brian Reece

Johnston High School

brian.reece@johnston.k12.ia.us

AP Computer Science Objectives:

▪ Design and implement computer-based solutions to problems in a variety of application areas

▪ Use and implement well-known algorithms and data structures

▪ Develop and select appropriate algorithms and data structures to solve problems

▪ Code fluently in an object-oriented paradigm using the programming language Java. Students are expected to be familiar with and be able to use standard Java library classes from the AP Java subset.

▪ Read and understand a large program consisting of several classes and interacting objects. Students should be able to read and understand a description of the design and development process leading to such a program

Specifications:

Your program should

▪ use the Random class (Built-in Java class) Review API online.

▪ use the Point class (user defined- you should already have this class)

▪ choose n random points from [0,1). n is defined by the user

▪ test whether the points are inside of the circle defined by the equation[pic].

▪ calculate the proportion of total points that lie inside of the circle.

▪ estimate [pic] (Hint: Geometric probability)

▪ Output the estimated value of [pic].

▪ include all necessary comments

Extra credit:

Alter your code so that the user can define the center and the radius of the circle.

import java.util.*;

import java.io.*;

public class PiEstimator {

/**

* Method main

*

*

* @param args

*

*/

public static void main(String[] args) {

File file = null;

long timeBefore = System.currentTimeMillis();

int repetitions =1000000;

int numberOfTrials = 100;

double OverallEstimate = 0;

if (args.length > 0) file = new File (args[0]);

if (file == null ) {

System.out.println ("Default: 10_DtextOutput.txt");

file = new File ("10_DtextOutput.txt");

}

try {

// Create a FileWriter stream to the file

FileWriter file_writer = new FileWriter (file);

// Put a buffered wrapper around it

BufferedWriter buf_writer = new BufferedWriter (file_writer);

// In turn, wrap the PrintWriter stream around this output stream

// and turn on the autoflush.

PrintWriter print_writer = new PrintWriter (buf_writer,true);

print_writer.println("1,000,000 points in 10-D space");

print_writer.println();

for(int j=0;j ................
................

In order to avoid copyright disputes, this page is only a partial summary.

Google Online Preview   Download