寻找不同类型的 divmod 函数答案

【问题标题】：Looking for a different kind of divmod function寻找不同类型的 divmod 函数
【发布时间】：2012-08-08 13:33:07
【问题描述】：

Python 的divmod 函数工作正常，这几乎是我想要的。但是，对于需要执行的操作，它对非整数的行为需要稍有不同。运行以下代码时，您可能会看到正在尝试执行的操作。

>>> function = divmod
>>> from math import pi
>>> function(pi * pi, pi) == (pi, 0)
False
>>>

如何在上面定义function，使得最终表达式的计算结果为True，而不是False？如果有人能弄清楚如何得到(pi, 0)而不是(3.0, 0.4448...) ，这就是答案。

编辑 1： 现在来看一个更复杂的例子，下面的代码应该产生[3, 2, 1, 3, 2, 1]。

>>> x = 1 * pi ** 5 + \
        2 * pi ** 4 + \
        3 * pi ** 3 + \
        1 * pi ** 2 + \
        2 * pi ** 1 + \
        3 * pi ** 0
>>> digits = []
>>> while x:
        x, y = function(x, pi)
        digits.append(y)


>>> digits
[0.3989191524449005, 0.2212554774328268, 2.309739581793931, 0.1504440784612413,
2.858407346410207, 1.0]
>>>

编辑 2： 以下显示的代码工作正常，只是它具有意外但有效的输出。

import math

def convert_dec_to_pi(number):
    digits = get_pi_digits(number)
    digits, remainder = correct_pi_digits(digits)
    return make_pi_string(digits, remainder)

def get_pi_digits(number):
    digits = []
    while number:
        number, digit = divmod(number, math.pi)
        digits.append(digit)
    digits.reverse()
    return digits

def correct_pi_digits(digits):
    last = len(digits) - 1
    for index, digit in enumerate(digits):
        if index < last and digit % 1 != 0:
            a, b = get_digit_options(digit, digits[index + 1])
            digits[index:index+2] = a if 0 <= a[1] < math.pi else b
    digit, remainder = divmod(digits[-1], 1)
    digits[-1] = digit
    return digits, remainder

def get_digit_options(digit, next_digit):
    a, b = math.floor(digit), math.ceil(digit)
    if a not in range(4):
        return (b, (digit - b) * math.pi + next_digit), None
    if b not in range(4):
        return (a, (digit - a) * math.pi + next_digit), None
    c, d = ((a, (digit - a) * math.pi + next_digit),
            (b, (digit - b) * math.pi + next_digit))
    return (c, d) if digit - a < 0.5 else (d, c)

def make_pi_string(digits, remainder):
    return '{} base \u03C0 + {} base 10'.format(
        ''.join(str(int(d)) for d in digits), remainder)

下面的函数可以用来逆运算并检查结果。

import re

def convert_pi_to_dec(string):
    match = re.search('^(\\d+) base \u03C0 \\+ (0\\.\\d+) base 10$', string)
    if not match:
        raise ValueError()
    digits, remainder = match.groups()
    return sum(int(x) * math.pi ** y for y, x in enumerate(reversed(digits))) \
           + float(remainder)

以下代码不会引发AssertionError，因此很明显一切正常。

for n in range(1, 36):
    value = convert_dec_to_pi(n)
    print(value)
    assert convert_pi_to_dec(value) == n

那么这让我想到了下面的例子。输出可以毫无问题地被转换回来，但人们会期待一些稍微不同的东西。

>>> convert_dec_to_pi(math.pi * math.pi)
'30 base π + 0.44482644031997864 base 10'
>>> convert_pi_to_dec(_) == math.pi * math.pi
True
>>>

字符串应该是100 base π + 0.0 base 10。此时输出是准确的，但不是“正确的”。

编辑 3：以下示例可能会提供一些额外的洞察力来了解我所追求的。在运行具有不同 π 幂的循环后，我希望所有输出的形式都是 10... base π + 0.0 base 10。结果与此不同，如下所示。

>>> for power in range(20):
    print(convert_dec_to_pi(math.pi ** power))


