Recipe 6.21. Calculating to Thousands of Digits


Recipe 6.21. Calculating π to Thousands of Digits

Problem

You want to impress people by showing how quickly Visual Basic 2005 can calculate π to a thousand or more decimal places. While you're at it, you might want to discover how to create multidigit mathematical functions using integer arrays of digits.

Solution

Sample code folder: Chapter 06\ CalculatePi

Create functions for basic mathematical operations (+, -, *, /) that operate on integer arrays of any reasonable size. Then demonstrate these functions by calculating π to many digits using one of the standard π-calculation algorithms.

Discussion

This recipe includes a module called PiCalculator that contains the functions needed to perform multidigit math, along with one to calculate π to any number of digits. The four main multidigit functions are named ArrayMult(), ArrayDivide(), ArrayAdd(), and ArraySub(). These are declared as Private to the module because they serve only as support routines to the FindPi() function, but you can change them to Public to experiment with them for other purposes. Other supporting functions include ArrayZero(), which sets all "digits" in an array to zeros, and ArcTangent(), which calls the other functions to calculate the arctangent of a multi-digit number.

The way the basic math functions work is similar to the way math is performed on paper by grade-schoolers: when two digits are added, any overflow is added into the next pair of digits, and so on. Calculating π to 500 decimal places requires a huge number of these small repetitive calculations, but that's what computers are really good at doing.

Here is the code to calculate π. It is based on the following calculation for π:

π/4 = (arctan 1/2)+ (arctan 1/3)

