バットとボールの問題を SageMath で解く

問題:バットとボールはセットで1ドル10セントします。バットはボールより1ドル高いです。ボールはいくらですか?

http://karapaia.com/archives/52047106.html

$ sage                                                             
┌────────────────────────────────────────────────────────────────────┐
│ SageMath version 7.2, Release Date: 2016-05-15                     │
│ Type "notebook()" for the browser-based notebook interface.        │
│ Type "help()" for help.                                            │
└────────────────────────────────────────────────────────────────────┘
sage: var('bat ball')
(bat, ball)
sage: f = [bat + ball == 110, bat == ball + 100]
sage: solve(f, bat, ball)
[[bat == 105, ball == 5]]

Gaucheで小町算にチャレンジ

Makoto HiroiさんのMemorandum(2013年1月5日)で小町算の問題が紹介されていました。

●パズルでプログラミング

パズルの世界では、1 から 9 までの数字を 1 個ずつすべて使った数字を「小町数」といいます。
たとえば、123456789 とか 321654987 のような数字です。「小町算」というものもあり、
たとえば 123 + 456 + 789 とか 321 * 654 + 987 のようなものです。

[問題] 小町算

1 から 9 までの数字を順番に並べ、間に + と - を補って三桁の値 (100 - 999) になる式を作ることにします。
100 になる式の一例を示します。

例:1 + 2 + 3 - 4 + 5 + 6 + 78 + 9 = 100

100 になる式は全部で 11 通りあります。それでは問題です。

1. 式の総数が最大になる値をすべて求めてください。
2. 解のない値で最小のものを求めてください。
3. 解のある値で最大のものを求めてください。 

Gauche(Scheme)でチャレンジしてみたいと思います。

はじめに、100になる式は全部で11通りあるとのことですので、それを確かめるプログラムを作成します。

(use util.match)

