|
| 1 | +def bailey_borwein_plouffe(digit_position: int, precision: int = 1000) -> str: |
| 2 | + """ |
| 3 | + Implement a popular pi-digit-extraction algorithm known as the |
| 4 | + Bailey-Borwein-Plouffe (BBP) formula to calculate the nth hex digit of pi. |
| 5 | + Wikipedia page: |
| 6 | + https://en.wikipedia.org/wiki/Bailey%E2%80%93Borwein%E2%80%93Plouffe_formula |
| 7 | + @param digit_position: a positive integer representing the position of the digit to extract. |
| 8 | + The digit immediately after the decimal point is located at position 1. |
| 9 | + @param precision: number of terms in the second summation to calculate. |
| 10 | + A higher number reduces the chance of an error but increases the runtime. |
| 11 | + @return: a hexadecimal digit representing the digit at the nth position |
| 12 | + in pi's decimal expansion. |
| 13 | + |
| 14 | + >>> "".join(bailey_borwein_plouffe(i) for i in range(1, 11)) |
| 15 | + '243f6a8885' |
| 16 | + >>> bailey_borwein_plouffe(5, 10000) |
| 17 | + '6' |
| 18 | + >>> bailey_borwein_plouffe(-10) |
| 19 | + Traceback (most recent call last): |
| 20 | + ... |
| 21 | + ValueError: Digit position must be a positive integer |
| 22 | + >>> bailey_borwein_plouffe(0) |
| 23 | + Traceback (most recent call last): |
| 24 | + ... |
| 25 | + ValueError: Digit position must be a positive integer |
| 26 | + >>> bailey_borwein_plouffe(1.7) |
| 27 | + Traceback (most recent call last): |
| 28 | + ... |
| 29 | + ValueError: Digit position must be a positive integer |
| 30 | + >>> bailey_borwein_plouffe(2, -10) |
| 31 | + Traceback (most recent call last): |
| 32 | + ... |
| 33 | + ValueError: Precision must be a nonnegative integer |
| 34 | + >>> bailey_borwein_plouffe(2, 1.6) |
| 35 | + Traceback (most recent call last): |
| 36 | + ... |
| 37 | + ValueError: Precision must be a nonnegative integer |
| 38 | + """ |
| 39 | + if (not isinstance(digit_position, int)) or (digit_position <= 0): |
| 40 | + raise ValueError("Digit position must be a positive integer") |
| 41 | + elif (not isinstance(precision, int)) or (precision < 0): |
| 42 | + raise ValueError("Please input a nonnegative integer for the precision") |
| 43 | + |
| 44 | + # compute an approximation of (16 ** (n - 1)) * pi whose fractional part is mostly accurate |
| 45 | + sum_result = ( |
| 46 | + 4 * _subsum(digit_position, 1, precision) |
| 47 | + - 2 * _subsum(digit_position, 4, precision) |
| 48 | + - _subsum(digit_position, 5, precision) |
| 49 | + - _subsum(digit_position, 6, precision) |
| 50 | + ) |
| 51 | + |
| 52 | + # return the first hex digit of the fractional part of the result |
| 53 | + return hex(int((sum_result % 1) * 16))[2:] |
| 54 | + |
| 55 | + |
| 56 | +def _subsum( |
| 57 | + digit_pos_to_extract: int, denominator_addend: int, precision: int |
| 58 | +) -> float: |
| 59 | + # only care about first digit of fractional part; don't need decimal |
| 60 | + """ |
| 61 | + Private helper function to implement the summation |
| 62 | + functionality. |
| 63 | + @param digit_pos_to_extract: digit position to extract |
| 64 | + @param denominator_addend: added to denominator of fractions in the formula |
| 65 | + @param precision: same as precision in main function |
| 66 | + @return: floating-point number whose integer part is not important |
| 67 | + """ |
| 68 | + sum = 0.0 |
| 69 | + for sum_index in range(digit_pos_to_extract + precision): |
| 70 | + denominator = 8 * sum_index + denominator_addend |
| 71 | + exponential_term = 0.0 |
| 72 | + if sum_index < digit_pos_to_extract: |
| 73 | + # if the exponential term is an integer and we mod it by the denominator before |
| 74 | + # dividing, only the integer part of the sum will change; the fractional part will not |
| 75 | + exponential_term = pow( |
| 76 | + 16, digit_pos_to_extract - 1 - sum_index, denominator |
| 77 | + ) |
| 78 | + else: |
| 79 | + exponential_term = pow(16, digit_pos_to_extract - 1 - sum_index) |
| 80 | + sum += exponential_term / denominator |
| 81 | + return sum |
| 82 | + |
| 83 | + |
| 84 | +if __name__ == "__main__": |
| 85 | + import doctest |
| 86 | + |
| 87 | + doctest.testmod() |
0 commit comments