Burnikel 和 Ziegler 算法 2 的 Eiffel 实现中的错误

Posted 2023-02-16

【中文标题】Burnikel 和 Ziegler 算法 2 的 Eiffel 实现中的错误

【英文标题】：Error in Eiffel implementation of Burnikel and Ziegler algorithm 2

【发布时间】：2017-04-21 21:11:49

【问题描述】：

我需要另一双眼睛来告诉我埃菲尔实现的 Burnikel 和 Ziegler 分区有什么问题，特别是“算法 2 - 3n/2n”。 Eiffel 特征如下所示。 “like Current”类型是 ARRAYED_LIST [NATURAL_8]。换句话说，该实现使用包含 8 位值的数字（即肢体），因此数字以 256 为基数。随后是失败调用的手动跟踪。 （抱歉，参数太大，但我无法用较短的值重现错误。）在这种情况下，执行遵循步骤 3b。

这就是问题所在。该算法似乎适用于第 5 步，其中余数“r”以比除数更多的数字结束。我相信错误出现在步骤 3b 中，可能是调用了“应该”提供一个值为“Beta^n - 1”的功能“ones”。 （也许我不明白 B&Z 的“Beta^n”符号。

这里是埃菲尔代码：

three_by_two_divide (a, a3, b: like Current): TUPLE [quot, rem: like Current]
        -- Called by `two_by_one_divide'.  It has similar structure as
        -- `div_three_halves_by_two_halfs', but the arguments to this
        -- function have type JJ_BIG_NATURAL instead of like `digit'.
        -- See Burnikel & Zieler, "Fast Recursive Division", pp 4-8,
        -- Algorithm 2.
    require
        n_not_odd: b.count >= div_limit and b.count \\ 2 = 0
        b_has_2n_digits: b.count = a3.count * 2
        a_has_2n_digits: a.count = a3.count * 2
    local
        n: INTEGER
        a1, a2: like Current
        b1, b2: like Current
        tup: TUPLE [quot, rem: like Current]
        q, q1, q2, r, r1: like Current
        c, d: like Current
    do
        n := b.count // 2
            -- 1) Split `a'
        a1 := new_sub_number (n + 1, a.count, a)
        a2 := new_sub_number (1, n.max (1), a)
            -- 2) Split `b'.
        b1 := new_sub_number (n + 1, b.count, b)
        b2 := new_sub_number (1, n.max (1), b)
            -- 3) Distinguish cases.
        if a1 < b1 then
                -- 3a) compute Q = floor ([A1,A2] / B1 with remainder.
            if b1.count < div_limit then
                tup := school_divide (a, b1)
            else
                tup := two_by_one_divide (a, b1)
            end
            q := tup.quot
            r1 := tup.rem
        else
                -- 3b) Q = beta^n - 1 and ...
            q := ones (n)
                -- ... R1 = [A1,A2] - [B1,0] + [0,B1] = [A1,A2] - QB1.
            r1 := a + b1
            if n > 1 then
                b1.shift_left (n)
            else
                b1.bit_shift_left (zero_digit.bit_count // 2)
            end
            r1.subtract (b1)
        end
            -- 4) D = Q * B2
        d := q * b2
            -- 5) R1 * B^n + A3 - D.  (The paper says "a4".)
        r1.shift_left (n)
        r := r1 + a3 - d
            -- 6) As long as R < 0, repeat
        from
        until not r.is_negative
        loop
            r := r + b
            q.decrement
        end
        check
            remainder_small_enough: r.count <= b.count
                -- because remainder must be less than divisor.
        end
        Result := [q, r]
    ensure
   --   n_digit_remainder: Result.rem.count = b.count // 2
        quotient_has_correct_count: Result.quot.count <= b.count // 2
    end

在跟踪中，箭头指向我认为不好的线，但我不知道如何处理它。这是跟踪：

three_by_two_divide (a = [227,26,41,95,169,93,135,110], 
                     a3 = [92,164,19,39],
                     b =  [161,167,158,41,164,0,0,0]) 

    n := b.count // 2 = 4
        -- 1) Split `a'.
    a1 := new_sub_number (n + 1, a.count, a) = [227,26,41,95]
    a2 := new_sub_number (1, n.max (1), a) = [169,93,135,110]
        -- 2) Split `b'.
    b1 := new_sub_number (n + 1, b.count, b) = [161,167,158,41]
    b2 := new_sub_number (1, n.max (1), b) = [164,0,0,0]
        -- 3b) Q = beta^n -1 and ...
--> q := ones (4) = [255,255,255,255]          <-- Is this the error?
        -- ... R1 = [A1,A2] - [B1,0] + [0,B1].
    r1 := a + b1 = [227,26,41,96,75,5,37,151]
    b1.shift_left (n) = [161,167,158,41,0,0,0,0]                        
    r1.subtract (b1) = [65,114,139,55,75,5,37,151]
    d := q * b2 = [163,255,255,255,92,0,0,0]
    r1.shift_left (n) = [227,25,135,184,172,220,37,151,0,0,0,0]   -- too big!
    r := r1 + a3 - d -= [227,25,135,184,8,220,37,152,0,164,19,39] -- too big!

我知道这很长，但感谢任何帮助。

【参考方案1】：

我建议在执行r1.shift_left (n) 之前检查r1 = [65,114,139,55,75,5,37,151] 是否仍然相同。有两种选择：

d := q * b2 会影响r1，但它不应该。很可能存在一些别名，即 r1 与其他一些更新的变量有别名，应该删除此别名。 r1 在d := q * b2 之后仍然相同。问题在于 shift_left 未能（重新）初始化某些数据或使用不应使用的全局数据。

【讨论】：

感谢 Alexander，没有别名，也没有全局数据。我想我已经通过避免调用“ones”来解决这个问题。 Burnikel & Ziegler 说，对于 A 除以 B，“让 A = Bβ^n，则得到一个 new_a := A - B 并记住商以 1 开头。此时，A 小于 B，所以按照 B&Z 的说法进行。我还有另一个问题，我将在一个新问题中发布。谢谢。

我真的看不出用r1.shift_left (n) 将[65,114,139,55,75,5,37,151] 左移4 位会得到[227,25,135,184,172,220,37,151,0,0,0,0]。如果你能解释这种行为，这将有助于理解实现。