25 September, 2014

[IVA] Chapter 2, Section 9, Exercise 13

Exercise 13 について自作コードで試してみた。(あれからIntではなくBigIntを使う様にしてみたり、計算順序の最適化を実施したりなど。sugerや斉次化はまだ未実装。) が、色々やばいですね、lexでの計算。grevlexなら簡単に終わるのに。

scala> import basic._
import basic._

  val f13_1_grevlex = Poly(Seq(
     (Mono(Seq(5,0,0))(grevlex), Q(1,1)),
     (Mono(Seq(0,4,0))(grevlex), Q(1,1)),
     (Mono(Seq(0,0,3))(grevlex), Q(1,1)),
     (Mono(Seq(0,0,0))(grevlex), Q(-1,1))
     ))(grevlex)                                  //> f13_1_grevlex  : basic.Poly = x^5 + y^4 + z^3 + -1
  val f13_2_grevlex = Poly(Seq(
     (Mono(Seq(3,0,0))(grevlex), Q(1,1)),
     (Mono(Seq(0,2,0))(grevlex), Q(1,1)),
     (Mono(Seq(0,0,2))(grevlex), Q(1,1)),
     (Mono(Seq(0,0,0))(grevlex), Q(-1,1))
     ))(grevlex)                                  //> f13_2_grevlex  : basic.Poly = x^3 + z^2 + y^2 + -1
  val buch3_base13_grevlex = Poly.reduce(Poly.buchberger(Seq(f13_1_grevlex, f13_2_grevlex))(grevlex))(grevlex)
                                                  //> buch3_base13_grevlex  : Seq[basic.Poly] = List(
                                                  //| x^3 + z^2 + y^2 + -1,  
                                                  //| y^4 + -x^2z^2 + -x^2y^2 + z^3 + x^2 + -1)