(define (komachi)
  (let loop ((expr '(1)) ; 計算式
             (rest '(2 3 4 5 6 7 8 9))) ; 残りの数字
    (cond ((null? rest)
           (when (and (integer? (car expr))
                      (= 100 (eval (reverse expr))))
             (print (reverse expr))))
          ((integer? (car expr))
           (loop (cons (+ (* 10 (car expr)) (car rest))
                       (cdr expr))
                 (cdr rest))
           (loop (cons '+ expr) rest)
           (loop (cons '- expr) rest))
          (else
           (loop (cons (car rest) expr) (cdr rest))))))

;;; 式を評価
(define (eval expr)
  (let loop ((expr (cdr expr))
             (acc  (car expr)))
    (match expr
      (() acc)
      (('+ x . z) (loop z (+ acc x)))
      (('- x . z) (loop z (- acc x))))))

実行結果です。

gosh> (komachi)
(123 + 45 - 67 + 8 - 9)
(123 + 4 - 5 + 67 - 89)
(123 - 45 - 67 + 89)
(123 - 4 - 5 - 6 - 7 + 8 - 9)
(12 + 3 + 4 + 5 - 6 - 7 + 89)
(12 + 3 - 4 + 5 + 67 + 8 + 9)
(12 - 3 - 4 + 5 - 6 + 7 + 89)
(1 + 23 - 4 + 56 + 7 + 8 + 9)
(1 + 23 - 4 + 5 + 6 + 78 - 9)
(1 + 2 + 34 - 5 + 67 - 8 + 9)
(1 + 2 + 3 - 4 + 5 + 6 + 78 + 9)
#<undef>

100になる式は全部で11通りありました。上のプログラムを元に以下の3つの問題を解いていきます。

  1. 式の総数が最大になる値をすべて求めてください。
  2. 解のない値で最小のものを求めてください。
  3. 解のある値で最大のものを求めてください。
(use util.match)
(use gauche.collection)
(use srfi-1) ; lset-difference

(define (komachi)
  (let ((ht (make-hash-table)))
    (let loop ((expr '(1)) ; 計算式
               (rest '(2 3 4 5 6 7 8 9))) ; 残りの数字
      (cond ((null? rest)
             (when (integer? (car expr))
               (let ((x (eval (reverse expr))))
                 (when (<= 100 x 999)
                   (hash-table-put! ht x (+ 1 (hash-table-get ht x 0)))))))
            ((integer? (car expr))
             (loop (cons (+ (* 10 (car expr)) (car rest))
                         (cdr expr))
                   (cdr rest))
             (loop (cons '+ expr) rest)
             (loop (cons '- expr) rest))
            (else
             (loop (cons (car rest) expr) (cdr rest)))))
    ht))

;;; 式を評価
(define (eval expr)
  (let loop ((expr (cdr expr))
             (acc  (car expr)))
    (match expr
      (() acc)
      (('+ x . z) (loop z (+ acc x)))
      (('- x . z) (loop z (- acc x))))))

(define (solve1) ; 問題1
  (let* ((alist (hash-table->alist (komachi)))
         (nmax (cdr (find-max alist :key cdr))))
    (map car (filter (^x (= nmax (cdr x))) alist))))

(define (solve2) ; 問題2
  (let ((alist (hash-table->alist (komachi))))
    (find-min (lset-difference = (iota 900 100) (map car alist)))))

(define (solve3) ; 問題3
  (let ((alist (hash-table->alist (komachi))))
    (car (find-max alist :key car))))

問題 1 の「全て求める」というのが結構くせ者ですね……。上のプログラムでは、まず最大値を見つけてから、再度その値を持つ要素を検索しています。1度きりの探索で答えを出せたら良いのですがうまく書けませんでした。

問題 2 は、存在しない値を求めるのにlset-differenceを使っています。最小値だけでなく、存在しない値を全て求めてしまっているのはあまり良くないかも知れません。

問題 3 は、単純にfind-maxで最も大きい値を求めています。

実行結果です。

gosh> (solve1)
(117 108 126)
gosh> (solve2)
160
gosh> (solve3)
972

参考

プログラミングGauche

プログラミングGauche

追記

  • 2013年1月12日のMemorandumに解答編が掲載されていました(2013-01-14)。

RubyのBinDataが便利!

これはすごい。。忘れないようにメモ。Rubyだとこんな書き方が出来るんですね。

深さ優先探索によるトポロジカルソート

トポロジカルソートのアルゴリズムは、『プログラマのうちあけ話―続・プログラム設計の着想』を読んで学びました。この本では次の方法により、トポロジカルソートします。

入力辺がないノードをキュー q に入れる
while (q に要素が含まれる):
  node = q.pop()
  print(node)  
  nodeの出力辺を全て削除する
  入力辺がないノードをキュー q に入れる(ただし、printしたnodeは除く)

トポロジカルソート - Wikipedia をみますと、上記アルゴリズムの他に、深さ優先探索による方法でトポロジカルソートを行うアルゴリズムが書かれてあります。

深さ優先探索によるトポロジカルソートのプログラムを Lua で書いてみました。

入力ファイルとして、Problem 79 - Project Eulerのkeylog.txtを使います。

tsort コマンド

tsort (Unix) - Wikipedia, the free encyclopedia というトポロジカルソートを行うコマンドがあるようです。この tsort コマンドをつかって、トポロジカルソートしてみました。

% awk '{ split($1, a, ""); print a[1], a[2], a[2], a[3] }' < keylog.txt | tsort
7
3
1
6
2
8
9
0

参考

池袋バイナリ勉強会(9)に参加しました

12月2日(日)に開催された池袋バイナリ勉強会(9)に参加しました。PDP-11のバイナリ解析に取り組んでいます。a.outの逆アセンブル、およびインタプリタの作成に取り組んでいます。

目標: hello.c を実行するインタプリタの作成

% cat hello.c
main() {
  write(1, "hello\n", 6);
}
% v6cc hello.c
% v6run a.out
hello

mov 命令の実装

mov 命令は、前回の勉強会で勉強したことを、そのまま素直に実装しました。

jsr 命令の実装(実装途中)

jsr 命令は、『pdp40/11 processor handbook』(pdp11-40.pdf)を読んでもさっぱり分かりませんでした。そこで、まずは以下のサンプルで jsr の仕組みを見ていきました。

% cat jsr.s
.globl _printf
.globl _exit

mov $10, r2
mov r2, -(sp)
mov $format, -(sp)
jsr pc, _printf
add $4, sp

mov $0, -(sp)
jsr pc, _exit

format: <value:%d\n\0>

v6asでアセンブルして、v6runで実行します。

% v6as jsr.s
% v6run a.out
000c,0000,0008,0000,0000,0000,sp=0004,pc=0006: invalid operand

エラー(invalid operand)が出ました。何かがおかしいです。逆アセンブルしてみます。

% v6as jsr.s
% v6strip a.out
% pdp11-aout-objdump -d a.out

a.out:     ファイル形式 a.out-pdp11


セクション .text の逆アセンブル:

00000000 <.text>:
   0:	15c2 0008      	mov	$10, r2
   4:	10a6           	mov	r2, -(sp)
   6:	15e6 001a      	mov	$32, -(sp)
   a:	09f7 fff2      	jsr	pc, 0x0     # <---- 0x0 は、おかしい
   e:	65c6 0004      	add	$4, sp
  12:	15e6 0000      	mov	$0, -(sp)
  16:	09f7 ffe6      	jsr	pc, 0x0
  1a:	6176 756c      	add	r5, 72554(sp)
  1e:	3a65           	bit	*-(r1), -(r5)
  20:	6425           	add	(r0)+, -(r5)
  22:	000a           	.word	12

_printfのアドレスが 0x0 になっているのがおかしいようです。この場合、v6ld a.out /lib/libc.a を実行します。

% v6as jsr.s
% v6ld a.out /lib/libc.a
% v6run a.out
value:8

今度はうまく行きました!ちなみに、このa.outに対して再び逆アセンブルを行うと、先ほどの_printfのアドレスが0x0から0x24に置き換わっており、さらに0x24以降の命令列を眺めると、それは v6src/s5/printf.s ファイルの中身がそのまま格納されているのが確認出来ます。

jsr.s における jsr 命令を見る限り、どうやら jsr 命令では _printf を呼び出しているようです。DEC PDP-11 Subroutines には、jsr 命令について、次のように解説されています。

        jsr  register,target   - push register; register <- PC; PC <- target
        rts  register          - PC <- register; pop register

If the specified register is r7 (the PC), then the PC is pushed and popped.

先ほどの jsr.s をもう一度見てみます。

% cat jsr.s
.globl _printf
.globl _exit

mov $10, r2
mov r2, -(sp)
mov $format, -(sp)
jsr pc, _printf
add $4, sp

mov $0, -(sp)
jsr pc, _exit

format: <value:%d\n\0>

つまり、jsr pc, _printf が行っていることは、

  • pc をプッシュ
  • pc に pc の値を代入
  • pc に _printf の値を代入

となります。プログラムカウンタ(pc)にはジャンプ先(_printf)を格納して、処理が終わると(ジャンプ前に保存しておいた)元の場所に戻る仕組みのようです。

jsr 命令の対として rts 命令があり、rts 命令が元の位置に戻るための命令になります。C言語のreturnのイメージです。

(補足): ちなみに、v6strip後にv6ldを実行すると以下のエラーが表示されました。

% v6as jsr.s
% v6strip a.out
% v6ld a.out /lib/libc.a
a.out: No relocation bits
a.out: Relocation error

bne 命令を用いたサンプル: loop.s

% cat loop.s 
.globl _printf
.globl _exit

    mov $10, r2
loop:
    mov r2, -(sp)
    mov $format, -(sp)
    jsr pc, _printf
    add $4, sp
    dec r2
    bne loop

    mov $0, -(sp)
    jsr pc, _exit

format: <%d\n\0>

実行結果です。

% v6as loop.s && v6ld a.out /lib/libc.a
% v6run a.out
8
7
6
5
4
3
2
1

v6runのverbose modeとsyscall mode

v6runコマンドは、-v, -s オプションをとるようになっていて、それぞれのオプションを指定するとsyscallと逆アセンブル結果を出力してくれるようになります。

% v6run
usage: v6run [-r V6ROOT] [-v/-s] cmd [args ...]
    -v: verbose mode (output syscall and disassemble)
    -s: syscall mode (output syscall)

実行例です。

% v6as jsr.s && v6ld a.out /lib/libc.a   
% v6run -s a.out
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001
v0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001
a0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001
l0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001
u0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001
e0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,0008,0000,fff4,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001
:0001,0000,ff6d,ffff,0001,ff60,sp=ff5e,pc=0226: sys indir; 029c
0001,0000,ff6d,ffff,0001,ff60,sp=ff5e,pc=0226: sys write; 02b2; 0001
80001,0000,ff6d,ffff,fff6,ff64,sp=ff62,pc=0226: sys indir; 029c
0001,0000,ff6d,ffff,fff6,ff64,sp=ff62,pc=0226: sys write; 02b2; 0001

0000,ffe8,0008,0000,0000,fff0,sp=fff0,pc=024e: sys exit

ところどころに、v a l u e : 8 の文字が出力されているのが確認できます。sys writeが文字を出力しているようです。

N, Z, V, C フラグ

『pdp11/40 processor handbook』(pdp11-40.dpf)を読むと、たとえば mov の説明には、N, Z, V, C フラグに関する記述があります。

Condition Codes: N: set if (src) < 0; cleared otherwise
                 Z: set if (src) = 0; cleared otherwise
                 V: cleard
                 C: not affected

mov 命令を実行した後、その処理結果に応じて N, Z, V, C の値が変わります。これらのフラグの値がどのように変更されるかは、それぞれ命令ごとに異なります。たとえば、上記の mov 命令の場合、V フラグは必ずクリアされることになります。

はじめてのOSコードリーディング ~UNIX V6で学ぶカーネルのしくみ (Software Design plus)』が出版されます。

以下の本が出版されることを勉強会で教えていただきました。とても楽しみです。

はじめてのOSコードリーディング ~UNIX V6で学ぶカーネルのしくみ (Software Design plus)

はじめてのOSコードリーディング ~UNIX V6で学ぶカーネルのしくみ (Software Design plus)

勉強会後にやったこと

7shiさんのv6runソースをGNU GLOBALで読む。

id:n7shiさんのv6runGNU GLOBALで読めるようにしました。

% gtags -v 
[Sat Dec 15 21:02:15 JST 2012] Gtags started.
 Using default configuration.
[Sat Dec 15 21:02:15 JST 2012] Creating 'GTAGS' and 'GRTAGS'.
 [1] extracting tags of AOut.cpp
 [2] extracting tags of AOut.h
 [3] extracting tags of main.cpp
 [4] extracting tags of Operand.cpp
 [5] extracting tags of Operand.h
 [6] extracting tags of utils.cpp
 [7] extracting tags of utils.h
 [8] extracting tags of VM.cpp
 [9] extracting tags of VM.h
 [10] extracting tags of VM.inst.cpp
 [11] extracting tags of VM.signal.cpp
 [12] extracting tags of VM.sys.cpp
[Sat Dec 15 21:02:15 JST 2012] Done.
% htags -ansx

参考:

v6runのverbose modeでN, Z, V, Cフラグの値を出力する

v6runのソースを改造して、N, Z, V, Cの値を出力するようにしました。VM.cppのdebug()メソッドを書き換えます。

void VM::debug(const std::string &msg)
{
    fprintf(stderr,
            "%04x,%04x,%04x,%04x,%04x,%04x,sp=%04x,pc=%04x: %s\n",
            r[0], r[1], r[2], r[3], r[4], r[5], r[6], prevPC, msg.c_str());

    fprintf(stderr, "N:%d Z:%d V:%d C:%d\n", N, Z, V, C); 
}

さきほどの loop.s プログラムを用いて、N, Z, V, C のフラグの変化を確認してみます。注目ポイントは、bne 命令を実行するときの Z フラグの値です。

サンプルプログラム(ループ回数を3回に変更しました):

% cat loop.s
.globl _printf
.globl _exit

    mov $3, r2
loop:
    mov r2, -(sp)
    mov $format, -(sp)
    jsr pc, _printf
    add $4, sp
    dec r2
    bne loop

    mov $0, -(sp)
    jsr pc, _exit

format: <%d\n\0>

アセンブル結果です。(先頭一部のみ)

% v6as loop.s && v6ld a.out /lib/libc.a 
% v6strip a.out
% pdp11-aout-objdump -d a.out

a.out:     ファイル形式 a.out-pdp11


セクション .text の逆アセンブル:

00000000 <.text>:
   0:	15c2 0003      	mov	$3, r2
   4:	10a6           	mov	r2, -(sp)
   6:	15e6 001e      	mov	$36, -(sp)
   a:	09f7 0014      	jsr	pc, 0x22
   e:	65c6 0004      	add	$4, sp
  12:	0ac2           	dec	r2
  14:	02f7           	bne	0x4       # <--- ここに注目
  16:	15e6 0000      	mov	$0, -(sp)
  1a:	09f7 0226      	jsr	pc, 0x244
<以降省略>

実行結果です。

% v6run2 -v a.out
...<省略>...
0000,ffe8,0002,0000,0000,0000,sp=fff6,pc=0014: bne
N:0 Z:0 V:0 C:0
...<省略>...
0000,ffe8,0001,0000,0000,0000,sp=fff6,pc=0014: bne
N:0 Z:0 V:0 C:0
...<省略>...
0000,ffe8,0000,0000,0000,0000,sp=fff6,pc=0014: bne
N:0 Z:1 V:0 C:0
...<省略>...

『pdp11/40 processor handbook』(pdp11-40.pdf)によると、bne 命令は Z = 0 の場合に pc の値を書き換えるとあります。上のv6run2による実行結果をみると、pc=0014 の bne 命令が 3 回実行されており、最初の 2 回は Z = 0 でループを繰り返しているのが確認出来ました。

Githubにプロジェクトを登録しました

https://github.com/noriok/pdp11

更新履歴

(2012-12-23) Githubのプロジェクトへのリンクを追加しました。

SRM563 div2 easy FoxAndHandleEasy

文字列Sがあり、その文字列Sの任意の位置に同じ文字列Sを挿入する。そうして作られる文字列をSのexpansionと呼ぶとする。文字列S, Tが与えられたとき、TがSのexpansionかどうかを求める問題。

コンテスト中に書いたコードは以下。std::stringの文字列の削除のやり方がすぐに分からず、Dashで調べながら小さなプログラムを書いて確認しようかなと思ったけど逆に時間がかかりそうだったので止めました。

#include <cstdio>
#include <cstdlib>
#include <ctime>
#include <iostream>
#include <string>
#include <vector>
using namespace std;

class FoxAndHandleEasy {
public:
    string isPossible( string S, string T ) {
        const int n = S.length();
        if (n*2 != T.length()) return "No";

        for (int i = 0; i < (int)T.length(); i++) {
            if (i+n >= T.length()) break;

            bool ok = true;
            for (int j = 0; j < n; j++) {
                if (S[j] != T[j+i]) {
                    ok = false;
                    break;
                }
            }

            if (!ok) continue;
            int p = 0;
            for (int j = 0; j < (int)T.length(); j++) {
                if (i <= j && j < i+n) continue;

                if (S[p++] != T[j]) {
                    ok = false;
                    break;
                }
            }

            if (ok) return "Yes";
        }
        return "No";
    }
};

Challenge phaseで他の人のコードを読むと、単純に文字列Sを検索して、見つかったらそれを削除して残りがSに等しいかを確認しているコードが殆どでした。ミスしてそうなコードは見当たらず、そのままChallenge phaseは終わりました。System Testが始まり、多くの人がこの問題を落としていて何事かと思ったら、以下のようなケースで引っかかっていました。

S = aba
T = ababaa

最初に見つかる aba を削除すると残る文字列は baa になり、Sと一致しませんが、2番目に見つかる aba を削除すると残る文字列は aba となり、Sと等しくなります。つまり、削除する文字列の位置が重要なんですね。これは気づきませんでした(気づかずに、Tに含まれる全てのSを検索するようなコードを書いていました)。

コンテスト後に書き直したコード:

class FoxAndHandleEasy {
public:
    string isPossible( string S, string T ) {
        string::size_type p = 0;
        while (string::npos != (p = T.find(S, p))) {
            string x = T.substr(0, p) + T.substr(p+S.length());
            if (x == S) return "Yes";
            p++;
        }
        return "No";
    }
};