DouKaku? 4783 shiro: 引き算のみで。 Z80マシン語で遊...(最大公約数（除算禁止）)

Comment detail

引き算のみで。

Z80マシン語で遊んでいたころ、除算を使わないでもこれでgcdが計算できることに気づいて大発見をしたような気分になりました。後でよく考えたらユークリッドの互除法を非効率にやっているだけでした。

隣り合うフィボナッチ数の場合は実は減算法でもあまり効率が悪くなりません。常に差がひとつ前のフィボナッチ数になるわけですから。

gosh> (time (gcd-nodiv (fib 2000) (fib 1999)))
;(time (gcd-nodiv (fib 2000) (fib 1999)))
; real   0.002
; user   0.000
; sys    0.000
1

ループ回数は互除法でやった場合と同じになります。

gosh> (gcd-nodiv/cnt (fib 2000) (fib 1999) 1)
1
1998

実はこの引数の場合は、組み込みのgcd (互除法使用) の方が遅かったりします。bignumの除算が重いんじゃないかな。

gosh> (time (gcd (fib 2000) (fib 1999)))
;(time (gcd (fib 2000) (fib 1999)))
; real   0.014
; user   0.010
; sys    0.000
1

(define (gcd-nodiv m n)
  (cond [(= n m) n]
        [(or (= m 1) (= n 1)) 1]
        [(> m n) (gcd-nodiv n (- m n))]
        [else    (gcd-nodiv m (- n m))]))

;;
;; 以下は測定用
;;
;;ループ回数のカウント
(define (gcd-nodiv/cnt m n cnt)
  (cond [(= n m) (values n cnt)]
        [(or (= m 1) (= n 1)) (values 1 cnt)]
        [(> m n) (gcd-nodiv/cnt n (- m n) (+ cnt 1))]
        [else    (gcd-nodiv/cnt m (- n m) (+ cnt 1))]))

;;メモワイズ版fib
(define fib
  (let1 memo (make-hash-table 'eqv?)
    (lambda (n)
      (cond [(<= n 2) 1]
            [(hash-table-get memo n #f)]
            [else (let1 val (+ (fib (- n 1)) (fib (- n 2)))
                    (hash-table-put! memo n val)
                    val)]))))

ちょっと補足。

> ループ回数は互除法でやった場合と同じになります。 

というのはあくまで引数が隣り合うフィボナッチ数である場合の話です。

出題者です.

>隣り合うフィボナッチ数の場合は実は減算法でもあまり効率が悪くなりません。常に差がひとつ前のフィボナッチ数になるわけですから。

全くもってその通りです. 1024bitと512bitの素数とかでやった方が良かったようでorz

この回答が破綻するのは一方が他方よりはるかに大きな場合ですね。m >> n で m = q*n + rのとき n >> r になるような組み合わせを続けてゆくと減算のみによる方法では悲惨なことに。

減算を繰り返して割り算の代わりにする方法だと差が大きいときに問題になる、ということで 「文字列の長さをみて小さい側を10**n倍する」 という荒技を使ってみました。

# Screencast-o-maticで撮りながら作ったのにアップロード中にブラウザが落ちて消滅…orz

def gcd(x, y):
    if x == y: return x
    if x < y: y, x = x, y

    lx = len(str(x))
    ly = len(str(y))
    diff = lx - ly

    z = abs(y * 10 ** diff - x)
    z2 = abs(y * 10 ** (diff + 1) - x)
    z = min(z, z2)
    if diff > 0:
        z2 = abs(y * 10 ** (diff - 1) - x)
        z = min(z, z2)

    # zはxとyと公約数を共有する
    if z < y:
        return gcd(y, z)
    return gcd(z, y)