scala> val lex = LexOrder(List('x,'y,'z))
lex: basic.LexOrder = LexOrder(List('x, 'y, 'z))

scala> val f13_1_lex = Poly(Seq(
     |      (Mono(Seq(5,0,0))(lex), Q(1,1)),
     |      (Mono(Seq(0,4,0))(lex), Q(1,1)),
     |      (Mono(Seq(0,0,3))(lex), Q(1,1)),
     |      (Mono(Seq(0,0,0))(lex), Q(-1,1))
     |      ))(lex) 
f13_1_lex: basic.Poly = x^5 + y^4 + z^3 + -1

scala> val f13_2_lex = Poly(Seq(
     |      (Mono(Seq(3,0,0))(lex), Q(1,1)),
     |      (Mono(Seq(0,2,0))(lex), Q(1,1)),
     |      (Mono(Seq(0,0,2))(lex), Q(1,1)),
     |      (Mono(Seq(0,0,0))(lex), Q(-1,1))
     |      ))(lex)
f13_2_lex: basic.Poly = x^3 + y^2 + z^2 + -1

scala> val buch3_base13_lex = Poly.reduce(Poly.buchberger(Seq(f13_1_lex, f13_2_lex))(lex))(lex)
buch3_base13_lex: Seq[basic.Poly] = List(x^3 + y^2 + z^2 + -1, x^2y^2 + x^2z^2 + -x^2 + -y^4 + -z^3 + 1, xy^4 + xz^3 + -x + y^4 + 2y^2z^2 + -2y^2 + z^4 + -2z^2 + 1, xy^2z + -xy^2 + 273/1408xz^12 + 1011/1408xz^11 + -227/1408xz^10 + -9665/4224xz^9 + 211/704xz^8 + 85/32xz^7 + -749/528xz^6 + 1/2xz^4 + -1/2xz^3 + -xz + x + -81/1408y^14 + -459/1408y^12z^2 + -9/176y^12z + 321/704y^12 + -819/1408y^10z^4 + -273/704y^10z^3 + 1323/704y^10z^2 + 39/176y^10z + -4373/4224y^10 + 273/1408y^8z^6 + -501/352y^8z^5 + 2329/1408y^8z^4 + 8551/4224y^8z^3 + -5871/1408y^8z^2 + 13/88y^8z + 395/352y^8 + -819/704y^6z^7 + 8373/1408y^6z^6 + 3213/704y^6z^5 + -3341/176y^6z^4 + -7555/2112y^6z^3 + 1219/66y^6z^2 + 53/88y^6z + -23341/4224y^6 + 273/704y^4z^9 + 9189/1408y^4z^8 + 773/352y^4z^7 + -33823/1056y^4z^6 + -17/3y^4z^5…
scala> buch3_base13_lex.length
res0: Int = 7

scala> buch3_base13_lex(0)
res2: basic.Poly = x^3 + y^2 + z^2 + -1

scala> buch3_base13_lex(1)
res3: basic.Poly = x^2y^2 + x^2z^2 + -x^2 + -y^4 + -z^3 + 1

scala> buch3_base13_lex(2)
res4: basic.Poly = xy^4 + xz^3 + -x + y^4 + 2y^2z^2 + -2y^2 + z^4 + -2z^2 + 1

scala> buch3_base13_lex(4)
res5: basic.Poly = xz^11 + 4xz^10 + xz^9 + -10xz^8 + -4xz^7 + 8xz^6 + 13145211549469524099276049970825174332001934982726555617/49464457656842324243991404386798887698363025825966149632y^22 + 336469772395565067430728152218633828000517403338752050447/98928915313684648487982808773597775396726051651932299264y^20z^2 + 41246821025682957247777391063405260354876124127586033329/98928915313684648487982808773597775396726051651932299264y^20z + -186688514037649654164494156802783046785181037030518263195/49464457656842324243991404386798887698363025825966149632y^20 + 3748115348128909106405657431258725581362590226678884005481/197857830627369296975965617547195550793452103303864598528y^18z^4 + 532878767777869557858761862556643001207912197601994175895/98928915313684648487982808773597775396726051651932299264...

scala> buch3_base13_lex(5)
res6: basic.Poly = y^12 + -y^10 + 3y^8z^3 + -5y^8z^2 + 2y^8 + -10y^6z^4 + 20y^6z^2 + -10y^6 + -7y^4z^6 + 30y^4z^4 + -6y^4z^3 + -30y^4z^2 + 13y^4 + -5y^2z^8 + 20y^2z^6 + -30y^2z^4 + 20y^2z^2 + -5y^2 + -z^10 + z^9 + 5z^8 + -13z^6 + 10z^4 + 3z^3 + -5z^2

scala> buch3_base13_lex(6)
res7: basic.Poly = x^2z^4 + x^2z^3 + -2x^2z^2 + -1/6xz^10 + -2/3xz^9 + -1/4xz^8 + 7/6xz^7 + -1/12xz^6 + -xz^5 + xz^4 + 1/2y^10z^2 + 1/3y^10z + -5/12y^10 + -1/6y^8z^4 + 1/3y^8z^3 + 1/4y^8z^2 + 1/2y^8 + y^6z^5 + -5/3y^6z^4 + -17/6y^6z^3 + 17/6y^6z^2 + 2/3y^6z + -1/3y^6 + -1/3y^4z^7 + -7/2y^4z^6 + -3y^4z^5 + 131/12y^4z^4 + 23/4y^4z^3 + -13y^4z^2 + -8/3y^4z + 29/6y^4 + -17/6y^2z^8 + -4y^2z^7 + 119/12y^2z^6 + 17/2y^2z^5 + -16y^2z^4 + -14/3y^2z^3 + 12y^2z^2 + 5/3y^2z + -55/12y^2 + -7/6z^10 + -7/6z^9 + 55/12z^8 + 21/4z^7 + -31/4z^6 + -15/2z^5 + 83/12z^4 + 41/12z^3 + -31/12z^2

scala>

21 September, 2014

[IVA] Chapter 2, Section 9, Exercise 1, 6, 11

Exersise 1.

$S=(c_1, \cdots, c_s)$, $T=(d_1, \cdots, d_s)$ は $F=(f_1, \cdots, f_s)$ の syzygy ゆえ、
$$\sum_{i=1}^s c_i \mathrm{LT}(f_i) = \sum_{i=1}^s d_i \mathrm{LT}(f_i) = 0$$.

$$\sum_{i=1}^s (c_i + d_i)  \mathrm{LT}(f_i) = \sum_{i=1}^s c_i \mathrm{LT}(f_i) + \sum_{i=1}^s d_i \mathrm{LT}(f_i) = 0$$,
$$\sum_{i=1}^s (g c_i)  \mathrm{LT}(f_i) = g \sum_{i=1}^s c_i \mathrm{LT}(f_i) = 0$$,
より
$S+T=(c_1 + d_1, \cdots, c_s + d_s)$, $g \cdot S = (g c_1, \cdots, g c_s)$ もまた syzygy となる。

Exercise 6.

$S_j$ が $S(G)$ の multideg $\gamma_j$ の homogeneous syzygy であるとき、$S_j \cdot G$ の multideg $< \gamma_j$ を示せ。

$G = (g_1, \cdots, g_s)$ とする。

$S_j \cdot G = \sum_{i=1}^s c_i x^{\alpha(i)} g_i$ であるが、
$$\mathrm{multideg}(S_j \cdot G) $$
$$= \mathrm{multideg}(\sum_{i=1}^s c_i x^{\alpha(i)} g_i) \leq \max_{i=1}^s \mathrm{multideg}(x^{\alpha(i)} g_i) $$
$$= \max_{i=1}^s \mathrm{multideg}(x^{\alpha(i)} \mathrm{LT}(g_i))$$

ここで等号成立は $\sum_{i=1}^s x^{\alpha(i)} \mathrm{LT}(g_i) \neq 0$ の場合であるが、

$S_j \subset S(G)$ であるから $\sum_{i=1}^s c_i x^{\alpha(i)} \mathrm{LT}(g_i) = 0$

より等号は成立しない。よって $\mathrm{multideg}(S_j \cdot G) < \max_{i=1}^s \mathrm{multideg}(x^{\alpha(i)} g_i)$ であるが、
homogeneous syzygy の定義より

$S_j = (c_1 x^{\alpha(1)}, \cdots, c_s x^{\alpha(s)})$ where $c_i \in k$ かつ $c_i \neq 0$ ならば $\alpha(i) + \mathrm{multideg}(g_i) = \gamma_j$

であるから、
 $\mathrm{multideg}(S_j \cdot G) < \gamma_j$ となる。

Exersise 11.

元の Buchberger では除算を 17 回、syzygy を使用したものは 4 回で済んでいる。

用いたコードは https://github.com/tmiya/scala_Buchberger にあります。

   // buchberger
  val buch_f1 = Poly(Seq(
     (Mono(Seq(3,0,0))(grlex), Q(1,1)),
     (Mono(Seq(1,1,0))(grlex), Q(-2,1))
     ))(grlex)                                    //> buch_f1  : basic.Poly = x^3 + -2xy
  val buch_f2 = Poly(Seq(
     (Mono(Seq(2,1,0))(grlex), Q(1,1)),
     (Mono(Seq(0,2,0))(grlex), Q(-2,1)),
     (Mono(Seq(1,0,0))(grlex), Q(1,1))
     ))(grlex)                                    //> buch_f2  : basic.Poly = x^2y + -2y^2 + x
  val buch_base1 = Poly.buchberger(Seq(buch_f1, buch_f2))(grlex)
                                                  //> divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| buch_base1  : Seq[basic.Poly] = ArrayBuffer(x^3 + -2xy, x^2y + -2y^2 + x, -
                                                  //| x^2, -2xy, -2y^2 + x)
  val buch_base2 = Poly.reduce(buch_base1)(grlex) //> buch_base2  : Seq[basic.Poly] = List(x^2, xy, y^2 + -1/2x)
  
  val buch2_base1 = Poly.buchberger2(Seq(buch_f1, buch_f2))(grlex)
                                                  //> divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| divmod
                                                  //| buch2_base1  : Seq[basic.Poly] = ArrayBuffer(x^3 + -2xy, x^2y + -2y^2 + x, 
                                                  //| -x^2, -2y^2 + x, -2xy)

08 September, 2014

ProofSummit 2014 & JSSST Coq Tutorial

名古屋で開催された ProofSummit 2014 と、同じく名古屋大学で開催されている日本ソフトウェア科学会のチュートリアル:定理証明支援系Coq入門 に参加しました。

ProofSummit は各種の定理証明系(今年は Coq, Agda, Mizarの話がありました)とその応用について年に一回(大体8月末〜9月初に)開催しているイベントで、2010年に開催された Coq庵 から数えると今年で5回目です。「Coq を使ってみた」的な初心者の発表から、定理証明系を使って計算機科学の研究をしているような研究者による発表など幅広い内容の発表があります。今年の発表のかなりの部分に付いては、上記リンクから発表資料が参照出来ると思います。
ProofSummit 今年は定理証明イベントと一緒に開催する関数型言語イベントが中止になったり、東京で ScalaMatsuri が開催されていたりなどで、ちょっと参加人数が減りましたが、来年もまた開催したいし、来年こそはなんかネタを用意して参加したいです。

Coq チュートリアルは、産総研でCoqを使って研究されている Affeldt さんによる、Coq を使う意義/Coqについて/Coqの拡張であるssreflectのチュートリアル、でした。聞きに行った主眼はssreflectのチュートリアルだったのですが、前半の定理証明系やCoqについての話の部分も非常に面白かったです。発表資料は Affeldt さんのページ から取得出来ます。
Coq チュートリアルは、100人以上が参加の、非常に盛況なチュートリアルでした。

その前の NII Shonan Coq セミナーと合わせて、この夏は Coq が熱い夏でした。今年のイベントとしてはあと、12/3-5に「高信頼な理論と実装のための定理証明および定理証明器」というのが九州で開催されます。

[IVA] Chapter 2, Section 8, Exercise 1,7

Mac版のAsirで計算しました

Exersise 1

$ /Applications/cfep.app/OpenXM/bin/openxm asir
[2189] G=gr([-x^3+y, x^2*y-z], [x,y,z], 0);
[-z*x+y^2,y*x^2-z,-x^3+y]
[2190] p_nf(x*y^3-z^2+y^5-z^3, G, [x,y,z], 0);
0

Exersise 7

$ /Applications/cfep.app/OpenXM/bin/openxm asir
[1891] load("gr");
[1998] load("noro_pd.rr");
[2182] B=[u*t-x, 1-u-y, u+t-u*t-z];
[-x+t*u,-y-u+1,-z+(-t+1)*u+t]
[2183] V=[u,t,x,y,z];
[u,t,x,y,z]
[2184] G=gr(B,V,2);
[y*x+y^2+(z-2)*y-z+1,-x-y-z+t+1,-y-u+1]
[2185] noro_pd.elimination(G,[x,y,z]);
[y*x+y^2+(z-2)*y-z+1]

01 September, 2014

NII Shonan School on Coq 参加報告


Coq というのは定理証明支援系のソフトウェアです。関数型言語に興味のある方の中にはカリー=ハワード対応という言葉をご存知の方もいるのではないかと思いますが、これはプログラム⇔証明、型⇔命題という対応がある、というものです。Coq は依存型という、通常の静的型付け言語 (Haskell とか Scala など) より強力な型システムを用いることで、述語論理 (「全てのεに対して、あるδが存在して、どうのこうの」の様な量化子を含んだ式を表現出来る論理) の証明に対応するプログラムを書き、そしてその証明が正しいかを型検査することで機械的に確認出来ます。

じゃあ Coq を使うと何が良いのかというと、まず数学的な定理の証明を機械検証可能な形で証明出来ます。最近だとケプラー予想の機械的証明が話題になりましたね。他にも四色問題などが定理証明系の上で証明されています。どちらも数学者の人力では証明の検証が難しい問題です。
それとは別に、ソフトウェアが満たすべき仕様を述語論理で記述して、プログラムを実装し、実装が仕様を満たすことを証明出来ます。ユニットテストなどのテストでは網羅的でなかったりなどの問題がありますがCoqを使えばそのような問題はありません。
では、そんな素晴らしいツールがなぜ世の中で広く使われていないのかというと...単純に証明を書くのは難しいからです。Coqで証明を書くにはCoqの技術を磨く必要があります。
そしてその機会が、今回のセミナーという訳です。

セミナーは、葉山のセミナーハウスにて合宿形式5日間コース(半日の鎌倉見物含む)で開催され、内容はCoqの開発元のINRIAから訪れた講師によるCoqチュートリアルでした。
朝から夕食の時間までCoqに関する授業あるいは課題による実習、その後も解けなかった問題を寝るまで解く、というハードな毎日でした。私自身はここんとこあまりCoqで何かを証明してなかったので、多少錆び付いていたCoqのスキルを取り戻す良い機会でした。

残念ながら講義資料は参加者以外に公開されていない様です(パスワードで保護)。ただ、Coq については名古屋のProofCafeや、東京だとSF読み進捗ダメです会議 #9 #readcoqart #Coq(こちらは私も参加しています)といった勉強会が定期開催されているので、興味があれば参加されるとCoqについて学べると思います。

[IVA] Chapter 2, Section 7, Exercise 7,13

Exercise 1.

計算問題故パス

Exersise 7.

monomial order を固定し、ideal $I$ のminimal Groebner base $G$, $\tilde{G}$ について、

a. LT($G$) = LT($\tilde{G}$) を示せ。

$G=\{g_1, \cdots, g_t\}$, $\tilde{G} = \{\tilde{g}_1, \cdots, \tilde{g}_{t'}\}$ とする。monomial order に従って、$g_1 \geq g_2 \geq \cdots$, $\tilde{g}_1 \geq \tilde{g}_2 \geq \cdots$ と仮定して良い。

$G$, $\tilde{G}$ はグレブナ基底だから $\langle \mathrm{LT}(G) \rangle = \langle \mathrm{LT}(\tilde{G}) \rangle = \langle \mathrm{LT}(I) rangle$.

$\mathrm{multideg}(\mathrm{LT}(g_1)) = \alpha_1$,  $\mathrm{multideg}(\mathrm{LT}(\tilde{g}_1)) = \tilde{\alpha}_1$ とする。

$$LT(\tilde{g}_1) = \sum_i f_i \mathrm{LT}(g_i)$$
と書けるから、$\tilde{\alpha}_1 \geq \alpha_1$ 逆も同様に言えるから $\alpha_1 \geq \tilde{\alpha}_1$ よって $\tilde{\alpha}_1 = \alpha_1$ である。LCは全て$1$ゆえ$\tilde{\mathrm{LT}(g_1)} = \mathrm{LT}(g1)$. $\mathrm{LT}(G)-\mathrm{LT}(g_1)$, $\mathrm{LT}(\tilde{G})-\mathrm{LT}(\tilde{g})_1$ について同様に議論でき、$\mathrm{LT}(G)=\mathrm{LT}(\tilde{G})$.

b. $G$ と $\tilde{G}$ が同じ数の要素を含む事を示せ。

a. より LT($G$) と LT($\tilde{G}$) は同じ数の元を持つ事が判る。仮に $| \tilde{G} | \;\gt \; | G |$  とすると、$\tilde{G}$ には LT が等しい事なる元が存在することになるが、それは minimal Groebner basis であることに矛盾する。よって$| \tilde{G} | = | G |$ .

Exercise 13.

$F$ がグレブナ基底で無い場合は、$F=(f_1, \cdots,f_s)$ を順序付きタプルと考える必要がある。
$\bar{f}^F = 0$ ならば $f = \sum_i h_i f_i$ と書ける。このとき、$F' = F \cup \tilde{f}$ とすると、$f = \sum_i h_i f_i + 0 \cdot \tilde{f}$ ゆえ $\bar{f}^{F'} = 0$.