LoginSignup
2
2

More than 1 year has passed since last update.

あのアルゴリズムはどこ?[Scheme(Gauche)版]

Last updated at Posted at 2021-11-24

0.はじめに

本家ッ!

あのアルゴリズムはどこ? Pythonを使用してAtCoderの緑色や水色を目指す方に、30以上のアルゴリズムスニペットと100問以上の問題(ACコード付き)を紹介!

説明ッ!

  • アルゴリズムのスニペット(?)の一覧はScheme版も欲しいなと思ったから。
    • Project EulerとAtCoderにScheme(Gauche)で参加しているから。
    • AtCoder最近再開して基本パターンを押さえておく重要性を感じたから。
    • あわよくばACLのScheme版を作るためのステップアップにしたいから。
  • Schemeの勉強をしたいから。
    • WebにSchemeの記事が少なすぎて他言語に比べて勉強しづらいと感じているので、ここをアウトプット&レスポンスの場にしたいから。
      • 特にAtCoder+Schemeとなると絶望的に記事が少ない(CommonLispも……)。
      • 戯れにPythonと比較してみたが比較するとかいうレベルではなくてマジ草。
    • Schemeの表現力、数値計算向きなところがプロコンでも活きると思っているから。
      • Pythonで書きづらくてもSchemeなら書きやすかったり、逆のケースもあると思うのでそういった違いを自分でも確認したいから。
    • コレみたいな感じでScheme版が欲しいと思ったから。

自分はSchemeぢからレベル3~5のあたりを低い意識でフラフラしている感じなので、間違い/より良い書き方があるかもしれませんが、そういったご指摘大歓迎ですので見つけられた方はコメントください。

なお、個人の勉強目的で書かれているので元のアルゴリズムをこの記事だけで理解しようとするのは難しいと思います。本家やリンク先の記事でアルゴリズムを理解してからか、もともと理解している人が実装を確認する、Schemeでの書き方を楽しむ目的で読むのが良いのではないでしょうか。

……実はまだ途中だけどまあいいや、公開!

環境ッ!

基本的に本記事では、

  • Windows 10
  • Gauche 0.9.9

を開発環境としています。Gaucheのビルトイン機能の、R*RS標準から外れた関数を使っている可能性がありますので他処理系で使う場合はご注意ください。

Gaucheのドキュメント参照にはできる限りオンラインマニュアルにも記載されているリンクを使用していますが、抜け漏れありましたらこちらもご指摘いただけますと幸いです。

http://practical-scheme.net/gauche/man/?l=ja&p=トピック名

余談ですがこの検索URLはブラウザの検索エンジン設定に登録しておくと便利です。

1.目次

(※本家のコピペです)
2.除算の切り上げ、切り下げ
3.約数列挙
4.素因数分解
5.エラトステネスの篩
6.BIT全探索
7.Union-Find
8.クラスカル法(最小全域木)
9.リスト(配列)の中身を連結・結合して文字列に
10.多次元配列のソート
11.defaultdict(辞書)
12.カウンター
13.ModInt
14.逆元
15.階乗
16.コンビネーション
17.ある値の累乗の何度も使用する場合
18.最大公約数
19.中国剰余定理
20.DFS
21.二分探索(bisect)
22.二分探索(その他)
23.最長増加部分列(LIS)
24.累積和(合わせて、いもす法)
25.短めな順列の作成
26.二次元累積和
27.ダイクストラ法
28.ワーシャルフロイド法
29.セグメント木(セグ木)
30.転倒数
31.桁DP
32.ダブリング
33.最長共通部分列
34.BitDP
35.10進数を2進数、8進数、16進数に変換
36.10進数とn進数の変換
37.正規表現
38.キュー
39.優先度付きキュー
40.しゃくとり法

2.除算の切り上げ、切り下げ

Scheme(Gauche)では除算の切り上げ、切り下げで関数呼び出しが遅いとかモジュールインポートしなければならないとかそういったことはあまり考えなくて良く、floorceilingを呼んでやればよい。

(floor (/ 4 3))
(ceiling (/ 4 3))
(truncate (/ 4 3))
(round (/ 4 3))

ただSchemeには正確数、不正確数という概念があるので、いわゆるIntを返すことを常に期待しているのなら->exactのついた関数を使ったほうがよい。

(floor->exact 2.3)
(ceiling->exact 2.3)
(truncate->exact 2.3)
(round->exact 2.3)

久しぶりにc++のコードを読んだらInt型の/演算子適用結果は小数点以下が切り捨てになるという事実に遭遇してびっくりしたが、言語によって「常識」は異なる……。

3.約数列挙、4.素因数分解

これはもうそのままの関数がmath.primeモジュールにあって、

(use math.prime)
(naive-factorize 840)

と書ける。約数だけにしたいなら重複排除の関数があるのでそれを呼べばよい。

(use math.prime)
(delete-duplicates (naive-factorize 840))

計算方法も同じようなので計算量も変わらないだろう。再実装しても良いがProject Eulerをやっていると素数周りの計算はもうおなか一杯になってしまうのでパス……。何か数学的な制約を使って計算時間を短縮できる問題があったらその時に書こうかな。

余談だが楔数の問題は、2から「Nを6で割った数」を超えない素数までの素数リストを作って、3重ループの再帰で回し、Nを超えた端から探索を打ち切るようにすれば効率が良さそうだ。combinationsを一瞬考えたが効率悪いだろうな……。MojaCoderはSchemeでの提出ができないので悲しい。

5.エラトステネスの篩

昔Project Eulerで書いたコードがあったので、今の知識で書き直してみた。