1 base π + 0.0 base 10
10 base π + 0.0 base 10
30 base π + 0.44482644031997864 base 10
231 base π + 0.8422899173517213 base 10
2312 base π + 0.6461318165449161 base 10
23122 base π + 0.029882968108176033 base 10
231220 base π + 0.0938801130760924 base 10
2312130 base π + 0.7397595138779653 base 10
23121302 base π + 0.3240230542211062 base 10
231213021 base π + 0.017948446735832846 base 10
2312130210 base π + 0.05638670840988885 base 10
23121302100 base π + 0.17714406890720072 base 10
231213021000 base π + 0.5565145054551264 base 10
2312130133130 base π + 0.6366321966964654 base 10
23121301331302 base π + 3.9032618162071486e-05 base 10
231213013313020 base π + 0.00012262302157861615 base 10
2312130133123211 base π + 0.24905356925301847 base 10
23121301331232110 base π + 0.7824248909895828 base 10
231213013312321102 base π + 0.4580601707952492 base 10
2312130133123211021 base π + 0.4390387422112354 base 10
>>> convert_pi_to_dec('2312130133123211021 base π + 0.4390387422112354 base 10')
2791563949.5978436
>>> convert_pi_to_dec('10000000000000000000 base π + 0.0 base 10')
2791563949.5978436
>>>

还显示了最后两个字符串如何等效，但输出应该是第二个字符串的形式。我发现 10000000000000000000 base π 和 2312130133123211021 base π 之间的区别是 0.4390387422112354 base 10 很有趣，但这种区别对表示有很大影响。输出应该如下所示。

1 base π + 0.0 base 10
10 base π + 0.0 base 10
100 base π + 0.0 base 10
1000 base π + 0.0 base 10
10000 base π + 0.0 base 10
100000 base π + 0.0 base 10
1000000 base π + 0.0 base 10
10000000 base π + 0.0 base 10
100000000 base π + 0.0 base 10
1000000000 base π + 0.0 base 10
10000000000 base π + 0.0 base 10
100000000000 base π + 0.0 base 10
1000000000000 base π + 0.0 base 10
10000000000000 base π + 0.0 base 10
100000000000000 base π + 0.0 base 10
1000000000000000 base π + 0.0 base 10
10000000000000000 base π + 0.0 base 10
100000000000000000 base π + 0.0 base 10
1000000000000000000 base π + 0.0 base 10
10000000000000000000 base π + 0.0 base 10

我有什么遗漏的吗，有没有解决这个问题的方法，或者这应该被认为是愚蠢的差事？

【问题讨论】：

所以基本上，您希望余数始终为 0？那为什么要使用divmod？只需使用除法！
半严肃的回答：使用分段函数func = lambda x,y: (pi,0) if x== pi*pi and y == pi else divmod(x,y)。如果这还不够，请显示更多示例输入/输出。很难从一个数据点推导出方程。
请查看编辑。您的函数适用于第一个示例。然而，这是我追求的第二个应用程序。

标签： python floating-point division modulus

【解决方案1】：

您正在寻找一种算法来确定浮点数的non-integer base 表示。

维基百科由于 Rényi 和 Frougny 描述了一个贪心算法；这是一个实现的尝试：

from math import log, floor
def expansion(x, b):
    k = int(floor(log(x) / log(b)))
    d, r = divmod(x / float(b ** k), 1)
    digits = [int(d)]
    for _ in range(k):
        d, r = divmod(b * r, 1)
        digits.append(int(d))
    def rest(b, d, r):
        while r:
            d, r = divmod(b * r, 1)
            yield int(d)
    return digits, rest(b, d, r)

这给出了字典顺序的初始扩展；你可以通过一点点摆弄来获得字典顺序的终端扩展：

def expansion(x, b, greedy=True):
    if not greedy:
        m = (floor(b) / (b - 1)) - 1
    k = int(floor(log(x) / log(b)))
    d, r = divmod(x / float(b ** k), 1)
    if not greedy and r < m:
        d, r = d - 1, r + 1
    digits = [int(d)]
    for _ in range(k):
        d, r = divmod(b * r, 1)
        if not greedy and r < m:
            d, r = d - 1, r + 1
        digits.append(int(d))
    def rest(d, r):
        while r:
            d, r = divmod(b * r, 1)
            if not greedy and r < m:
                d, r = d - 1, r + 1
            yield int(d)
    return digits, rest(d, r)

不幸的是，这仍然行不通，因为 OP 的扩展在第一个数字上是非贪婪的，但在最后一个数字上是贪婪的。

【讨论】：

我不确定贪心算法是否会给 OP 带来预期的结果。 ISTM 第一个贪心数字是 2，代表不会终止。
这很复杂，因为非整数基表示通常不是唯一的。就像以10为底的1也可以表示为.9999...，示例数字，以123123为底的数字与2000000几乎相同（尽管后者省略了一些小数位） .贪心算法更喜欢后一种形式，我不确定有没有办法获得替代版本。
示例编号很好。执行convert_dec_to_pi(1 * pi ** 5 + 2 * pi ** 4 + 3 * pi ** 3 + 1 * pi ** 2 + 2 * pi ** 1 + 3 * pi ** 0 + 0.0000000000001) 返回'123123 base π + 3.597122599785507e-14 base 10'，这已经足够接近了。然而，令人失望的是，运行 convert_dec_to_pi(pi * pi) 并没有返回 '100 base π + 0.0 base 10'。

