Much of the underlying math is on Wikipedia. First, Continued_fraction: Best rational within an interval gets us to look at how many terms of the continued fraction expansion of left
and right
match. The previous section, Continued_fraction: Best rational approximations, helps us to choose the single value to tack onto the end of the continued fraction expansion. Admittedly these Wikipedia sections could be better written. Admittedly, a primary or secondary source might be better than Wikipedia here.
Here is some code that implements the idea. Please try it out to see whether / where it does not do what you would hope.
from fractions import Fraction
import math
import re
def best_fraction(target, left_open=False, right_open=True):
"""
best_fraction computes a low-denominator fraction that
approximates the input `target` string.
This is much in the spirit of
fractions.Fraction.limit_denominator() but, instead of using a
default or user-specified maximum denominator, the code assumes
that the input `target` is a string representing a rounded decimal
value. It infers `left` as the lowest value that `target` could
have been rounded from, and `right` as the highest value that
`target` could have been rounded from. By default, the subroutine
then returns a Fraction that approximates `target` with low
denominator, from within the half-open interval [left, right).
Whether the interval from left to right should be closed, open, or
half-open the other way can be changed by changing the values of
`left_open` and `right_open`.
The implementation uses a continued fraction expansion of the
`left` and `right` values. These will match for a number of terms
but then will have a first term where they do not match. The
returned value comes from the matching part of continued fraction
expansion plus a single additional term. The single additional
term is chosen so that the fraction will be nearest the `target`
value.
For example, "3.14" could be assumed to be rounded from any value
in [3.135, 3.145). best_fraction("3.14") will return Fraction(22,
7) as the best fraction. In contrast,
Fraction("3.14").limit_denominator() returns Fraction(157, 50).
"""
# Check types of inputs
if not isinstance(target, str):
raise ValueError(f"target ({target}) must be a str")
if not isinstance(left_open, bool):
raise ValueError(
f"left_open ({left_open}) must be True or False"
)
if not isinstance(right_open, bool):
raise ValueError(
f"right_open ({right_open}) must be True or False"
)
# If supplied string is malformed, this will throw an exception
target_fraction = Fraction(target)
# If `target` is positive then we can figure the right end of the
# rounding interval by attaching a "5" suffix to the mantissa. If
# the mantissa does not include a decimal point then append a
# decimal point before appending the "5".
right = re.sub(
r"^([+-]?[0-9_]*)(\.([0-9_]*))?([eE][+-]?[0-9_]+)?$",
r"\g<1>.\g<3>5\g<4>",
target,
)
if right == target:
raise ValueError(
f"target ({target}) not recognized as a decimal number"
)
right_fraction = Fraction(right)
# The left side of the rounding interval is equidistant from the
# target, but on the other side.
left_fraction = 2 * target_fraction - right_fraction
# If `target` is negative, we have left and right swapped
if left_fraction > right_fraction:
(left_fraction, right_fraction) = (
right_fraction,
left_fraction,
)
# We now know the rounding interval
return best_fraction_from_interval(
target_fraction,
left_fraction,
right_fraction,
left_open,
right_open,
)
def best_fraction_from_interval(
target, left, right, left_open=False, right_open=True
):
"""
best_fraction_from_interval computes a low-denominator fraction
that approximates the input `target`; much in the spirit of
fractions.Fraction.limit_denominator().
Specifically, this subroutine assumes that `target` is a rounded
decimal value, `left` is the lowest value that `target` could have
been rounded from, and `right` is the highest value that `target`
could have been rounded from. By default, the subroutine then
returns a Fraction that approximates `target` with low
denominator, from within the half-open interval [left, right).
Whether the interval from left to right should be closed, open, or
half-open the other way can be changed by changing the values of
`left_open` and `right_open`.
The implementation uses a continued fraction expansion of the
`left` and `right` values. These will match for a number of terms
but then will have a first term where they do not match. The
returned value comes from the matching part of continued fraction
expansion plus a single additional term. The single additional
term is chosen so that the fraction will be nearest the `target`
value.
For example, "3.14" could be assumed to be rounded from any value
in [3.135, 3.145). best_fraction_from_interval(Fraction("3.14"),
Fraction("3.135"), Fraction("3.145")) will return Fraction(22, 7)
as the best fraction.
"""
# Check types of inputs
if not isinstance(target, Fraction):
raise ValueError(
f"target value ({target}) must be a Fraction"
)
if not isinstance(left, Fraction):
raise ValueError(f"left value ({left}) must be a Fraction")
if not isinstance(right, Fraction):
raise ValueError(f"right value ({right}) must be a Fraction")
if not target > left:
raise ValueError(
f"target value ({target}) must be greater than left value ({left})"
)
if not target < right:
raise ValueError(
f"target value ({target}) must be less than right value ({right})"
)
if not isinstance(left_open, bool):
raise ValueError(
f"left_open ({left_open}) must be True or False"
)
if not isinstance(right_open, bool):
raise ValueError(
f"right_open ({right_open}) must be True or False"
)
# Define useful helper function
def update_fraction(endpoint, endpoint_plus, endpoint_floor):
endpoint = (
1 / (endpoint - endpoint_floor)
if endpoint != endpoint_floor
else math.inf
)
endpoint_plus = not endpoint_plus
endpoint_floor = (
math.inf
if endpoint == math.inf
else math.floor(endpoint)
if endpoint_plus
else math.ceil(endpoint) - 1
)
return (endpoint, endpoint_plus, endpoint_floor)
# Is each endpoint effectively plus epsilon (as opposed to minus
# epsilon)?
left_plus = left_open
right_plus = not right_open
# To get to the next iteration of the continued fraction
# expansion, we need to know the floor of each endpoint. Note
# that an exact integer has a "floor" that depends upon whether
# our value is implicitly +epsilon or -epsilon.
left_floor = (
math.floor(left) if left_plus else math.ceil(left) - 1
)
right_floor = (
math.floor(right) if right_plus else math.ceil(right) - 1
)
# The start up convergents for a continued fraction expansion
(prev_numer, prev_denom) = (0, 1)
(curr_numer, curr_denom) = (1, 0)
while left_floor == right_floor:
# Compute one more convergent and retain one previous
# convergent
(curr_numer, prev_numer) = (
left_floor * curr_numer + prev_numer,
curr_numer,
)
(curr_denom, prev_denom) = (
left_floor * curr_denom + prev_denom,
curr_denom,
)
# Update the values so that we can compute the next
# convergents
(left, left_plus, left_floor) = update_fraction(
left, left_plus, left_floor
)
(right, right_plus, right_floor) = update_fraction(
right, right_plus, right_floor
)
# We know the continued fraction so far, but we have reached the
# point where the continued fraction expanion for `left` and
# `right` are no longer the same. We need to choose the last term
# within the interval and nearest the target. In the general case
# there will be a nearest to the left and a nearest to the right.
best_floor_unrounded = Fraction(
prev_numer * target.denominator
- prev_denom * target.numerator,
curr_denom * target.numerator
- curr_numer * target.denominator,
)
best_floor_left = max(
math.floor(best_floor_unrounded),
min(left_floor, right_floor) + 1,
)
best_floor_right = min(
math.ceil(best_floor_unrounded), max(left_floor, right_floor)
)
answer_left = Fraction(
best_floor_left * curr_numer + prev_numer,
best_floor_left * curr_denom + prev_denom,
)
answer_right = Fraction(
best_floor_right * curr_numer + prev_numer,
best_floor_right * curr_denom + prev_denom,
)
if abs(answer_left - target) <= abs(answer_right - target):
return answer_left
return answer_right