(define (sieve n)
  (let [(limit (sqrt n))]
    (let loop [(p-list (iota (- n 1) 2))]
      (cond [(null? p-list) '()]
            [(< limit (car p-list)) p-list]
            [else (cons (car p-list)
                        (loop (filter-map (^(x) (if (= 0 (modulo x (car p-list))) #f x))
                                          (cdr p-list))))]))))

(※なお上記の「^」はGaucheにおけるlamndaの簡易記法。またSchemeのBooleanは「#t」と「#f」。)
リストだとインデクスで悩まなくて良い。

実際に解いた問題

6.BIT全探索

0から(- (expt 2 n) 1)までの数を2進数で表現したい。表現方法を考えたとき、とにかく0と1が並んだリストにしたいので、一旦整数のリストを作ってからその整数を文字リストに分解してやるのが良さそうだ。
というわけでこうなる。

(define (make-bit-list bit-len)
  (map (^(l) (if (= bit-len (length l))
               l
               (append (make-list (- bit-len (length l)) #\0) l)))
       (map (^(x) (string->list (number->string x 2)))
            (iota (expt 2 bit-len) 0))))

動作例

gosh> (make-bit-list 3)
gosh> ((#\0 #\0 #\0) (#\0 #\0 #\1) (#\0 #\1 #\0) (#\0 #\1 #\1) (#\1 #\0 #\0) (#\1 #\0 #\1) (#\1 #\1 #\0) (#\1 #\1 #\1))

さらにPython de アルゴリズム(bit全探索)の例題を解くなら、金額のリストとfoldしてやって、金額の合計を計算し、300円以下になるものだけを残してカウントしてやれよい。

(length
  (filter-map (^(x) (if (< 300 (apply + x)) #f x))
              (map (^(bits) (fold (^(x y lis) (if (eqv? #\1 y)
                                                (cons x lis)
                                                (cons 0 lis)))
                                  '()
                                  '(100 200 300)     
                                  bits))
                   (make-bit-list 3))))
  => 5

今回は整数のリストにしてから合計を計算するため、ビット列中0に当たった場合は整数0を返すようにしたが、問題によってビット=0に当たった場合の挙動を後の計算で便利なように変えてやれば良い('()で返すとか#fで返すとか)。

どこまでメモリ使用量を気にしてやれば良いかわからないが、最初にビット列をすべて求めるのが効率悪かったら、ビット列を逐次計算するパターンの方が良いだろう。

;;; 数え上げながら逐次bit表現を作って適用するパターン
(define (convert-to-bit-list n bit-len)
  (let [(bit-list (string->list (number->string n 2)))]
    (if (= bit-len (length bit-list))
      bit-list 
      (append (make-list (- bit-len (length bit-list)) #\0) bit-list))))

動作例

gosh> (convert-to-bit-list 5 3)
gosh> (#\1 #\0 #\1)

'(100 200 300) から選ぶ/選ばないの組み合わせ全パターンのうち300以下になるものの数。

(length
  (filter-map (^(x) (if (< 300 (apply + x)) #f x))
              (let* [(bit-len 3)
                     (limit (expt 2 bit-len))]
                (let loop [(n 0)]
                  (cond [(<= limit n) '()]
                        [else (cons (fold (^(x y lis) (if (eqv? #\1 y)
                                                        (cons x lis)
                                                        (cons 0 lis)))
                                          '()
                                          '(100 200 300)     
                                          (convert-to-bit-list n bit-len))
                                    (loop (+ n 1)))])))))
  => 5

7.Union-Find

ABC177-D問題の@chokudaiさんの解説が一番わかりやすかったのでこれを参考にSchemeで実装した。下のコード例もこの問題の条件に沿うようにしている。
木構造とは言っても限定的な用途なので、データ構造自体はvecotrや配列で事足りる。valueを親ノードへのポインタ(index)として利用するわけだ。rootのvalueはツリーのサイズとして利用するので、ポインタとの区別のためマイナス符号をつけ、最初は-1で埋め、以降unionするたび逐次|value|がツリーのサイズになるように更新してく。

(define *uf-tree* (make-vector 5 -1))

(define (uf-find uf-tree x)
  (cond [(> 0 (vector-ref uf-tree x)) x]
        [else (let [(ret (uf-find uf-tree (vector-ref uf-tree x)))]
                (vector-set! uf-tree x ret)
                ret)]))

(define (uf-union uf-tree x y)
  (let* [(x-root (uf-find uf-tree x))
         (y-root (uf-find uf-tree y))
         (size-x (vector-ref uf-tree x-root))
         (size-y (vector-ref uf-tree y-root))]
    (cond [(= x-root y-root) #f]
          [(> x-root y-root)
                (vector-set! uf-tree y-root (+ size-x size-y))
                (vector-set! uf-tree x-root y-root)
                #t]
          [else (vector-set! uf-tree x-root (+ size-x size-y))
                (vector-set! uf-tree y-root x-root)
                #t])))

(define (uf-size uf-tree x)
  (- (vector-ref uf-tree (uf-find uf-tree x))))

試しに動かしてみるとこうなる。

gosh> *uf-tree*
gosh> #(-1 -1 -1 -1 -1)
gosh> (uf-union *uf-tree* 0 1)
gosh> #t
gosh> (uf-union *uf-tree* 1 2)
gosh> #t
gosh> *uf-tree*
gosh> #(-3 0 0 -1 -1)
gosh> (uf-size *uf-tree* 2)
gosh> 3

全グループについて要素をまとめたリストを返す関数はこう書ける。

(define (all-members-list uf-tree)
  (let [(limit (vector-length uf-tree))
        (root&elem (make-hash-table))]
    (let loop [(n 0)]
      (cond [(<= limit n) (hash-table-values root&elem)]
            [(memv (vector-ref uf-tree n) (hash-table-keys root&elem)) 
             (hash-table-push! root&elem (vector-ref uf-tree n) n)
             (loop (+ n 1))]
            [(= n (uf-find uf-tree n))
             (hash-table-push! root&elem n n)
             (loop (+ n 1))]
            [else
              (hash-table-push! root&elem (vector-ref uf-tree n) n)
              (loop (+ n 1))]))))

動作例

gosh> (all-members-list *uf-tree*)
gosh> ((2 1 0) (4) (3))

findされた結果、直接rootを親に持つようになったノードに対してまたfindするのは効率が悪いように思う(ほんの少しの差かもしれないが……)。ので、

  1. ハッシュテーブルのkeysつまりすでに判明しているroot要素をvalueに持っていないかチェックし持っていればpushしている。
  2. findした結果rootだと判明したノードは自分をkeyとしてpushする。
  3. またfindされるごとにvectorが書き換わるので、最後のelse節ではnのvalueは必ずrootノードを指していることが保証されるので素直にpushしてやる。

members関数が欲しい場合はリストではなくハッシュテーブル自体を返すようにすれば使いやすいように思う。ただ、all-membersもmembersもベクタの全要素を走査しないと結果が出せないのどちらを使っても計算量自体はそんなに変わらないように思える。

余談だが、「-」や「&」を関数名や変数名に自在に使えるあたり他言語プログラマは気持ち悪さを覚えるのではないだろうか。Lisp/Schemeを書いているとそのうちに気にならなくなって、むしろ()を虹色ハイライトしてくれないことを気持ち悪いと思うようになるのだ……。カッコってこんなに綺麗なのにな……(じっとVimを見る)。

8.クラスカル法(最小全域木)

ABC 065-D問題の解説がわかりやすかった。全域木(全ての要素が網羅されている木)でかつ辺の重みの合計を最小にするにはどうすればよいか、という話。
記憶の片隅に置いておかなければならないのは、実際に最小全域木で辿れる個々の2点間の道筋の重みは評価しないということ。つまりどこかよくわからないけどどっかしらに最小重みの辺でつながっていれば良しとされるわけだ。
なので、2点間の辺の重みが小さいものから順に評価し「つながっているもの」リストに入ってなければリストに入れる(=道を造る)、入っていれば入れない(=道は造らない)
、を繰り返していけばよい。この「つながっているもの」リストにUnion-Findを使うというわけ。辺の重みのリストアップ方法は問題によって異なるだろうが、リストを作ったら昇順に評価する、メモとしてUnion-Findを使うと効率が良い、という構造は変わらない……のだろうたぶん。

コードの見通しを良くしたいのでUF木をクロージャから生成できるようにまとめてみた。

;;; Union-Findをclosureにまとめたもの
(define (make-union-find k)
  (define (uf-find uf-tree x)
    (cond [(> 0 (vector-ref uf-tree x)) x]
          [else (let [(ret (uf-find uf-tree (vector-ref uf-tree x)))]
                  (vector-set! uf-tree x ret)
                  ret)]))
  (define (uf-union uf-tree x y)
    (let* [(x-root (uf-find uf-tree x))
           (y-root (uf-find uf-tree y))
           (size-x (vector-ref uf-tree x-root))
           (size-y (vector-ref uf-tree y-root))]
      (cond [(= x-root y-root) #f]
            [(> x-root y-root)
             (vector-set! uf-tree y-root (+ size-x size-y))
             (vector-set! uf-tree x-root y-root)
             #t]
            [else (vector-set! uf-tree x-root (+ size-x size-y))
                  (vector-set! uf-tree y-root x-root)
                  #t])))
  (define (uf-size uf-tree x)
    (- (vector-ref uf-tree (uf-find uf-tree x))))
  (let [(uf (make-vector k -1))]
    (values (^(x) (uf-find uf x))
            (^(x y) (uf-union uf x y))
            (^(x) (uf-size uf x)))))

ABC-065-D問題を解くのは下記のコード。問題の読み込みや前処理がごちゃごちゃしてて長ったらしいが本質的には最後の3行だと思う。UF木をメモ兼判定機として使う。

(define (answer o)
  (format #t "~A~%" o))

; (x座標 座標ID)または(y座標 座標ID)のリストを入力として、座標でソートし、
; (隣り合う2点間の距離 座標ID-1 座標ID-2)のリストを返す
(define (convert-to-distance-info lst)
  (let loop [(l (sort lst))]
    (cond [(null? l) '()]
          [(null? (cdr l)) '()]
          [else (cons (list (- (caadr l) (caar l)) (cadar l) (cadadr l))
                      (loop (cdr l)))])))

(answer
  (let [(N (read))] (read-line) ; read-lineしてEOLを捨てる
    (let-values [((uf-find uf-union uf-size) (make-union-find N))] ; UF木の生成
      ; 問題の座標群の読み込み。分けて書くのも長くなるのでnamed letの初期値としている。
      (let loop [(distances (let loop [(n 0) (x-list '()) (y-list '())]
                              (cond [(<= N n) (sort (append (convert-to-distance-info x-list)
                                                            (convert-to-distance-info y-list)))]
                                    [else (let* [(x (read)) (y (read))] (read-line)
                                            (loop (+ n 1) (cons (list x n) x-list)(cons (list y n) y-list)))])))
                 (ans 0)]
        ; 問題を解くアルゴリズムの表層部分。
        ; 2点間の距離とIDのリストを作ったら、UF木でチェックしながら道を造っていけばよい。
        (cond [(null? distances) ans]
              [(uf-union (cadar distances) (caddar distances)) (loop (cdr distances) (+ ans (caar distances)))]
              [else (loop (cdr distances) ans)])))))

Lisp/Schemeに慣れてない人は眩暈を起こすだろうな……。

9.リスト(配列)の中身を連結・結合して文字列に

基本なので簡単に。

gosh> (string->list "Alice")
gosh> (#\A #\l #\i #\c #\e)
gosh> (list->string '(#\A #\l #\i #\c #\e))
gosh> "Alice"

文字列はさておき数字の場合。Project Eulerでは整数を各桁に分解して1つずつ評価する問題が頻出するので、n桁の整数をn要素の整数リストに分解する関数を作って使いまわしている。

; 数字を文字リストに変換する
(define (number->list n)
  (string->list (number->string n)))  
; 文字リストを数字に変換する
(define (list->number lis)
  (string->number (list->string lis)))
; 桁を逆順に読んだ数字に変換する
(define (reverse-number n)
  (list->number (reverse (number->list n)))) 

; 数を1桁数字リストに変換する
(define (integer->intlist n)
  (map (^(c) (- (char->integer c) 48)) (number->list n)))
; 1桁数字リストを数に変換する
(define (intlist->integer lis)
  (list->number (map (^(x) (integer->char (+ x 48))) lis)))
gosh> (number->list 12345)
gosh> (#\1 #\2 #\3 #\4 #\5)
gosh> (list->number '(#\1 #\2 #\3 #\4 #\5))
gosh> 12345
gosh> (reverse-number 12345)
gosh> 54321
gosh> (integer->intlist 12345)
gosh> (1 2 3 4 5)
gosh> (intlist->integer '(1 2 3 4 5))
gosh> 12345

1桁数字リストをまとめて数字にする関数は、1桁ずつ10^nを掛けて合計する手法も取りうる。サンプルコードとして見かけるのは下記のような書き方なのだが、今回実行時間を計測してみたらstring->numberを使った書き方の方が圧倒的に早かった。標準関数最強ですね。

; 別Ver
(define (intlist->integer2 lis)
  (let loop ((digit (- (length lis) 1))
             (lis lis)
             (ans 0))
    (if (null? lis) ans
      (loop (- digit 1) (cdr lis) (+ ans (* (car lis) (expt 10 digit)))))))

; もう一つ別Ver
(define (intlist->integer3 lis)
  (let loop ((lis lis)
             (ans 0))
    (if (null? lis) ans
      (loop (cdr lis) (+ (* ans 10) (car lis))))))
gosh> (define *long-intlist*
  (fold append '() (map (^(x) (iota x 1)) (iota 100 9 0))))

gosh> (time-this 1000 (lambda () (intlist->integer *long-intlist*)))
gosh> #<time-result 1000 times/  0.256 real/  0.297 user/  0.000 sys>
gosh> (time-this 1000 (lambda () (intlist->integer2 *long-intlist*)))
gosh> #<time-result 1000 times/ 11.007 real/ 13.594 user/  0.766 sys>
gosh> (time-this 1000 (lambda () (intlist->integer3 *long-intlist*)))
gosh> #<time-result 1000 times/  1.370 real/  2.578 user/  0.140 sys>

さらに余談だが、Schemeの世界には(正確数であれば)桁あふれというものはない。メモリの続く限り永遠に大きな数を扱うことができる。例えば上記のlong-intlistからは900桁の整数が作れるが、何も考えずに実行して問題ない。

gosh> (intlist->integer *long-intlist*)
gosh> 123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789123456789

他言語で桁あふれを気にしなければならないシーンで、ノーケアでコードを書くことに集中できるのは大きなアドバンテージではないだろうか(ただしメモリ制限を超えない限り……)。

10.多次元配列のソート

(define *sample-list* '((3 3 3) (1 4 4) (2 2 5) (1 2 3)))

; 普通に昇順ソート
(sort *sample-list*)
; 普通に昇順ソートしてreverse取って降順にする
(reverse (sort *sample-list*))
; 比較器は変えずに第1要素だけでソートする
(sort *sample-list* default-comparator car)
; 比較関数を与えることで降順にもできる
(sort *sample-list* > car)
; 第2要素だけで降順ソートする
(sort *sample-list* > cadr)
; 第1要素を削って第2要素以降だけでソートする
(sort *sample-list* default-comparator cdr)
; 第3要素と第1要素をこの優先度でソートする
(sort *sample-list* default-comparator (^ (l) (list (caddr l) (car l))))

第1キーは昇順、第2キーは降順、というように昇順降順を切り替えるスマートな書き方は思いつかなかった。sortを2回適用するので良いだろう。

11.defaultdict(辞書)

Gaucheではdictionary型としてはhash-tableが使える
Pythonのdefaultdictの利点はkeyが無いのにvalueを参照しようとしたときにエラーが出ない、ということらしい。Scheme/Gaucheではインクリメント+代入をあまり使わないので、相当する処理を考えると以下のようなことになるだろうか。

; 文字列中に出現する文字をhash-tableを使ってカウントする
(define *ht* (make-hash-table))

(let loop [(lis (string->list "mygkgoqqllqiqeewywrtuljoujehmzebhojbzvltzmhgyqddjpdkkhudysygdeydsqpnbnomfgwwtiwlqneadrbgqmwrpcqptxli"))]
  (cond [(null? lis) #t]
        [else (hash-table-put! *ht*
                               (car lis)
                               (+ (hash-table-get *ht* (car lis) 0) 1))
              (loop (cdr lis))]))

hash-table-get関数はkeyに対応するvalueが無かった場合デフォルトではエラーを返すが、:optional defaultを与えればエラーの代わりにdefaultを返すことになっている。今回は初期値として使いたいのでdefaultに0を与えている。

ただ、単なる数え上げなら既存の値の参照などせずともhash-table-push!関数でとりあえず詰め込んでしまって、あとからlengthで要素数をカウントするのが簡単ではないかという気もする。

(let loop [(lis (string->list "mygkgoqqllqiqeewywrtuljoujehmzebhojbzvltzmhgyqddjpdkkhudysygdeydsqpnbnomfgwwtiwlqneadrbgqmwrpcqptxli"))]
  (cond [(null? lis) #t]
        [else (hash-table-push! *ht* (car lis) (car lis))
              (loop (cdr lis))]))
  => #t
(map (^(x) (list x (length (hash-table-get *ht* x))))
     (hash-table-keys *ht*))
  => ((#\p 4) (#\y 6) (#\g 6) (#\b 4) (#\k 3) (#\t 4) (#\x 1) (#\f 1) (#\o 4) (#\s 2) (#\j 4) (#\a 1) (#\n 3) (#\e 6) (#\w 6) (#\i 3) (#\r 3) (#\d 7) (#\v 1) (#\m 5) (#\z 3) (#\q 9) (#\h 4) (#\c 1) (#\l 6) (#\u 3))

この辺は好みもあるだろう。

12.カウンター

要素の数を数える方法はSchemeにもある。

(define *sample-list* '(#\a #\a #\a #\a #\b #\c #\c))

(length *sample-list*)
(count (^(x) (eqv? x #\a)) *sample-list*)

PythonのCounterクラスというのはデータ構造自体はhash-tableのようなのでこのように作ってやる。

; keyが文字、valueが文字の出現個数のhash-tableを返す
(define (make-counter-hash-table lis)
  (let [(ht (make-hash-table 'eqv?))]
    (map (^(x) (hash-table-put! ht x (count (^(y) (eqv? y x)) lis)))
         (delete-duplicates lis))
    ht))

実行例

; 実行例
gosh> (define *counter-ht*
  (make-counter-hash-table *sample-list*))
gosh> *counter-ht*
gosh> (hash-table-get *counter-ht* #\a 0)
gosh> 4
gosh> (hash-table-get *counter-ht* #\x 0)
gosh> 0

hash-tableなので以下のような手続きが使える。

; keys
gosh> (hash-table-keys *counter-ht*)
gosh> (#\c #\b #\a)
; values
gosh> (hash-table-values *counter-ht*)
gosh> (2 1 4)
; alistを生成
gosh> (hash-table->alist *counter-ht*)
gosh> ((#\a . 4) (#\b . 1) (#\c . 2))
; listのlistを生成
gosh> (hash-table-map *counter-ht* list)
gosh> ((#\a 4) (#\b 1) (#\c 2))

Scheme/Gaucheにはmost_commonメソッドなどというものはないが、データがリスト化できればsort関数で降順ソートしてやればよい。

; alistでソート
(sort (hash-table->alist *counter-ht*) > cdr)
; listのlistでソート
(sort (hash-table-map *counter-ht* list) > cadr)

あとはcaarするなりlastするなりmap carするなりすればよい。

条件を満たす要素の個数をカウントするには、filterとlengthを組み合わせればよい。

(length (filter (^(x) (< x 0)) (iota 11 -5)))

ジェネレータ式に相当するものがわからなかったが、foldで条件判定と総和を一度に処理してやればよいのではないだろうか。

(fold (^(x y) (if (< x 0) (+ 1 y) (+ 0 y))) 0 (iota 11 -5))

条件を反転させたい場合……notくらいSchemeにもあるヨ。なおnotは関数なのでカッコを外すことはできない(Lisperには当たり前でも他言語話者にとっては当たり前ではないらしい)。

(length (filter (^(x) (not (< x 0))) (iota 11 -5)))

文字列の最後がeで終わるかどうか、を判定するのはこう書ける。

(eqv? #\e (last (string->list "hogehoge")))

あとはlambda式にしてfilterに渡すなりなんなり好きにすればよい。

Scheme/Gaucheでは文字列の置換はstring-replaceではなくて、regexp-replace-allを使う。例えば文字列から「,.」を除去したい場合はこうなる。

(regexp-replace-all #/\.|,/ "government of the people, by the people, for the people." "")

string-splitは普通。

gosh> (string-split "government of the people by the people for the people" " ")
gosh> ("government" "of" "the" "people" "by" "the" "people" "for" "the" "people")

なお、make-counter-hash-table関数は等価性の判定にeqv?を使っているので文字列の評価も可能で、変更無しでそのまま文字列の出現回数をカウントするために使える。

13.ModInt、14.逆元

元記事で紹介されていた意図だと、計算時間の短縮と桁あふれを防ぐためらしい。Scehemeでは後者は気にする必要が無さそうだが前者は気にしたほうが良いことは間違いない。Modulo計算をじっくりやったのは素数夜曲を読んだ時以来か。Project Eulerでもあまりお目にかからない。

Pythonでの実装はModIntクラスを実装して四則演算をオーバーライドしている。こういうときオブジェクト指向な言語はラクだなと思う。Schemeも関数のオーバーライドくらい朝飯前だが、Schemeぢからが足りないせいかクラスで隠蔽するのと同じようなことを実現する方法が思いつかない。ある環境下で四則演算は全てModulo計算に置き換えてしまったほうがいいなら、letで書き換えてしまってもいいかもしれない。

実際に実装する際は元記事のさらに大元のこの記事が役に立った。

加減乗算までは特に文句ないだろう。

(define (m+ m x y)
  (let [(x (if (< m x) (modulo x m) x))
        (y (if (< m y) (modulo y m) y))]
    (modulo (+ x y) m)))

(define (m- m x y)
  (let [(x (if (< m x) (modulo x m) x))
        (y (if (< m y) (modulo y m) y))]
    (modulo (- x y) m)))

(define (m* m x y)
  (let [(x (if (< m x) (modulo x m) x))
        (y (if (< m y) (modulo y m) y))]
    (modulo (* x y) m)))

除算の前に逆元を取る計算だが、大体下記のようになるだろう。modが素数に限定される場合のアルゴリズムとして勉強がてらmodinv-prime関数も実装してみた。標準関数にexpt-modがあるのでありがたく使う。modinvは拡張Euclidの互除法を使うことでmodが素数でない場合でも計算できるようにしてある。

(define (modinv-prime m x)
  (expt-mod x (- m 2) m))

(define (ext-gcd a b)
  (cond [(= b 0) (values a 1 0)]
        [else (let-values [((d x y) (ext-gcd b (remainder a b)))]
                (values d y (- x (* (quotient a b) y))))]))

(define (modinv m x)
  (let-values [((d x y) (ext-gcd x m))]
    (modulo x m)))

計算結果を確かめてみると、確かにmodが素数でない場合は結果が違っているし、法modとaが互いに素でない場合もおかしな結果になっている。

gosh> (iota 12 1)
gosh> (1 2 3 4 5 6 7 8 9 10 11 12)
gosh> (map (^(x) (modinv-prime 12 x)) (iota 12 1))
gosh> (1 4 9 4 1 0 1 4 9 4 1 0)
gosh> (map (^(x) (modinv 12 x)) (iota 12 1))
gosh> (1 1 1 1 5 1 7 11 11 11 11 0)
gosh> (map (^(x) (m* 12 x (modinv-prime 12 x))) (iota 12 1))
gosh> (1 8 3 4 5 0 7 8 9 4 11 0)
gosh> (map (^(x) (m* 12 x (modinv 12 x))) (iota 12 1))
gosh> (1 2 3 4 1 6 1 4 3 2 1 0)

というわけなので、うっかり使ってしまわないように逆元が取れない場合は#fを返すようにしてみた。

; 修正版
(define (modinv m x)
  (cond [(= 0 (modulo x m)) #f]
        [(= 1 (gcd m x)) (let-values [((d x y) (ext-gcd x m))]
                           (modulo x m))]
        [else #f]))
gosh> (map (^(x) (modinv 12 x)) (iota 12 1))
gosh> (1 #f #f #f 5 #f 7 #f #f #f 11 #f)
gosh> (map (^(x) (if (modinv 12 x) (m* 12 x (modinv 12 x)) #f)) (iota 12 1))
gosh> (1 #f #f #f 1 #f 1 #f #f #f 1 #f)

この時除算はこう書ける。

(define (m/ m x y)
  (let [(x (if (< m x) (modulo x m) x))
        (inv (modinv m y))]
    (if inv 
      (modulo (* x inv) m)
      #f)))
gosh> (m/ 12 11 7)
gosh> 5
gosh> (m/ 12 11 2)
gosh> #f

累乗はそれこそexpt-modで良いし、二項係数はメモ化してやれば良いように思う。

冒頭で述べた「let環境下で関数をオーバーライドする」方法は下記。

(define (modint-functions mod)
  (define (m+ m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (y (if (< m y) (modulo y m) y))]
      (modulo (+ x y) m)))

  (define (m- m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (y (if (< m y) (modulo y m) y))]
      (modulo (- x y) m)))

  (define (m* m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (y (if (< m y) (modulo y m) y))]
      (modulo (* x y) m)))

  (define (ext-gcd a b)
    (cond [(= b 0) (values a 1 0)]
          [else (let-values [((d x y) (ext-gcd b (remainder a b)))]
                  (values d y (- x (* (quotient a b) y))))]))

  (define (modinv m x)
    (cond [(= 0 (modulo x m)) #f]
          [(= 1 (gcd m x)) (let-values [((d x y) (ext-gcd x m))]
                             (modulo x m))]
          [else #f]))

  (define (m/ m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (inv (modinv m y))]
      (if inv 
        (modulo (* x inv) m)
        #f)))

  (values (^ (x y) (m+ mod x y))
          (^ (x y) (m- mod x y))
          (^ (x y) (m* mod x y))
          (^ (x y) (m/ mod x y))
          (^ (x y) (expt-mod x y mod))))


; 加減乗除関数には2引数しか取れなくなることに注意
(let-values [((+ - * / expt) (modint-functions 1000000007))]
  (* 111111111 (* 123456789 987654321)))
;  => 769682799 

通常の+-*/のようにn引数を取るように改良することもできるはずなんだけど、Schemeぢからが足りてないので宿題にします。

15.階乗

それぞれ逐次計算方式、メモ化再帰、Modulo取った乗算でのメモ化再帰Ver。

; 毎回素直に再帰計算
(define (simple-factorial n)
  (let loop [(n n)]
    (cond [(<= n 1) 1]
          [else (* n (loop (- n 1)))])))

; 計算結果を捨てずにメモする。永続的に使いまわしたいのでクロージャにしている。
(define (make-memoized-factorial)
  (let [(cache (hash-table-r7 'eqv? 0 1))]
    (^ (x)
       (let loop [(n x)]
         (let [(ans (hash-table-get cache n #f))]
           (cond [ans ans]
                 [(<= n 0) 1]
                 [else (let [(ans (* n (loop (- n 1))))]
                         (hash-table-put! cache n ans)
                         ans)]))))))

; 乗算をmodを法としたmoduloを取る算法で計算する。こちらもメモ化してクロージャ化。
(define (make-memoized-factorial-mod mod)
  (define (m* m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (y (if (< m y) (modulo y m) y))]
      (modulo (* x y) m)))

  (let [(cache (hash-table-r7 'eqv? 0 1))
        (m* (^(x y) (m* mod x y))) ]
    (^ (x)
       (let loop [(n x)]
         (let [(ans (hash-table-get cache n #f))]
           (cond [ans ans]
                 [(<= n 0) 1]
                 [else (let [(ans (m* n (loop (- n 1))))]
                         (hash-table-put! cache n ans)
                         ans)]))))))

クロージャは使うときにdefineする。

実行結果。メモ化してないと2000程度でもうヤバい。メモ化すると160ms。さらにModulo計算使うと12ms。メモ化してるから繰り返し評価しても仕方ないのでこれで終わり。

gosh> (define memoized-factorial (make-memoized-factorial))
gosh> (define memoized-factorial-mod (make-memoized-factorial-mod (+ 7 (expt 10 9))))
gosh> (use gauche.time)
gosh> (time-this 1 (^() (map simple-factorial (iota 2001 1))))
gosh> #<time-result 1 times/  2.843 real/  4.750 user/  0.235 sys>
gosh> (time-this 1 (^() (map memoized-factorial (iota 10001 1))))
gosh> #<time-result 1 times/  0.160 real/  0.188 user/  0.047 sys>
gosh> (time-this 1 (^() (map memoized-factorial-mod (iota 10001 1))))
gosh> #<time-result 1 times/  0.012 real/  0.015 user/  0.000 sys>

16.コンビネーション

二項係数nCrはProject Eulerでも頻出のテーマで、問題によって、つまりどの項がいつどういったタイミングで必要になるのか、によって取りうる選択肢が変わってくるように思う。

元記事では階乗と階乗の逆元をメモ化しておき利用する戦略が取られているが、どうせメモ化するのならnCrの値自体をメモ化してしまって、nCrとnC(r-1)の漸化式に持ち込むか、もしくはパスカルの三角形を使ったほうがいいんじゃないだろうか。計算負荷を考えたらメモ化しているとはいえ加算の方が乗算より速いだろう。

(define (make-memoized-pascal-triangle)
  (let [(cache (hash-table-r7 'equal? (list 0 0) 1 (list 1 0) 1 (list 0 1) 1))]
    (^ (n r)
       (let loop [(n n) (r r)]
         (let [(ans (hash-table-get cache (list n r) #f))]
           (cond [ans ans]
                 [(or (< n 0) (< r 0)) #f]
                 [(= r 0) 1]
                 [(= n r) 1]
                 [else (let [(ans (+ (loop (- n 1) (- r 1)) (loop (- n 1) r)))]
                         (hash-table-put! cache (list n r) ans)
                         ans)]))))))

階乗と同様にModulo版。

(define (make-memoized-pascal-triangle-mod mod)
  (define (m+ m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (y (if (< m y) (modulo y m) y))]
      (modulo (+ x y) m)))

  (let [(cache (hash-table-r7 'equal? (list 0 0) 1 (list 1 0) 1 (list 0 1) 1))
        (m+ (^(x y) (m+ mod x y)))]
    (^ (n r)
       (let loop [(n n) (r r)]
         (let [(ans (hash-table-get cache (list n r) #f))]
           (cond [ans ans]
                 [(or (< n 0) (< r 0)) #f]
                 [(= r 0) 1]
                 [(= n r) 1]
                 [else (let [(ans (m+ (loop (- n 1) (- r 1)) (loop (- n 1) r)))]
                         (hash-table-put! cache (list n r) ans)
                         ans)]))))))

ところがn=100万程度でメモリを食いつぶし帰ってこなくなってしまった。

(define pascal-triangle-mod (make-memoized-pascal-triangle-mod (+ 7 (expt 10 9))))
(pascal-triangle-mod 1000000 123456)
;  => ...メモリ食いすぎ&時間かかりすぎ!

メモリ使用量はO(n^2)に対して、階乗&階乗逆元版の方はO(2n)なのでその点で有利なのかな。というわけでそちらも書いてみる。

; 階乗&階乗の逆元をメモ化する改良版
(define (make-memoized-combination-mod mod)
  ; Mod計算の関数を読み込み
  (define (m* m x y)
    (let [(x (if (< m x) (modulo x m) x))
          (y (if (< m y) (modulo y m) y))]
      (modulo (* x y) m)))

  (define (ext-gcd a b)
    (cond [(= b 0) (values a 1 0)]
          [else (let-values [((d x y) (ext-gcd b (remainder a b)))]
                  (values d y (- x (* (quotient a b) y))))]))

  (define (modinv m x)
    (cond [(= 0 (modulo x m)) #f]
          [(= 1 (gcd m x)) (let-values [((d x y) (ext-gcd x m))]
                             (modulo x m))]
          [else #f]))

  ; メモ化階乗関数
  (define (make-memoized-factorial-mod mod)
    (let [(cache (hash-table-r7 'eqv? 0 1))
          (m* (^(x y) (m* mod x y)))]
      (^ (x)
         (let loop [(n x)]
           (let [(ans (hash-table-get cache n #f))]
             (cond [ans ans]
                   [(<= n 0) 1]
                   [else (let [(ans (m* n (loop (- n 1))))]
                           (hash-table-put! cache n ans)
                           ans)]))))))

  ; メモ化階乗逆元関数。メモを共用したいので、上の階乗関数クロージャをfact-procに渡す。
  (define (make-memoized-factinv-mod fact-proc mod)
    (let [(cache (make-hash-table))]
      (^ (x)
         (let [(ans (hash-table-get cache x #f))]
           (cond [ans ans]
                 [else (let [(ans (modinv mod (fact-proc x)))]
                         (hash-table-put! cache x ans)
                         ans)])))))

  ; 二項係数の計算式本体
  (let* [(m* (^(x y) (m* mod x y)))
         (fact (make-memoized-factorial-mod mod))
         (fact-inv (make-memoized-factinv-mod fact mod))]
    (^ (n r)
       (cond [(or (< n 0) (< r 0)) #f]
             [(= r 0) 1]
             [(= n r) 1]
             [else
               (m* (fact n) (m* (fact-inv r) (fact-inv (- n r))))]))))

Mod計算の関数も詰め込んでいるのでいろいろ長くなった。階乗の逆元もメモ化したほうがいいのかちょっと疑問(逆元くらい階乗のメモから逐次計算してもバチは当たるまい)だが、まあいいや。

gosh> (use gauche.time)
gosh> (define combination-mod (make-memoized-combination-mod (+ 7 (expt 10 9))))
gosh> (time-this 1 (^() (combination-mod 1000000 123456)))
gosh> #<time-result 1 times/  1.238 real/  2.890 user/  0.063 sys>

十分早くなったと思う。同じnとrなら2回目からはほぼコスト0で呼び出せる。
注意点としては、modとして与えた数よりも大きなnを扱う場合は逆元が計算できずエラーが出る可能性がある。この辺は問題の制約条件と相談になるかな。

; 階乗&階乗の逆元をメモ化する改良版*Modulo無し*
(define (make-memoized-combination)
  ; メモ化階乗関数
  (define (make-memoized-factorial)
    (let [(cache (hash-table-r7 'eqv? 0 1))]
      (^ (x)
         (let loop [(n x)]
           (let [(ans (hash-table-get cache n #f))]
             (cond [ans ans]
                   [(<= n 0) 1]
                   [else (let [(ans (* n (loop (- n 1))))]
                           (hash-table-put! cache n ans)
                           ans)]))))))

  ; メモ化階乗逆元関数。メモを共用したいので、上の階乗関数クロージャをfact-procに渡す。
  (define (make-memoized-factinv fact-proc)
    (let [(cache (make-hash-table))]
      (^ (x)
         (let [(ans (hash-table-get cache x #f))]
           (cond [ans ans]
                 [else (let [(ans (/ 1 (fact-proc x)))]
                         (hash-table-put! cache x ans)
                         ans)])))))

  ; 二項係数の計算式本体
  (let* [(fact (make-memoized-factorial))
         (fact-inv (make-memoized-factinv fact))]
    (^ (n r)
       (cond [(or (< n 0) (< r 0)) #f]
             [(= r 0) 1]
             [(= n r) 1]
             [else
               (* (fact n) (fact-inv r) (fact-inv (- n r)))]))))

Modulo計算無し版も作ってみたが間違い探しレベル。ただ扱う数が莫大な桁数になってしまうのでメモが破裂する。n=10000くらいなら計算可能。なるほど、だからModulo計算が必要なのね……。

17.ある値の累乗の何度も使用する場合

expt-modをメモ化したらいいと思う。割愛(実際の問題を解いたらまた更新するかも)。

18.最大公約数

Scheme/Gaucheの場合gcdもlcmもn引数に対応している。

gosh> (gcd (* 2 3 5) (* 13 3 5) (* 107 3 5))
gosh> 15
gosh> (lcm (* 2 3 5) (* 13 3 5) (* 107 3 5))
gosh> 41730

19.中国剰余定理

いろいろ見て回ったが、結局けんちょんさんの記事に流れ着きこれを参考にした。そしてこのあたりからだんだんPythonのコードよりC++のコードを見る時間の方が増えてきた……。

(define (ext-gcd a b)
  (cond [(= b 0) (values a 1 0)]
        [else (let-values [((d x y) (ext-gcd b (remainder a b)))]
                (values d y (- x (* (quotient a b) y))))]))

; 2元の場合
(define (chinese-remainder-theorem b1 m1 b2 m2)
  (let-values [((d p q) (ext-gcd m1 m2))]
    (cond [(not (= 0 (modulo (- b2 b1) d))) (value 0 -1)]
          [else (let* [(m (lcm m1 m2))
                       (tmp (modulo (* (quotient (- b2 b1) d) p) (quotient m2 d)))
                       (r (modulo (+ b1 (* m1 tmp)) m))]
                  (values r m))])))

(chinese-remainder-theorem 2 3 3 5)
;  => 8
;     15

; n元の場合
(define (chinese-remainder-theorem blist mlist)
  (let loop [(r 0) (M 1) (b blist) (m mlist)]
    (cond [(null? b) (values (modulo r M) M)]
          [else (let-values [((d p q) (ext-gcd M (car m)))]
                  (cond [(not (= 0 (modulo (- (car b) r) d))) (value 0 -1)]
                        [else (let* [(tmp (modulo (* (quotient (- (car b) r) d) p) (quotient (car m) d)))]
                                (loop (+ r (* M tmp)) (* M (quotient (car m) d)) (cdr b) (cdr m)))]))])))

(chinese-remainder-theorem '(2 3 2) '(3 5 7))
;  => 23
;     105

余談だが、どこだったかに「中国剰余定理は巨大な人数を数えるための方法だった」と書かれていて、始皇帝が円周率を計算する話を思い出した。中国ェ……。

20.DFS

とりまc++のコードを元のエッセンスを残しつつ書き直すとこうなる。

(let [(N 10) (M 2)]
  (let loop [(a '()) (i 0)]
    (cond [(<= N (length a)) (print (reverse a))]
          [(<= (- M 1) i) (loop (cons i a) 0)]
          [else (loop (cons i a) 0)
                (loop a (+ i 1))])))
;  => (0 0 0 0 0 0 0 0 0 0)
;     (0 0 0 0 0 0 0 0 0 1)
;     (0 0 0 0 0 0 0 0 1 0)
;     (0 0 0 0 0 0 0 0 1 1)
;     (0 0 0 0 0 0 0 1 0 0)
;     ...
;     (1 1 1 1 1 1 1 0 1 1)
;     (1 1 1 1 1 1 1 1 0 0)
;     (1 1 1 1 1 1 1 1 0 1)
;     (1 1 1 1 1 1 1 1 1 0)
;     (1 1 1 1 1 1 1 1 1 1)

元のc++コードでは再帰関数の中でさらにfor文を使っていて幾分気持ち悪さを感じるが、これなら1関数でループしている風情になってスッキリ気持ちが良い(偏見)。

「6.BIT全探索」と似てるけどちょっと違う。

21.二分探索(bisect)

GaucheにはPythonのbisectに対応するものは無いので、tree-mapを使って実装する。tree-mapはhash-table同様keyとvalueを格納できるが、keyの順序が保たれる。ので、ABC-143-D問題ではまず棒を長い順にソートしvectorに入れ、棒の長さをkey、vectorのindexをvalueにしたtree-mapを作成する。

ACコードは下記。Python版のように300msでは終わらず、大体1500msほどかかっている。ギリギリ。

(define (answer o)
  (format #t "~A~%" o))

(answer
  (let* ((N (read))) (read-line)
    (let [(bar-list (let loop [(n N) (ans '())]
                      (cond ((= n 0) (reverse (sort ans)))
                            (else (loop (- n 1) (cons (read) ans))))))]
      (let [(v (list->vector bar-list))
            (tm (make-tree-map))]
        (vector-for-each-with-index (^(i x) (tree-map-put! tm x i)) v)
        (let loop [(a 0) (b 1) (ans 0)]
          (cond
            [(>= a (- N 2)) ans]
            [(>= b (- N 1)) (loop (+ a 1) (+ a 2) ans)]
            [else (let* [(c (tree-map-successor-value tm (- (vector-ref v a) (vector-ref v b))))
                         (count (if c (max (- c b) 0) 0))]
                    (loop a (+ b 1) (+ count ans)))]))))))

vector-for-each-with-indexが使える。

余談だが、3重ループで処理する方法だとTLEしてしまった

さらに余談だが、GaucheではN個の数字列が1行で与えられたとき、文字列として読み込んでsplitしてstring->numberするのと、1個ずつreadして数字を読み込むのとでは速度にそう大差はないようだ。実際にACが出るコードも試してみたがあまり変わらない

22.二分探索(その他)

その他というかめぐる式二分探索がテーマ。元記事のリンクをたどっていって、結局けんちょんさんの記事にたどり着いたのでこれを参考に実装した。

; isOK関数の代わり
(define (greater-equal? v index key)
  (>= (vector-ref v index) key))

(define (vector-binary-search v is-ok? key)
  (define (medium a b) (quotient (+ a b) 2))
  (let loop [(ng -1) (ok (vector-length v)) (mid (medium -1 (vector-length v)))]
    (cond [(<= (abs (- ok ng)) 1) ok]
          [(is-ok? v mid key) (loop ng mid (medium ng mid))]
          [else (loop mid ok (medium mid ok))])))

Schemeにもwhileはあるが、whileを使うと冗長になるというか見通しが悪くなるのでnamed ledでループを書いている。あとisOK関数は引数として渡せるようにしてみた。

使うときはこうする。

(define int-vector (list->vector '(1 14 32 51 51 51 243 419 750 910)))
(vector-binary-search int-vector greater-equal?  50)
;  => 3

一般化したのでこうすることも可能。

(define (greater? v index key)
  (> (vector-ref v index) key))
(vector-binary-search int-vector greater? 50)
;  => 3
(vector-binary-search int-vector greater? 51)
;  => 6
(vector-binary-search int-vector greater? 52)
;  => 6

ただ、元記事にも書いてある通り、vectorは昇順ソートされていることが前提だし、昇順ソートされたvectorを渡した場合は<(小なり)のときに#tを返す評価関数は受け付けない。これも元記事にあったように、昇順か降順かに合わせてngとokの初期値を逆転させてやったり、使う評価関数に気を使ってやらないといけない。いっかなめぐる関数といえども、いかなる状況でも使える、というわけにはいかないようだ。

かといってGaucheのtree-mapのtree-map-floor関数などが直観的にわかりやすいかというとそうでもない。なんとなくfloorとceilingが逆ではという気もする。慣れの問題なのかもしれないが……。

関数 意味 解釈というかイメージ
tree-map-floor x <= probe を満たす最大のx probeから見て足元(今いるところ=probeを含む)
tree-map-ceiling x >= probe を満たす最小のx probeから見て頭がつっかえた天井(頭がつっかえているのでprobeを含む……だいぶ苦しいか)
tree-map-predecessor x < probe を満たす最大のx probeから見て前の人
tree-map-successor x > probe を満たす最小のx probeから見て後ろの人

23.最長増加部分列(LIS)

本家元記事が参照していたこちらの記事から例を取り、
https://qiita.com/python_walker/items/d1e2be789f6e7a0851e5
実装はこちらのGitHubを参考にした。
https://tjkendev.github.io/procon-library/python/dp/lis.html

ベクタに対して二分探索を取りつつ、更新していく。

(define *inf* (expt 10 10))

(define *dp-vector* (make-vector 10 *inf*))
(vector-set! *dp-vector* 0 -1)
; (vector-ref *dp-vector* 0)
;  => -1

; keyを超えない場合#t
(define (less-than? v index key)
  (< (vector-ref v index) key))

; ベクタを対象に二分探索をする
(define (vector-binary-search v is-ok? key)
  (define (medium a b) (quotient (+ a b) 2))
  ; vが昇順ソートされていて、かつ小さい側を答えにしたい場合はokを-1、ngをlength+1に設定する
  (let loop [(ng (vector-length v)) (ok -1) (mid (medium -1 (vector-length v)))]
    (cond [(<= (abs (- ok ng)) 1) ok]
          [(is-ok? v mid key) (loop ng mid (medium ng mid))]
          [else (loop mid ok (medium mid ok))])))

; 回答本体
(let loop [(l '(1 3 5 2 4 6))]
  (cond [(null? l) (vector-binary-search *dp-vector* less-than? *inf*)]
        [else (let [(idx (vector-binary-search *dp-vector* less-than? (car l)))]
                (vector-set! *dp-vector*
                             (+ 1 idx)  
                             (min (car l) (vector-ref *dp-vector* (+ 1 idx)))))
              (loop (cdr l))]))
;  => 4
*dp-vector*
;  => #(-1 1 2 4 6 10000000000 10000000000 10000000000 10000000000 10000000000)

実際に解いた問題

……ところで、c++よく知らないんですけどGitHubの実装例の

dp[idx] = min(a, dp[idx])

って

dp[idx++] = min(a, dp[idx++])

の間違いなんじゃ……???

24.累積和(合わせて、いもす法)

え!! ワンライナーで累積和を!?

出来らあっ!

(fold-right (^ (x l) (cons x (map (^(a) (+ a x)) l))) '() (iota 10 1))
;  => (1 3 6 10 15 21 28 36 45 55)
(cdr (reverse (fold (^(x l) (cons (+ x (car l)) l)) '(0) (iota 10 1))))
;  => (1 3 6 10 15 21 28 36 45 55)

冗談はさておき1つ目はmapを使っているためn^2/2の計算量だし、2つ目はreverseしてknilに与えている0を削るなど無駄が多い。……と思っていたら、けんちょんさんの記事によれば、累積和の先頭S_0は0であったほうが良いようだ。

(reverse (fold (^(x l) (cons (+ x (car l)) l)) '(0) (iota 10 1)))
;  => (0 1 3 6 10 15 21 28 36 45 55)

素直に書くなら次のようになるだろう。

(define (accumulate l)
  (let loop [(l l) (tmp 0)]
    (cond [(null? l) '()]
          [else (cons (+ (car l) tmp) (loop (cdr l) (+ (car l) tmp)))])))

(accumulate (iota 10 1))
;  => (1 3 6 10 15 21 28 36 45 55)

このあたりからはけんちょんさんの記事を参考にして実装した。
https://qiita.com/drken/items/56a6b68edef8fc605821

累積和の計算方法にどれを使うかは好みだが、例えば3から6の部分数列の和はこう書ける。

(define *cum* (list->vector (reverse (fold (^(x l) (cons (+ x (car l)) l)) '(0) (iota 10 1)))))

(define (sub-sequence-sum head tail)
  (- (vector-ref *cum* tail) (vector-ref *cum* head)))

(sub-sequence-sum 3 6)
;  => 15

ABC-084-D問題、もしくはエラトステネスの篩再び

これも上述のけんちょんさの記事にある実装例を参考にして、エラトステネスの篩を作る部分はこう書ける。

(define (sieve-list n)
  (let [(limit (sqrt n))]
    (let loop [(p-list (iota (- n 1) 2))]
      (cond [(null? p-list) '()]
            [(< limit (car p-list)) p-list]
            [(= 0 (car p-list)) (cons (car p-list) (loop (cdr p-list)))]
            [else (cons (car p-list)
                        (loop (filter-map (^(x) (if (= 0 (modulo x (car p-list))) 0 x))
                                          (cdr p-list))))]))))

(list->vector (cons 0 (cons 1 (sieve-list 100))))
;  => #(0 1 2 3 0 5 0 7 0 0 0 11 0 13 0 0 0 17 0 19 0 0 0 23 0 0 0 0 0 29 0 31 0 0 0 0 0 37 0 0 0 41 0 43 0 0 0 47 0 0 0 0 0 53 0 0 0 0 0 59 0 61 0 0 0 0 0 67 0 0 0 71 0 73 0 0 0 0 0 79 0 0 0 83 0 0 0 0 0 89 0 0 0 0 0 0 0 97 0 0 0)

なおGaucheのmath.primeモジュールには素数の無限遅延シーケンスprimesという便利なものがあって、これを利用するとこんな風にも書ける。

(use math.prime)
(define (sieve-list n)
  (let loop [(prime-list *primes*) (k 0)]
    (cond [(> k n) '()]
          ; [(= k 1) (cons 1 (loop prime-list (+ k 1)))]
          [(= k (car prime-list)) (cons (car prime-list) (loop (cdr prime-list) (+ k 1)))]
          [else (cons 0 (loop prime-list (+ k 1)))])))

(list->vector (sieve-list 100))
;  => #(0 0 2 3 0 5 0 7 0 0 0 11 0 13 0 0 0 17 0 19 0 0 0 23 0 0 0 0 0 29 0 31 0 0 0 0 0 37 0 0 0 41 0 43 0 0 0 47 0 0 0 0 0 53 0 0 0 0 0 59 0 61 0 0 0 0 0 67 0 0 0 71 0 73 0 0 0 0 0 79 0 0 0 83 0 0 0 0 0 89 0 0 0 0 0 0 0 97 0 0)

※1が素数リストに載っているのはちょっと違和感が拭えないのでindex=1は0としている。
一つ性能比較もしておこう。

(use gauche.time)

; 素朴な実装による篩関数をdefineしたあとに実行
(time-this 1 (lambda () (sieve-list (expt 10 5))))
;  => #<time-result 1 times/  1.404 real/  1.625 user/  0.031 sys>

; *prime*を使った実装による篩関数をdefineしたあとに実行
(reset-primes)
(time-this 1 (lambda () (sieve-list (expt 10 5))))
;  => #<time-result 1 times/  0.016 real/  0.016 user/  0.000 sys>

prime版の方だけなんだかめちゃくちゃ速いぞ?
(※実はprimeを使った関数の方は、モジュール読み込みや無限遅延シーケンスのリセット関数などの実行時間も含めて計測してみたのだがそれでも早い)
これはあれかな、素朴な実装の方は毎回listにmapを適用しているから遅いのかな? というわけで破壊的代入ありのvector版。

(define (sieve-vector n)
  (let [(limit (sqrt n))
        (ans (list->vector (iota (+ n 1) 0)))]
    (let loop [(k 2)]
      (cond [(< limit k) ans]
            [(= 0 (vector-ref ans k)) (loop (+ k 1))]
            [else (let loop2 [(a 2)]
                    (cond [(< n (* a k)) #t]
                          [else (vector-set! ans (* a k) 0)
                                (loop2 (+ a 1))]))
                  (loop (+ k 1))]))))

(sieve-vector 100)
;  => #(0 1 2 3 0 5 0 7 0 0 0 11 0 13 0 0 0 17 0 19 0 0 0 23 0 0 0 0 0 29 0 31 0 0 0 0 0 37 0 0 0 41 0 43 0 0 0 47 0 0 0 0 0 53 0 0 0 0 0 59 0 61 0 0 0 0 0 67 0 0 0 71 0 73 0 0 0 0 0 79 0 0 0 83 0 0 0 0 0 89 0 0 0 0 0 0 0 97 0 0 0)
(time-this 1 (lambda () (sieve-vector (expt 10 5))))
;  => #<time-result 1 times/  0.016 real/  0.015 user/  0.000 sys>

おお、速くなった。
あとはこの篩を使って「2017に似た数」かどうかを判定し、累積和をやっていけばよい。

いもす法

実際に解いた問題

indexの扱いが面倒くさい。問題に最適化すると個々の関数が汎用性を失ってしまうので、使うときに気を付けなければならない。あと普通に1桁数字リストを数に変換すると先頭の0が消えてしまうので、単純に1桁数字を文字に変換し結合するだけの関数もこさえてみた。

; 1桁数字リストを文字列に連結する
(define (intlist->string lis)
  (list->string (map (^(x) (integer->char (+ x 48))) lis)))

しゃくとり法

しゃくとり法については同じくけんちょんさんのこちらの記事を参考にした。
https://qiita.com/drken/items/ecd1a472d3a0e7db8dce

まず最もシンプルなAOJ Course The Number of Windowsの問題は下記のようになる。AOJはSchemeでの提出を受け付けていないので実装のみ。

(let [(n 12)
      (x 25)
      (a '(4 6 7 8 1 2 110 2 4 12 3 9))]
  (let loop [(left 0) (right 0) (v (list->vector a)) (sum 0) (res 0)]
    (cond [(>= left n) res]
          [(and (< right n) (<= (+ sum (vector-ref v right)) x))
           (loop left (+ right 1) v (+ sum (vector-ref v right)) res)]
          [(= left right) (loop (+ left 1) (+ right 1) v sum (+ res (- right left)))]
          [else (loop (+ left 1) right v (- sum (vector-ref v left)) (+ res (- right left)))])))

loopが1個しかないのに手続き型の書き方では2重、3重ループになっているところがLispの面白いところだナ……。

25.短めな順列の作成

もちろん順列・組み合わせ、累乗集合もある。

(use util.combinations)

; 順列
(permutations '(1 1 2))
;  => ((1 1 2) (1 2 1) (1 1 2) (1 2 1) (2 1 1) (2 1 1))
; 重複する要素は除く
(permutations* '(1 1 2))
;  => ((1 1 2) (1 2 1) (2 1 1))

; 残念だが文字列を渡してもエラーになる
(permutations "Qiita")
;  => ERROR
; リストにしてから渡す
(permutations* (string->list "Qiita"))
;  => ((#\Q #\i #\i #\t #\a) (#\Q #\i #\i #\a #\t) ... (#\a #\t #\i #\i #\Q))

; 組み合わせ
(combinations '(1 1 2 3) 2)
;  => ((1 1) (1 2) (1 3) (1 2) (1 3) (2 3))
; 重複する要素は除く
(combinations* '(1 1 2 3) 2)
;  => ((1 1) (1 2) (1 3) (2 3))

; 累乗集合(全てのサブセット)
(power-set '(1 1 2 3))
;  => (() (1) (1) (2) (3) (1 1) (1 2) (1 3) (1 2) (1 3) (2 3) (1 1 2) (1 1 3) (1 2 3) (1 2 3) (1 1 2 3))
; 重複する要素は除く
(power-set* '(1 1 2 3))
;  => (() (1) (2) (3) (1 1) (1 2) (1 3) (2 3) (1 1 2) (1 1 3) (1 2 3) (1 1 2 3))

26.二次元累積和

今回も結局けんちょんさんの記事にお世話になった。
https://qiita.com/drken/items/56a6b68edef8fc605821#4-%E4%BA%8C%E6%AC%A1%E5%85%83%E7%B4%AF%E7%A9%8D%E5%92%8C

(use gauche.array)

; 2次元累積和配列を作成する
(define (make-accumulate-array h w arry)
  (let [(s (make-array (shape 0 (+ h 1) 0 (+ w 1)) 0))]
    (dotimes (row h)
      (dotimes (col w)
        (array-set! s (+ row 1) (+ col 1)
                    (-  (+ (array-ref s row (+ col 1))
                           (array-ref s (+ row 1) col)
                           (array-ref arry row col))
                        (array-ref s row col)))))
    s))

; 2次元累積和配列から目当ての範囲の総和を得る
; 範囲指定は[x1 x2)[y1 y2)の半開区間
(define (two-dimension-accumulate arry x1 y1 x2 y2)
  (-  (+ (array-ref arry x2 y2)
         (array-ref arry x1 y1))
      (array-ref arry x1 y2) 
      (array-ref arry x2 y1)))

(let* [(H 4)
         (W 5)]
    ; (read-line)
    (let [(a (make-array (shape 0 H 0 W) 0))]
      (let loop [(row 0) (col 0) (lis '((1 8 7 3 2)
                                        (9 1 3 4 6)
                                        (3 5 8 1 4)
                                        (2 7 3 2 5)))]
        (cond [(null? lis) #t]
              [(null? (car lis)) (loop (+ row 1) 0 (cdr lis))]
              [else (array-set! a row col (caar lis))
                    (loop row (+ col 1) (cons (cdar lis) (cdr lis)))])) 
      ; 二次元累積和
      (let [(s (make-accumulate-array H W a))]
        ; 結果
        (list
          (two-dimension-accumulate s 1 2 3 5)
          (two-dimension-accumulate s 0 1 2 3)
          (two-dimension-accumulate s 0 0 4 5)))))
;  => (26 19 84) 

範囲指定の引数を「x1 x2 y1 y2」から「x1 y1 x2 y2」にしているけどこの辺は好み。実際の問題を解いたときに不都合だったら変えるかも。

それにしても配列を使うためにモジュールをuse(他言語でいうところのimport)しないといけないのはSchemeらしい。

27.ダイクストラ法

本家も含めていろいろなサイトを参照したが、いまいちわかりにくい。ここまで例題を解きながら実装することで理解してきたのだが、AtCoderにダイクストラ法の単純な問題が見つからず苦労した。いろいろな人の実装例を見ても、ちょっとずつ違うので余計混乱する。結局Wikipediaのダイクストラ法のページを最も参考にした。

; ダイクストラ法の優先度付きキュー(ヒープ)を利用したコード。
; 例題として下記のものを参考にしている
; https://ja.wikipedia.org/wiki/%E3%83%80%E3%82%A4%E3%82%AF%E3%82%B9%E3%83%88%E3%83%A9%E6%B3%95
; エッジの表現やヒープに詰める情報は(node番号 . コスト)としている。
(use data.heap)

(let* [(INF (expt 10 10))
       (N 6)
       (M 9)]
  (let* [(edges (make-hash-table))
         (costs (make-vector N INF))
         (nodes-in-searching (make-binary-heap :key cdr))
         (prev (make-vector N -1))]
    ; 例題は無向グラフなのでエッジは両方向登録
    (map (^(lis) (begin
                   (hash-table-push! edges (- (car lis) 1) (cons (- (cadr lis) 1) (caddr lis)))
                   (hash-table-push! edges (- (cadr lis) 1) (cons (- (car lis) 1) (caddr lis)))))
         ; 次のエッジ情報はWikipediaの例題の通り。
         '((1 2 7) (1 3 9) (1 6 14) (2 3 19) (2 4 15) (3 4 11) (3 6 2) (4 5 6) (5 6 9)))
    ; 初期化
    (vector-set! costs 0 0)
    (vector-map-with-index (^(i x) (binary-heap-push! nodes-in-searching (cons i x)))
                           costs)
    ; 本計算
    (let loop []
      (cond [(binary-heap-empty? nodes-in-searching)
             (values costs prev)]
            [else
              (let [(current-node (binary-heap-pop-min! nodes-in-searching))]
                    (map (^(next-node)
                           (let [(to (car next-node))
                                 (to-cost (+ (cdr current-node) (cdr next-node)))]
                             (if (< to-cost (vector-ref costs to))
                               (begin
                                 (vector-set! costs to to-cost)
                                 (vector-set! prev to (car current-node))
                                 (binary-heap-push! nodes-in-searching (cons to to-cost))))))
                         (hash-table-ref edges (car current-node)))
                    (loop))]))))
;  => #(0 7 9 20 20 11)
;     #(-1 0 0 2 5 2)

ところでヒープがGahceのデフォルトライブラリに存在したのでありがたかった。SLIBからひっぱってこなければならないかなと思っていた……。

ゴールが最初から決まっている場合は、ゴールノードの最短経路が求まったら計算を打ち止めしたり、他にも枝刈りの手法がありそう。特にヒープの初期化の際、始点ノードのコストだけ入れても処理はできるんだけど、上記コードでは愚直に実装してある。これは始点からたどり着けないノード対策かなと思ったけど、例えそうだったとしても最小コストの更新が起こりえないので無意味だ。

気になる場合は初期化のコード

    (vector-map-with-index (^(i x) (binary-heap-push! nodes-in-searching (cons i x)))
                           costs)

を次のように書き換えておけば良い(ただしindex=0が始点ノードの場合)。

(binary-heap-push! nodes-in-searching (cons 0 0))

28.ワーシャルフロイド法

本家と同じくこちらを参考にした。
https://qiita.com/okaryo/items/8e6cd73f8a676b7a5d75

(use gauche.array)

(let* [(INF (expt 10 10))
       (edges (make-array (shape 0 4 0 4) INF))]
  ; 初期化
  (map (^(lis) (begin
                 (array-set! edges (car lis) (cadr lis) (caddr lis))
                 (array-set! edges (cadr lis) (car lis) (caddr lis))))
       '((0 1 5) (0 2 3) (0 3 8) (1 3 2) (2 3 20)))
  ; 始点と終点が同じエッジは0埋めする
  (array-for-each-index edges (^(i j) (if (= i j) (array-set! edges i j 0))))
  ; 計算本体
  (let loop [(k 0) (i 0) (j 0)]
    (cond [(>= k 4) edges]
          [(>= i 4) (loop (+ k 1) 0 0)]
          [(>= j 4) (loop k (+ i 1) 0)]
          [else 
            (array-set! edges i j (min (array-ref edges i j)
                                       (+ (array-ref edges i k) (array-ref edges k j))))
            (loop k i (+ j 1))])))
;  => #,(<array> (0 4 0 4) 0 5 3 7 5 0 8 2 3 8 0 10 7 2 10 0)

総当たりしてるのに結果が正しいのか直観的に理解できないアルゴリズムって難しいなと思う。

29.セグメント木(セグ木)

元記事も参照しているスライドが役に立った。
https://www.slideshare.net/iwiwi/ss-3578491

; セグメント木のクロージャ
; nは2の累乗
(define (make-segtree n)
  (let* [(INF (expt 10 10))
         (dat (make-vector (- (* 2 n) 1) INF))]
    ; ツリーの幅
    (define (size) n)
    ; 初期化
    (define (initialize segtree lis)
      (vector-map-with-index (^(i x) (vector-set! segtree (+ i (- n 1)) x)) (list->vector lis)) 
      (let loop [(i 0)]
        (cond [(and (<= (- n 1) i) (< i (* 2 n))) (vector-ref segtree i)]
              [else (vector-set! segtree i (min (loop (+ (* i 2) 1)) (loop (+ (* i 2) 2))))
                    (vector-ref segtree i)])))
    ; 検索
    (define (reference segtree i)
      (let [(i (+ i (- n 1)))]
        (vector-ref segtree i)))
    ; 実データ部分を取得
    (define (reference-all segtree)
      (vector-copy segtree (- n 1) (- (* 2 n) 2)))
    ; 更新
    (define (update segtree i x)
      (let [(i (+ i (- n 1)))]
        (vector-set! segtree i x)
        (let loop [(i i)]
          (cond [(<= i 0) #t]
                [else (let [(i (quotient (- i 1) 2))]
                        (vector-set! segtree i
                                     (min (vector-ref segtree (+ (* i 2) 1))
                                          (vector-ref segtree (+ (* i 2) 2))))
                        (loop i))]))))
    ; 区間最小を探し出す
    (define (query segtree a b k l r)
      (cond [(or (<= r a) (<= b l)) INF]
            [(and (<= a l) (<= r b)) (vector-ref segtree k)]
            [else (min (query segtree a b (+ (* 2 k) 1) l (quotient (+ l r) 2))
                       (query segtree a b (+ (* 2 k) 2) (quotient (+ l r) 2) r))]))
    (values
      (^(lis) (initialize dat lis))
      (^(i) (reference dat i))
      (^() (reference-all dat))
      (^(i x) (update dat i x))
      (^(a b) (query dat a b 0 0 n)))))

updateとqueryがあれば事足りるようだが、確認のためにreferenceなど用意した。
実際に使うとしたらこんな感じ。これもqueryに渡す区間は半開区間。

(let-values [((segtree-init segtree-ref segtree-ref-all segtree-update! segtree-query) (make-segtree (expt 2 3)))]
  (segtree-init '(5 3 7 9 1 4 6 2))
  (values 
    (segtree-ref-all)
    (segtree-query 0 1)
    (segtree-query 0 2)
    (segtree-query 2 4)
    (segtree-query 7 8)
    (segtree-query 0 8)
    (segtree-query 8 9)))
;  => #(5 3 7 9 1 4 6)
;     5
;     3
;     7
;     2
;     1
;     10000000000

とりあえずの実装なので後々書き直すかも。

30.転倒数

さあついに大台だー! 転倒数の前にBIT。
実装にはこのページが役に立った。
https://algo-logic.info/binary-indexed-tree/

; BIT
(define (make-bit n)
  (let* [(dat (make-vector (+ n 1) 0))]
    ; 1-indexedとしてi番目の要素にxを加算する 
    (define (bit-add bit i x)
      (let loop [(i i)]
        (cond [(> i n) #t]
              [else (vector-set! bit i (+ (vector-ref bit i) x))
                    (loop (+ i (logand i (- i))))])))
    ; 1-indexedとして1からi番目の要素までの合計を計算する
    (define (bit-sum bit i)
      (let loop [(i i)]
        (cond [(<= i 0) 0]
              [else (+ (vector-ref bit i) 
                       (loop (- i (logand i (- i)))))])))
    (values
      (^(i x) (bit-add dat i x))
      (^(i) (bit-sum dat i)))))

(let-values [((bit-add bit-sum) (make-bit 8))]
  (bit-add 1 1)
  (bit-add 3 2)
  (bit-add 5 1)
  (bit-sum 7))
;  => 4

転倒数の仕組みについてはこのスライドが良さそう。
https://www.slideshare.net/hcpc_hokudai/binary-indexed-tree

31.桁DP

ひとまずけんちょんさんの記事を参考に、ABC-117-D問題Schemeで解いてみた

1桁ごとに状態管理をしつつ走査していくイメージか。

その後効率的な書き方がわからなくて結局元記事が参考にしている記事を同じく参考にした。長くなりすぎたので別記事に分けた。

32.ダブリング

これも長くなったので別記事にまとめた。

33.最長共通部分列(Common Longest Subsequence)

34.bitDP

ひとまず例題を1問AC通したが記事としては書き途中。

35.10進数を2進数、8進数、16進数に変換

とても簡単。

(define i 255)

(number->string i 2)
(number->string i 8)
(number->string i 16)

実行結果。

gosh> i
gosh> "11111111"
gosh> "377"
gosh> "ff"

逆関数ももちろんある。

(string->number "11111111" 2)
(string->number "377" 8)
(string->number "ff" 16)

実行結果。

gosh> 255
gosh> 255
gosh> 255

数xを2進数表記した場合に現れる1の数をpopulation countと言うらしい。調べたらSRFI-60にlogcountという関数があった。ついでにinteger-lengthについても紹介しておく。

(logcount (string->number "10100111" 2))
;  => 5
(integer-length (string->number "10100111" 2))
;  => 8

36.10進数とn進数の変換

nを基数がbaseの数とみなして、10進数に変換する関数。

(define (change-base-to-10 n base)
  (let loop [(lis (reverse (map (^x (string->number (list->string (list x))))
                                (string->list (number->string n)))))
             (dig-base 1)]
    (cond [(null? lis) 0]
          [else (+ (* (car lis) dig-base)
                   (loop (cdr lis) (* dig-base base)))])))

例題として挙がっているABC-192-D-Base nを解いてみた。

最初はそもそも二分探索を適用しないとTLEしてしまうことや、エッジケースであるM=(expt 10 18)の時でWAが出てしまうなど単純なミスでつまづいてしまっていた。修行が足りない……。

37.正規表現

2
2
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
2
2