ゲーデルの不完全性定理の証明のアレをRacketで書いてみる (5)

「10.8.3 列」

定義6 n番目の要素

(define (CanDivideByPower x n k)
  (CanDivide x (expt (prime n x) k)))

(define (elm x n)
  (Min k ≦ x (and (CanDivideByPower x n k)
                  (not (CanDivideByPower x n (+ k 1))))))

(P n)ではなく(prime n x)を使っているので、歯抜けでも問題はないはず

> (elm 504 3)
1

おｋ

定義7 列の長さ

(define (len x)
  (Min k ≦ x (and (> (prime k x) 0)
                  (= (prime (+ k 1) x) 0))))

定義8 列の連結

⇒が出てきました

「ならば」とよみますけどifそのままってわけじゃないですね
A⇒Bは(if A B #t)か
まあこれでもいいけど毎回#tを書くのもアレだから定義しましょう

あ、これ関数で定義するとマズいやつ？
副作用はないけど余分な計算が走るのは避けたいので

(define-syntax-rule (⇒ x y)
  (or (not x) y))

これくらいならsyntax-parseでなくてもいいでしょう
さて定義

(define (M8 x y)
  (expt (P (+ (len x) (len y))) (+ x y)))

(define (** x y)
  (Min z ≦ (M8 x y)
       (and (∀ m ≦ (len x)
               (⇒ (<= 1 m) (= (elm z m) (elm x m))))
            (∀ n ≦ (len y)
               (⇒ (<= 1 n) (= (elm z (+ (len x) n)) (elm y n)))))))

テスト

> (** 8 4)
remainder: undefined for 0

ぐふ
CanDivideはただのreminderじゃなかった
0で割ってるときもエラーにせず#fを返すようにしてないと

(define (CanDivide x d)
  (and (not (zero? d))
       (= (remainder x d) 0)))

もう少し大きな数でテスト
[2 1]という列と[2 1]という列を連結してみます

> (time (** 8 4))   ; [3]という列と[2]という列の連結
cpu time: 3 real time: 3 gc time: 0
72
> (time (** 12 4))  ; [2 1]という列と[2]という列の連結
cpu time: 35571 real time: 35565 gc time: 1153
300
> (time (** 12 12)) ; [2 1]という列と[2 1]という列の連結

お返事がありません
うーん
ほとんど最小レベルの列の連結ですらここまでとは

どんなときもかならず1から順番に探すという作戦が無理あるんだよなー
せっかくマクロ書いたけどさー
どうしようかなー

ゲーデルの不完全性定理の証明のアレをRacketで書いてみる (4)

「10.8.2 整数論」の続き

（単位はミリ秒）

> (time (P 8))
cpu time: 12989 real time: 12992 gc time: 1441
19

これだとテストもままならなくなりそうなのでもうちょっと速くなるようにしてみますか

とりあえず簡単なところから

; もとの定義
;(define (CanDevide x d)
;  (∃ n ≦ x (= x (* d n))))

(define (CanDevide x d)
  (= (remainder x d) 0))

測定

> (time (P 8))
cpu time: 10410 real time: 10410 gc time: 1369
19

もうちょっと減るかと思ったんだけどなあ
こっちのほうが効くかなあ

; もとの定義
;(define (IsPrime x)
;  (and (> x 1)
;       (not (∃ d ≦ x (and (not (= d 1))
;                          (not (= d x))
;                          (CanDevide x d))))))

(define (IsPrime x)
  (and (> x 1)
       (let loop ((d 2))
         (cond ((> (* d d) x) #t)
               ((CanDevide x d) #f)
               (else (loop (+ d 1)))))))

どれ

> (time (P 8))
cpu time: 5008 real time: 5006 gc time: 1103
19

お、効いた

> (time (P 9))
cpu time: 92482 real time: 92457 gc time: 3484
23

前は300秒位かかってたからまあマシに
とはいえ分単位で待つわけには

やっぱここが無駄多すぎるんだな
ちゃんとプロファイリングしないとダメですよってことでしょうか

(define (M5 n)
  (+ (factorial n) 1))

(define (P n)
  (cond ((= n 0) 0)
        (else (Min p ≦ (M5 n) (and (< (P (- n 1)) p)
                                   (IsPrime p))))))

(P 8)を求めるとき正直にp=1から始めて
しかもそれぞれのpにつき(P 7)から(P 0)を計算し直してます
で、(P 7)を計算するときには正直に以下同文

ぱっと見(M5 n)まで計算するのかとまずびっくりしますが
実際のところpが(M5 n)まで達することはないので実は関係なし

いろいろやりかたはあると思いますがとりあえず見つけた素数は覚えておくことに

(define primes (make-hash))
; 0番目の素数は0という定義
(hash-set! primes 0 0)
(hash-set! primes 1 2)

(define (P n)
  (cond ((hash-ref primes n #f))
        (else (let loop ((k (+ (P (- n 1)) 1)))
                (cond ((IsPrime k)
                       (hash-set! primes n k)
                       k)
                      (else (loop (+ k 1))))))))

どれ

> (time (P 9))
cpu time: 0 real time: 0 gc time: 0
23
> (time (P 100))
cpu time: 1 real time: 0 gc time: 0
541
> (time (P 1000))
cpu time: 11 real time: 11 gc time: 0
7919
> (time (P 10000))
cpu time: 261 real time: 261 gc time: 0
104729

よしこれくらいなら先へ進んでももうしばらくは動かせるかな

このへんってただたまたま論理式の表現方法に素数指数表現を選んだから
こうなってるってだけで別に不完全性とは関係ないとこなんですよねえ
あ、モジュールに分けたほうがキレイかも？

コード的には中途半端な気もするけどこれくらいで先に進もう

ゲーデルの不完全性定理の証明のアレをRacketで書いてみる (3)

「10.8.2 整数論」

このへんは単なる数の計算

定義1 割り切れる？
関数名はできるだけ数学ガールに合わせていきます

(define (CanDevide x d)
  (∃ n ≦ x (= x (* d n))))

割り切れないときにx以下のnすべてについて試してしまってもったいないですが
とりあえずそういうことは言わないで進みます
言ってたらキリがないので！
超富豪です！

定義2 素数？

(define (IsPrime x)
  (and (> x 1)
       (not (∃ d ≦ x (and (not (= d 1))
                          (not (= d x))
                          (CanDevide x d))))))

言わないぞ言わないぞ・・・

定義3 n番目の素因数

(define (CanDevideByPrime x p)
  (and (CanDevide x p) (IsPrime p)))
(define (prime n x)
  (cond ((= n 0) 0)
        (else (Min p ≦ x (and (< (prime (- n 1) x) p)
                              (CanDevideByPrime x p))))))

例では2352=2^4×3^1×7^2が使われてます
5が抜けてますが歯抜けはアリってことですかね
ゲーデル数の定義ってそうなってるのかな

p.220

この自然数の列を一つの自然数にまとめたいとき、素数を小さい方から
並べた列〈2,3,5,7,11,13...〉を別途用意する。

「小さい方から並べた列」は「小さい方からすべて」なのか
「2,3,5,7」はただの一例なのか
微妙なところです

そういうコード、じゃないな関数定義になってるますしまあ気にしないで進みます

定義4 nの階乗

(define (factorial n)
  (cond ((= n 0) 1)
        (else (* n (factorial (- n 1))))))

ふつうです
再帰入門みたいなのにでてくるサンプルと同じ
ていうかプログラミングよりもこういう世界が先にあったんですよね

ifじゃなくてcondを使ってるのは、場合分けってなんとなくifより
condっぽいかなーという気分だけの話です

定義5 n番目の素数

(define (M5 n)
  (+ (factorial n) 1))

(define (P n)
  (cond ((= n 0) 0)
        (else (Min p ≦ (M5 n) (and (< (P (- n 1)) p)
                                   (IsPrime p))))))

すごい数の計算が行われてそうだ・・・

> (time (P 5))
cpu time: 7 real time: 7 gc time: 0
11
> (time (P 6))
cpu time: 42 real time: 42 gc time: 0
13
> (time (P 7))
cpu time: 683 real time: 682 gc time: 83
17
> (time (P 8))
cpu time: 12989 real time: 12992 gc time: 1441
19
> (time (P 9))
cpu time: 304833 real time: 305456 gc time: 37796
23

普通に処理時間が爆発してますね
(P 10)はちょっとやる気がしません
列が出てきたら苦しいだろうなとは思ってましたがその手前でもうこんな状態とは

「目で見る数学」

まんが数学史とかあったらいいんですけどね

娘が歴女になってしまいます - kb84tkhrのブログ

とまでは行きませんがブックオフで「目で見る数学」って本を見つけました

目で見る数学―美しい数・形の世界

 数字の変遷から0の発見あたりはちょっと歴史っぽいかも！
無限とかフラクタルとかカオスとかそのあたりまで出てきます
まずまず求めるものに近い感じ

ということで買って置いておいたら読んでました
数学者列伝のあたりもちょっと興味を引いてたみたいです

読む本、読まない本を見比べてみると
うーん見た目？
カラフルなことが大事っぽいんですよねえ特に表紙
中身は白黒の字ばっかりでも読むときは読むんですが
年齢的なところもあるのかな

続編もいっしょに置いてあったんですが数学というより物理に近い感じでした
ちょっと保留

いっしょに買ってきたこの本に食いついちゃったんでしまった失敗したかと
思いましたが読んでくれてよかった 

人生はあはれなり… 紫式部日記

子供向けじゃありませんがこういうのでもけっこう喜んで読んでます
微妙にダークな感じとかわかるのかな

ゲーデルの不完全性定理の証明のアレをRacketで書いてみる (2)

「10.8.1 装備を整える」の続き

同様に∃x≦M[...]とmin x≦M[...]を定義します

(define (∃≦ max f)
  (let loop ((x 1))
    (cond ((> x max) #f)
          ((f x) #t)
          (else (loop (+ x 1))))))

(define-syntax (∃ stx)
  (syntax-parse stx
    #:literals (≦)
    [(_ v:id ≦ max:expr body:expr)
     #'(∃≦ max (λ (v) body))]))
     
(define (Min≦ max f)
  (let loop ((x 1))
    (cond ((> x max) 0)
          ((f x) x)
          (else (loop (+ x 1))))))

(define-syntax (Min stx)
  (syntax-parse stx
    #:literals (≦)
    [(_ v:id ≦ max:expr body:expr)
     #'(Min≦ max (λ (v) body))]))

minだと名前がカブるのでMinにしました

Min≦やMinは(f x)を満たすxが見つからなかった場合は0を返します
x=0から開始すると見つからなかったのかx=0で見つかったのかわからないので
x=1から開始することにしました
∀や∃も同様に修正
カンが悪かった

ところで、∀・∃・Minはどう見てもそっくりですね
関数とマクロを作るマクロを作ればまとめられそうな気がします

(require (for-syntax syntax/parse))
(require (for-syntax racket/syntax))

(define ≦ #f)

(define-syntax (define-equipment stx)
  (syntax-parse stx
    ((_ name term notfound found)
     #:with fname (format-id stx "~a≦" #'name)
     #'(begin
         (define (fname max f)
           (let loop ((x 1))
             (cond ((> x max) (notfound x))
                   ((term (f x)) (found x))
                   (else (loop (+ x 1))))))
         (define-syntax (name stx)
           (syntax-parse stx
             #:literals (≦)
             [(_ v:id ≦ max:expr body:expr)
              #'(fname max (λ (v) body))]))))))

(define-equipment ∀ not (const #t) (const #f))
(define-equipment ∃ identity (const #f) (const #t))
(define-equipment Min identity (const 0) identity)

若干無理やり感ありますができたみたいです
∀と∃だけだったらもうちょっとシンプルになったと思うんですが、
Minもいっしょにしようとしたんでconstとか使うことに

けっこう調べた知識を総動員したかも
短くはなりましたがあとで見たらわかるかどうかちょっと心配

とりあえずこれくらいで先に進むことにします

定数の定義

; 定数

(define c0 1)
(define cf 3)
(define cnot 5)
(define cor 7)
(define call 9)
(define clp 11)
(define crp 13)

どんな名前にするのがいいのか悩みましたがまあこんなところで
（全角のゼロとかfとか使ってやろうかと思ったけどさすがにやめた）

それでは長旅に出かけましょう

ゲーデルの不完全性定理の証明のアレをRacketで書いてみる (1)

数学ガールのゲーデル巻を読んでると、「定義1」から「定義46」までなんだか
プログラムを読んでるみたいですよね？ですよね？

数学ガール/ゲーデルの不完全性定理 (数学ガールシリーズ 3)

しかも関数型っていうか関数そのものですね
なんか証明の中に定義が46個並んでるのは壮観な感じですが
プログラムの中に関数が46個っていうとなんてことない気もします

なーんてことを思ってたらじゃあちょっと書いてみるかという気に
書いたからどうなるってこともない気はするんですけど
ややこしくて理解が難しいのは多分その先の表現定理が出てきてからだし
どう考えたってマトモに動くはずもないし

でもまあやってみようかと
どこまで行けるかな
書くだけなら行けると思うけど

「10.8 《冬》証明可能性へ至る長い長い旅」

このあたりから
必要なら戻って

「10.8.1 装備を整える」

まずは∀x≦M[...]の定義から
こんな感じですかね
それっぽいunicode文字は使っていこうかと思ってます

(define (∀≦ max f)
  (let loop ((x 0))
    (cond ((> x max) #t)
          ((not (f x)) #f)
          (else (loop (+ x 1))))))

ちょっと動かしてみます
3行書いたら動かしてみるくらいじゃないと不安なんです

> (∀≦ 3 (λ (x) (< x 4)))
#t
> (∀≦ 3 (λ (x) (< x 3)))
#f

うん大丈夫

ところでxは0から始まるべきなのか1から始まるべきなのか
ちょっと自信がありませんがとりあえず0で進みます
辻褄が合わないところが出てきたら戻す

ここでsyntax-parse登場
ていうかここで使ってみようと思って先に勉強してたわけです

(require (for-syntax syntax/parse))

(define ≦ #f)

(define-syntax (∀ stx)
  (syntax-parse stx
    #:literals (≦)
    [(_ v:id ≦ max:expr body:expr)
     #'(∀≦ max (λ (v) body))]))

literalは何か定義しないとエラーになってしまいます
とりあえず定義されていれば値はなんでもいいようなので#fにして進みます

> (∀ x ≦ 3 (< x 4))
#t
> (∀ x ≦ 3 (< x 3))
#f

おー
ちょっとそれっぽくなりました
≦は完全に飾りですが
これからたくさん使うので、ラムダラムダ書かずに済むのはありがたいのでは

とりあえずここまでで満足！

RacketのMacroを調べてみる (11)

ではSyntax: Meta-Programming Helpersを見ていきます

1 Parsing and Specifying Syntax

• ここではsyntax-parseとdefine-syntax-classを説明していく
• syntax-parseはsyntax patternにしたがって自動的にエラーメッセージを生成する

1.1 Introduction

• syntax-parseとsyntax classを使って堅固なマクロの書き方を紹介する
• letを題材とする

最も単純には以下のように書けます

(define-syntax-rule (mylet ([var rhs] ...) body ...)
  ((lambda (var ...) body ...) rhs ...))

しかしちゃんと動くときはいいんですが

• 間違ったときはどこでどう間違ったのかわかりづらい
• 最悪、エラーにならずおかしな動きをする可能性がある
• マクロは構文を検証し、適切なメッセージを表示するべき
• 読みやすく維持しやすいコードで機械的にチェックされるのはマクロ作成者にも利点

このあと細かく段階を追ってひとつずつsyntax-parseの機能を紹介し、
それがどういうエラーをチェックしてどんなメッセージを出すか示しながら進みます

まずは基本

; syntax-parse関連をインポートする
(require (for-syntax syntax/parse))

(define-syntax (mylet stx)
  ; まずはsyntax-caseと同じような感じで書く
  (syntax-parse stx
    ; パターン変数をsyntax classでアノテートすることによって構文をチェックする
    ; デフォルトでは識別子を表すid、式を表すexpr（など？）が利用できる
    ; ここでは、varは識別子でrhsは式ですよ、ってこと
    ; `...+`は1個以上の繰り返し
    [(_ ([var:id rhs:expr] ...) body ...+)
     #'((lambda (var ...) body ...) rhs ...)]))

これでvarが識別子かどうかを調べてダメなら"mylet: expected identifier"みたいな
エラーを出してくれるようになります

次はsyntax classの定義です

(define-syntax (mylet stx)

  ; カスタムsyntax classの定義
  (define-syntax-class binding
    ; syntax-classの説明（エラーメッセージで使われる）
    #:description "binding pair"
    ; "binding"とは(識別子 式)という形をしたもの
    (pattern (var:id rhs:expr)))

  (syntax-parse stx
    ; カスタムsyntax classを使ってアノテーション
    [(_ (b:binding ...) body ...+)
     ; パターン変数.属性名という形で参照する 
     #'((lambda (b.var ...) body ...) b.rhs ...)]))

これで、"mylet: expected binding pair"のようなエラーメッセージが
出るようになります
これはなんだかいい感じがしてきました

さらに、letで割り当てる識別子が一意でないといけないという制約を入れるために
もうひとつsyntax classを定義します

(define-syntax (mylet stx)

  (define-syntax-class binding
    #:description "binding pair"
    (pattern (var:id rhs:expr)))

  (define-syntax-class distinct-bindings
    #:description "sequence of distinct binding pairs"
    ; bはbindingであり・・・
    (pattern (b:binding ...)
             ; さらに、一意でなければならない
             ; エラーメッセージも指定する
             #:fail-when (check-duplicate-identifier
                          (syntax->list #'(b.var ...)))
                          "duplicate variable name"
             ; 属性名をb.varと書く代わりにvarと書けるようにする
             #:with (var ...) #'(b.var ...)
             #:with (rhs ...) #'(b.rhs ...)))

  (syntax-parse stx
    ; bsはdistinct-bindings
    [(_ bs:distinct-bindings body ...+)
     ; 上で#:withしてないと、bs.b.varなどと書かなければいけない
     #'((lambda (bs.var ...) body ...) bs.rhs ...)]
    ; 名前付きletも使えるようにする
    [(_ loop:id bs:distinct-bindings body ...+)
     #'(letrec ([loop (lambda (bs.var ...) body ...)])
         (loop bs.rhs ...))]))

確かにsyntax-caseと比べて構文チェックはかなり強力になってる感じですね