【解决方案2】：

认识到浮点运算在定义上是不精确的。 pi*pi 之类的操作不能保证等于数学常数 π^2（就此而言，math.pi 仅与“可用精度”一样准确——这意味着它也不是正确的值）。因此，实际上不可能对浮点数进行操作，将它们视为实数。

一般的解决方案是检查与某个 epsilon 值的距离，但这有明显的局限性。您最好重新检查您的基本要求（为什么需要实数精度？）并尝试从不同的方向解决问题。

对于您描述的示例，为什么您需要实际使用 π 的值？你能把 π 的实际计算留到最后，只对系数进行运算吗？

例如，直接存储列表[3, 2, 1, 3, 2, 1]，并使用它们是系数的隐式契约进行操作和转换，然后定义如下内容：

toFloat(ls,mult):
  pow = 0
  ret = 0
  for coef in ls:
    ret += coef * mult**pow
    pow += 1
  return ret

作为打印前的最后一步。更好的是，您可以将这种行为封装在一个类中（我敢打赌有人曾经这样做过）并使__str__() 执行toFloat() 的行为，以便显示您的对象为您提供最精确的您可以获得的价值。

【讨论】：

我不认为提问者的问题与浮点精度有关（尽管这可能会在稍后出现，如果先解决概念问题）。
嗯，我认为期望 divmod(pi * pi, pi) == (pi, 0) 为 True 是对浮点运算的误解。第二个例子同样想要应用然后去除 π 的幂，就好像它们是可以随意添加和删除的精确值一样。致 OP - 如果这不是一个有用的答案，请为您的问题添加更多背景故事。
问题是，我没有系数。我的愿望是提取系数。第二次编辑应该有助于更好地解释情况。
divmod(pi * pi, pi) 在我的 Python 2 和 3 上是 (3.0, 0.44482644031997864)。不是精度问题。 "For floating point numbers the result is (q, a % b), where q is usually math.floor(a / b) ..."

【解决方案3】：

这个非常简单，似乎比 OP 更有效。我认为结果中的缺陷与精度有关：

import math
import struct
import os
from decimal import Decimal, getcontext

getcontext().prec = 1000

def digits_base_b(n, b):
    n = Decimal(n)
    b = Decimal(b)
    digits = {}
    while n >= b:
        exp = int(math.log(n, b))
        digit = int(n/b**exp)
        digits[exp] = digit
        n -= digit*b**exp
    return digits, n # n is less than b**1, idk how you want to handle those

def digits_2_str(digits, base):
    exps = sorted(digits, reverse=True)
    result = []
    format_spec = '%d*'+base+'^%d'
    for exp in exps:
        result.append(format_spec % (digits[exp], exp))
    return ' + '.join(result)

pi = Decimal(
'3.14159265358979323846264338327950288419716939937510'
'58209749445923078164062862089986280348253421170679'
'82148086513282306647093844609550582231725359408128'
'48111745028410270193852110555964462294895493038196'
'44288109756659334461284756482337867831652712019091'
'45648566923460348610454326648213393607260249141273'
'72458700660631558817488152092096282925409171536436'
'78925903600113305305488204665213841469519415116094'
'33057270365759591953092186117381932611793105118548'
'07446237996274956735188575272489122793818301194912'
'98336733624406566430860213949463952247371907021798'
'60943702770539217176293176752384674818467669405132'
'00056812714526356082778577134275778960917363717872'
'14684409012249534301465495853710507922796892589235'
'42019956112129021960864034418159813629774771309960'
'51870721134999999837297804995105973173281609631859'
'50244594553469083026425223082533446850352619311881'
'71010003137838752886587533208381420617177669147303'
'59825349042875546873115956286388235378759375195778'
'18577805321712268066130019278766111959092164201989'
)

if __name__ == '__main__':
    random_float = lambda: struct.unpack('d', os.urandom(8))[0]
    x = random_float()
    while x < pi: # some floats are no good, i've only tested with positives
        x = random_float()

    digits, leftover = digits_base_b(x, pi)
    print x, '='
    print digits_2_str(digits, u'\u03C0')

    for i in range(20):
        digits, leftover = digits_base_b(pi**i, pi)
        print float(pi**i), '=', digits_2_str(digits, u'\u03C0'), '+', float(leftover)

更新我从互联网上得到了 pi 的前一千位，并使用了十进制。十进制并少了一些错误，但仍然有几个。因此，我相信差异与精度有关。此外，随着精度的提高，计算所需的时间也会急剧增加。

【讨论】：