c++ - x(x-1)/2 = c 的快速整数解

Question

给定一个非负整数 c，我需要一个高效的算法来找到最大的整数 x 使得

x*(x-1)/2 <= c

同样，我需要一种高效且可靠准确的算法来计算:

x = floor((1 + sqrt(1 + 8*c))/2)        (1)

为了明确起见，我将这个问题标记为 C++，因此答案应该是用该语言编写的函数。您可以假设 c 是一个无符号的 32 位 int。

此外，如果您可以证明 (1)(或涉及浮点运算的等价表达式)总是给出正确的结果，那也是一个有效的答案，因为现代处理器上的浮点运算可能比整数算法更快。

Answer 1

如果您愿意假设 IEEE double 对所有运算(包括平方根)进行正确舍入，那么您编写的表达式(加上转换为 double)会给出所有输入的正确答案。

这是一个非正式的证明。由于 c 是一个 32 位无符号整数，正在转换为具有 53 位有效数的浮点类型，因此 1 + 8*(double)c 是精确的，并且sqrt(1 + 8*(double)c) 正确四舍五入。 1 + sqrt(1 + 8*(double)c) 精确到一个 ulp 以内，因为最后一项小于 2**((32 + 3)/2) = 2**17.5表示后一项最后一位的单位小于1，因此(1 + sqrt(1 + 8*(double)c ))/2 精确到一个 ulp，因为除以 2 是精确的。

最后一项业务是地板。这里的问题案例是当 (1 + sqrt(1 + 8*(double)c))/2 被四舍五入为整数时。当且仅当 sqrt(...) 舍入为奇数时才会发生这种情况。由于 sqrt 的参数是一个整数，对于正奇数 z，最坏的情况看起来像 sqrt(z**2 - 1)，我们绑定(bind)

z - sqrt(z**2 - 1) = z * (1 - sqrt(1 - 1/z**2)) >= 1/(2*z)

通过泰勒展开。由于 z 小于 2**17.5，因此与最接近整数的差距至少为 1/2**18.5 结果幅度小于 2**17.5，这意味着此错误不可能由正确舍入的 sqrt 引起。

采用Yakk的化简，我们可以写

(uint32_t)(0.5 + sqrt(0.25 + 2.0*c))

无需进一步检查。