SPOJ`CUBERT`中的错误答案

Question

我在SPOJ上找到了解决这个问题的错误答案.

该问题要求计算整数的立方根(可以长达150位),并输出截断最多10位小数的答案.
它还要求将答案模10中所有数字的总和计算为"校验和"值.

这是确切的问题陈述:

您的任务是计算给定正整数的立方根.我们不记得为什么我们需要这个,但它与公主,一个年轻的农民,亲吻和一半的王国(一个巨大的,我们可以向你保证)有一些共同之处.

编写一个程序来解决这个关键任务.

输入

输入以包含单个整数t <= 20(测试用例数)的行开头.t测试案例如下.

下一行包含最多150个十进制数字的大正整数.每个数字都在输入文件的单独行中.输入文件可能包含空行.数字可以在空格之前或之后,但没有行超过255个字符.

产量

对于输入文件中的每个数字,程序应输出一行,其中包含由单个空格分隔的两个值.第二个值是给定数字的立方根,在第10个小数位后截断(不是舍入!).第一个值是立方根的所有打印数字的校验和,计算为模10的打印数字之和.

例

输入:
5
1

8

1000

2 33076161

输出:
1 1.0000000000
2 2.0000000000
1 10.0000000000
0 1.2599210498
6 321.0000000000

这是我的解决方案:

from math import pow


def foo(num):
    num_cube_root = pow(num, 1.0 / 3)
    # First round upto 11 decimal places
    num_cube_root = "%.11f" % (num_cube_root)
    # Then remove the last decimal digit
    # to achieve a truncation of 10 decimal places
    num_cube_root = str(num_cube_root)[0:-1]

    num_cube_root_sum = 0
    for digit in num_cube_root:
        if digit != '.':
            num_cube_root_sum += int(digit)
    num_cube_root_sum %= 10

    return (num_cube_root_sum, num_cube_root)


def main():
    # Number of test cases
    t = int(input())
    while t:
        t -= 1
        num = input().strip()
        # If line empty, ignore
        if not num:
            t += 1
            continue

        num = int(num)
        ans = foo(num)
        print(str(ans[0]) + " " + ans[1])


if __name__ == '__main__':
    main()

它适用于示例案例:现场演示.

任何人都可以告诉这个解决方案有什么问题吗？

Answer 1

您的解决方案有两个问题,都与使用浮点运算有关.第一个问题是Python float只携带大约16位有效精度的十进制数字,所以只要你的答案需要超过16位有效数字左右(所以在点之前超过6位,之后超过10位),你就是获得正确的尾随数字的希望很小.第二个问题更微妙,甚至影响到小的价值n.这就是你的舍入到11位十进制数然后丢弃最后一位数的方法由于双舍入而遭受潜在错误.举个例子,拿n = 33.立方根n,到20位小数左右,是:

3.20753432999582648755...

如果在该点之后四舍五入到11位,那么你最终会得到

3.20753433000

现在丢弃最后一位数3.2075343300,这不是你想要的.问题是,小数到11位小数可能会影响到第11位数左边的数字.

那么你能做些什么来解决这个问题呢？好吧,你可以完全避免浮点数并将其减少为纯整数问题.我们需要一些整数n到10个小数位的立方根(将最后一个位置向下舍入).这相当于将立方根计算为10**30 * n最接近的整数,再次向下舍入,然后将结果除以10**10.所以这里的基本任务是计算任何给定整数的立方根的底限n.我无法找到关于计算整数多维数据集根目录的任何现有Stack Overflow答案(在Python中更少),所以我认为值得详细说明如何这样做.

计算整数的整数根证明是非常容易的(借助一点点数学).有各种可能的方法,但是一种既有效又易于实现的方法是使用Newton-Raphson方法的纯整数版本.在实数上,牛顿求解方程的方法x**3 = n近似于x立方根n,并迭代以返回改进的近似.所需的迭代是:

x_next = (2*x + n/x**2)/3

在实际情况中,您将重复迭代,直到达到某个所需的容差.事实证明,在整数,本质上是相同的迭代工作,并与合适的退出条件,它会给我们正好正确的答案(无需公差).整数情况下的迭代是:

a_next = (2*a + n//a**2)//3

(注意使用地板分割算子//代替/上面通常的真正除法运算符.)数学上,a_next正好是地板(2*a + n/a**2)/3.

这是基于此迭代的一些代码:

def icbrt_v1(n, initial_guess=None):
    """
    Given a positive integer n, find the floor of the cube root of n.

    Args:
        n : positive integer
        initial_guess : positive integer, optional. If given, this is an
            initial guess for the floor of the cube root. It must be greater
            than or equal to floor(cube_root(n)).

    Returns:
        The floor of the cube root of n, as an integer.
    """
    a = initial_guess if initial_guess is not None else n
    while True:
        d = n//a**2
        if a <= d:
            return a
        a = (2*a + d)//3

一些例子使用:

>>> icbrt_v1(100)
4
>>> icbrt_v1(1000000000)
1000
>>> large_int = 31415926535897932384626433
>>> icbrt_v1(large_int**3)
31415926535897932384626433
>>> icbrt_v1(large_int**3-1)
31415926535897932384626432

icbrt_v1我们很快就会解决一些烦恼和效率低下的问题.但首先,简要说明上述代码的工作原理.请注意,我们从一个初始猜测开始,假设它大于或等于立方根的底面.我们将展示这个属性是一个循环不变量:每次我们到达while循环的顶部时,a至少是这样floor(cbrt(n)).此外,每次迭代产生的值a严格小于旧的值,因此我们的迭代最终会收敛到floor(cbrt(n)).为了证明这些事实,请注意,当我们进入while循环时,有两种可能性:

案例1. a严格大于立方根n.然后a > n//a**2,代码进入下一次迭代.写a_next = (2*a + n//a**2)//3,然后我们有:

• a_next >= floor(cbrt(n)).这是从以下事实(2*a + n/a**2)/3是至少的立方根n,这反过来从如下AM-GM不等式适用于a,a和n/a**2:这三个量的几何平均数是准确的立方根n,所以算术平均值必须在至少是立方根n.因此,我们的循环不变量将保留用于下一次迭代.

• a_next < a:因为我们假定a比立方根较大,n/a**2 < a以及接下去(2a + n/a**2) / 3是小比a,从而获得一种floor((2a + n/a**2) / 3) < a.这保证了我们在每次迭代时都能在解决方案上取得进展.

情况2. a小于或等于立方根n.然后a <= floor(cbrt(n)),从上面建立的循环不变量我们也知道a >= floor(cbrt(n)).所以我们完成了:a是我们追求的价值.此时,while循环退出a <= n // a**2.

上面的代码有几个问题.首先,与初始猜测开始n是低效的:该代码将花费其最初的几个迭代(粗略地)分割的当前值a由3每一次,直到它进入溶液的附近.初始猜测(以及在Python中可以轻松计算的一个)的更好选择是使用超过立方根的第一个2的幂n.

initial_guess = 1 << -(-n.bit_length() // 3)

更好的是,如果n小到足以避免溢出,则使用浮点算法来提供初始猜测,例如:

initial_guess = int(round(n ** (1/3.)))

但这带来了我们的第二个问题:我们的算法的正确性要求初始猜测不小于实际的整数立方根,并且随着n变大,我们无法保证对于基于浮点数的initial_guess上述(尽管足够小)n, 我们可以).幸运的是,有一个非常简单的修复:对于任何正整数a,如果我们执行单次迭代,我们总是得到一个至少是一个值floor(cbrt(a))(使用我们上面使用的相同AM-GM参数).因此,我们所要做的就是在开始测试收敛之前至少执行一次迭代.

考虑到这一点,这是上述代码的更高效版本:

def icbrt(n):
    """
    Given a positive integer n, find the floor of the cube root of n.

    Args:
        n : positive integer

    Returns:
        The floor of the cube root of n, as an integer.
    """
    if n.bit_length() < 1024:  # float(n) safe from overflow
        a = int(round(n**(1/3.)))
        a = (2*a + n//a**2)//3  # Ensure a >= floor(cbrt(n)).
    else:
        a = 1 << -(-n.bit_length()//3)

    while True:
        d = n//a**2
        if a <= d:
            return a
        a = (2*a + d)//3

并且icbrt在手中,很容易将所有内容放在一起以计算立方根到十位小数.这里,为简单起见,我将结果输出为字符串,但您可以轻松地构造Decimal实例.

def cbrt_to_ten_places(n):
    """
    Compute the cube root of `n`, truncated to ten decimal places.

    Returns the answer as a string.
    """
    a = icbrt(n * 10**30)
    q, r = divmod(a, 10**10)
    return "{}.{:010d}".format(q, r)

示例输出:

>>> cbrt_to_ten_places(2)
'1.2599210498'
>>> cbrt_to_ten_places(8)
'2.0000000000'
>>> cbrt_to_ten_places(31415926535897932384626433)
'315536756.9301821867'
>>> cbrt_to_ten_places(31415926535897932384626433**3)
'31415926535897932384626433.0000000000'