Each part of the algorithm is performed manually, including the arctangent calculation:

 Module PiCalculator    Private NumberDigits As Integer    Public Function FindPi(ByVal digits As Integer) As String       ' ----- Calculate Pi to the specified number of digits,       '       based on the formula:       '          Pi/4 = arctan(1/2) + arctan(1/3)       Dim result As New System.Text.StringBuilder("PI=3.")       Dim digitIndex As Integer       Dim divFactor As Integer       ' ----- Build an array that will hold manual calculations.       NumberDigits = digits + 2       Dim targetValue(NumberDigits) As Integer       Dim sourceValue(NumberDigits) As Integer       ' ----   Perform the calculation.       divFactor = 2       ArcTangent(targetValue, sourceValue, divFactor)       divFactor = 3       ArcTangent(targetValue, sourceValue, divFactor)       ArrayMult(targetValue, 4)       ' ----- Return a string version of the calculation.       For digitIndex = 1 To NumberDigits - 3          result.Append(Chr(targetValue(digitIndex) + Asc("0"c)))       Next digitIndex       Return result.ToString    End Function    Private Sub ArrayMult(ByRef baseNumber() As Integer, _          ByRef multiplier As Integer)       ' ----- Multiply an array number by another number by hand.       ' The product remains in the array number.       Dim carry As Integer       Dim position As Integer       Dim holdDigit As Integer       ' ----- Multiple each base digit, from right to left.       For position = NumberDigits To 0 Step -1          ' ----- If the multiplication went past 9, carry the          ' tens value to the next column.          holdDigit = (baseNumber(position) * multiplier) + carry          carry = holdDigit \ 10          baseNumber(position) = holdDigit Mod 10       Next position    End Sub      Private Sub ArrayDivide(ByRef dividend() As Integer, ByRef divisor As Integer)       ' ----- Divide an array number by another number by hand.       '       The quotient remains in the array number.       Dim borrow As Integer       Dim position As Integer       Dim holdDigit As Integer       ' ----- Process division for each digit.       For position = 0 To NumberDigits          ' ----- If the division can't happen directly, borrow from          '       the previous position.          holdDigit = dividend(position) + borrow * 10          dividend(position) = holdDigit \ divisor          borrow = holdDigit Mod divisor       Next position    End Sub    Private Sub ArrayAdd(ByRef baseNumber() As Integer, ByRef addend() As Integer)       ' ----- Add two array numbers together.       '       The sum remains in the first array number.       Dim carry As Integer       Dim position As Integer       Dim holdDigit As Integer       ' ----- Add each digit from right to left.       For position = NumberDigits To 0 Step -1          ' ----- If the sum goes beyond 9, carry the tens          '       value to the next column.          holdDigit = baseNumber(position) + addend(position) + carry          carry = holdDigit \ 10          baseNumber(position) = holdDigit Mod 10       Next position    End Sub    Private Sub ArraySub(ByRef minuend() As Integer, ByRef subtrahend() As Integer)       ' ----- Subtract one array number from another.       '       The difference remains in the first array number.       Dim borrow As Integer       Dim position As Integer       Dim holdDigit As Integer       ' ---- Subtract the digits from right to left.       For position = NumberDigits To 0 Step -1          ' ----- If the subtraction would give a negative value          '       for a column, we will have to borrow.          holdDigit = minuend(position) - subtrahend(position) + 10          borrow = holdDigit \ 10          minuend(position) = holdDigit Mod 10          If (borrow = 0) Then minuend(position - 1) -= 1       Next position    End Sub    Private Function ArrayZero(ByRef baseNumber() As Integer) As Boolean       ' ----- Report whether an array number is all zero.       Dim position As Integer       ' ----- Examine each digit.       For position = 0 To NumberDigits          If (baseNumber(position) <> 0) Then             ' ----- The number is nonzero.             Return False          End If       Next position       ' ----- The number is zero.       Return True    End Function    Private Sub ArcTangent(ByRef targetValue() As Integer, _          ByRef sourceValue() As Integer, _          ByVal divFactor As Integer)       ' ----- Calculate an arctangent of a fraction,       '       1/divFactor. This routine   performs a modified       '       Maclaurin series to calculate the arctangent.       '       The base formula is:       '         arctan(x) = x - x^3/3 + x^5/5 -' x^7/7 + x^9/9 - …       '       where -1 < x < 1 (1/divFactor in this case).       Dim workingFactor As Integer       Dim incremental As Integer       ' ----- Figure out the "x" part, 1/divFactor.       sourceValue(0) = 1       incremental = 1       workingFactor = divFactor       ArrayDivide(sourceValue, workingFactor)       ' ----- Add "x" to the total.       ArrayAdd(targetValue, sourceValue)       Do          ' ----- Perform the "- (x^y)/y" part.          ArrayMult(sourceValue, incremental)          workingFactor = divFactor * divFactor          ArrayDivide(sourceValue, workingFactor)          incremental += 2          workingFactor = incremental          ArrayDivide(sourceValue, workingFactor)          ArraySub(targetValue, sourceValue)          ' ----- Perform the "+ (x^y)/y" part.          ArrayMult(sourceValue, incremental)          workingFactor = divFactor * divFactor          ArrayDivide(sourceValue, workingFactor)          incremental += 2          workingFactor = incremental          ArrayDivide(sourceValue, workingFactor)          ArrayAdd(targetValue, sourceValue)       Loop Until ArrayZero(sourceValue)    End Sub End Module 

To exercise these procedures, the following statement uses the FindPi() function to calculate π to 500 digits:

 MsgBox(FindPi(500)) 

You can change the 500 argument to obtain a different number of digits. However, even though the time required to calculate π to 500 or even 1,000 digits is fairly negligible, every time you double the count, the FindPi() function requires around four times as long to return the results. Try smaller counts first, moving up to larger counts when you have a good feel for just how long the calculation will take on your computer.

Figure 6-21 shows the first 500 digits of π as formatted by the FindPi() function. If you prefer to format the digits differently, say in groups of 10 digits or with occasional end-of-line characters, you might want to change FindPi() to return an array of digits. The calling code can then format the digits as desired. The "digits" in the array have values in the range 0 to 9, and they need to be converted to ASCII digits by adding the ASCII equivalent of "0" (zero) to their value before applying the Chr() conversion function.

Figure 6-21. Pi calculated to 500 decimal places


See Also

There are many places on the Web to see many digits of π and to learn of the different algorithms used for calculating π to even millions of decimal places. See, for example, http://www.exploratorium.edu/pi/Pi10-6.html.




Visual Basic 2005 Cookbook(c) Solutions for VB 2005 Programmers
Visual Basic 2005 Cookbook: Solutions for VB 2005 Programmers (Cookbooks (OReilly))
ISBN: 0596101775
EAN: 2147483647
Year: 2006
Pages: 400

flylib.com © 2008-2017.
If you may any questions please contact us: flylib@qtcs.net