ビット演算集大成の問題！？

ビット演算集大成の良問をPOJで見つけたので紹介。

実は、自分でも似たような問題を考えていてTopCoderに問題投稿しようかな・・と思ってましたが、似たような問題は既にありました。。

しかも、2次元バージョンでよりトリッキーな感じ。。

問題はこちら。

http://poj.org/problem?id=2965

簡単に日本語訳すると、

4*4マスに+と-の文字列がある。

（i, j）マスを選択すると、(i,j)マスが

+の場合は、-に

-の場合は、+に

反転変換できる。但し、(i,j)マスを選択すると、i行のマス及びj行のマスの値もすべて反転してしまう。

すべての文字を-にするためには、最低何回のマスを選択しないといけないか？

またそのときに選択するマスを列挙せよ。

んー、これ系の問題はメジャーなんですかね。。TopCoderでも似たような問題ありました。

ビット演算を使うと解けますが、実はこの問題はビット演算の使いどころが２種類あります。

①組み合わせをビットで表現する

②状態および処理をビット演算で計算し、処理を高速化

①は、どのマスを選択するかを0, 1のビットで表すという意味。

②は、4*4マスの状態を16bitのビットで保持します。さらに、あるマスを選択した場合の処理をXORを用いて実現します。

以下ソース。

int main() {
    int ch[16];
    REP(i,16) {
        ch[i] = 0;
        REP(k, 16)
            if (k/4 == i/4 || k%4 == i%4)
                ch[i] |= 1<<k;
    }

    string in;
    int ref = 0;

    REP(i, 4) {
        cin >> in;
        REP(j, 4)
            if (in[j] == '+')
                ref |= 1<<(4*i+j);
    }

    int ret = INF, bcomb = -1;
    REP(comb, 1<<16) {
        int ref_t = ref, cnt=0;
        REP(i, 16)
            if (comb >> i & 1)
                ref_t ^= ch[i], ++cnt;

        if (!ref_t && ret > cnt)
            ret = cnt, bcomb = comb;
    }
    cout << ret << endl;
    REP(i, 16)
        if (bcomb >> i & 1)
            cout << i/4+1 << " " << i%4+1 << endl;

    return 0;
}

クロージャとは何か？(2)

前回に引き続き、クロージャについて。

前回のブログで、クロージャとは”汎用的な機能を持つ関数の機能を何らかの変数で固定して使うような仕組み”と説明しました。

実際に例を見て、クロージャとは何かを更に詳しく調べてみます。

1,2,3, ...., n

という数列の和を求める関数を書いてみます。

def seriesSum(n):
  ret = 0
  for x in range(1, n+1):
      ret += x
  return ret

久々のPython。。

では、次に、

1, 4, 9, ...., n^2

という数列の和を求める関数を書いてみましょう。

def seriesSum(n):
  ret = 0
  for x in range(1, n+1):
      ret += x*x
  return ret

じゃあ、一般的に

f(1), f(2), ..., f(n)

という数列の和を求める関数は書けないものでしょうか？

クロージャを使えば書けます。

こんな感じ。

def genSeriesSum(_f):
  func = _f
  def seriesSum(n):
      ret = 0
      for x in range(1, n+1):
          ret += func(x)
      return ret
  return seriesSum

使うときは、どうするのでしょうか？

例えば、fはx -> x^2の写像とする場合は以下のようにして使います。

def square(x):
  return x*x

squareSum = genSeriesSum(square)
print squareSum(10)

では、説明してみます。

まず、genSeriesSum()について。

この関数は、

まず、関数オブジェクト_fを受け取ります。

そのあと、変数funcに_fを代入しています。

そして、関数内で定義された関数seriesSum()を関数オブジェクトとして返します。

ここでポイントは、seriesSum()の中のfuncという変数は、_fという値に固定されているということです。

つまり、genSeriesSum()は、機能を固定した（funcを_fにした）関数seriesSum()を生成するような関数であると言えます。

そう、これです。

「機能を固定した関数を生成する関数」

これが、クロージャです。

クロージャを使うと、こんなことも出来ます。

関数の微分値を求める関数。

def genDeriverable(_f, _dx):
  f = _f
  dx = _dx
  def deriverable(x):
      return (f(x + dx) - f(x)) / dx

  return deriverable

どうでしょうか？

クロージャ、すごいです。

個人的には、クロージャとは、関数ファクトリーという感覚を覚えました。

私の認識が間違ってたら、是非是非指摘してくださいーー。

クロージャとは何か？(1)

最近よく耳にするクロージャ。

JavaやC++の最新バージョンがクロージャをサポートするとか、しないとか。。

で、結局、クロージャって何なの？

いろいろ調べてみましたが、「なるほど！」と思わせるサイトがいくつかあったので、それを読んで得た知識を紹介します。（間違いなどありましたらご指摘ください。）

wikiをみると、以下のような解説があります。

In computer science, a closure is a first-class function with free variables that are bound in the lexical environment.　Such a function is said to be "closed over" its free variables. A closure is defined within the scope of its free variables, and the extent of those variables is at least as long as the lifetime of the closure itself.

日本語訳。（間違いあるかも。。）

クロージャとは自由変数を持つ第一級関数である。そして、その自由変数はそのスコープにおいて境界となる。このような関数はその自由変数に包括されたと言われる。

クロージャはその自由変数のスコープによって定義され、それらの変数の有効期限は、クロージャそのもののライフタイムと同じ程度である。

これだけ見ても？？な感じですが、具体的な例を見るとその意味が分かります。

クロージャの例や便利な使い方（pythonでの使用例）は次回紹介するとして、今回は、クロージャのポイントをまとめます。

①変数を持ち、かつ、オブジェクトのように扱うことのできる関数である。

②自由変数によってその機能がclosed（固定）される。

③自由変数は、クロージャそのものが存在している間は、記憶され続ける。

どうやら、汎用的な機能を持つ関数の機能を何らかの変数で固定して使うような仕組みを”クロージャ”と呼ぶのではないでしょうか？

逆に言うと、異なる変数で"closed"してあげると、異なる機能を持った関数を簡単に作成することができるということにもなります。

続きを読む..

包除原理とその実装

日曜のSRM、no-contestとは。。悲しすぎる。。

本番は解けなかったHardの問題。包除原理を使うことは分かっていたけど、実装が間に合わなかった。

今日は、簡単な包除原理の復習。英語では、Inclusion-exclusion principle.

そのまんま。

この原理を用いると、集合A1, A2, A3, ... Anがあるときに、すべての集合の和集合の要素数を求めることができます。

一番簡単な例で、A1とA2の2つしか集合が無い場合を考えます。これは簡単。

|A1| + |A2| - |A1 and A2|

ですね。ここで、|x|は集合xの要素数です。

じゃ、A1, A2, A3と集合が3つになったら？？

これはベン図を書いて、少し考えると分かります。

|A1| + |A2| + |A3| - |A1 and A2| - |A1 and A3| - |A2 and A3| + |A1 and A2 and A3|

です。

じゃあ、集合がn個あるときは？を説明したのがInclusion-exclusion principleです。

いろいろ解説サイトがありますが、wikiの説明が1番分かりやすいです。

http://en.wikipedia.org/wiki/Inclusion%E2%80%93exclusion_principle

証明も載ってます。

それでは、実装してみましょう。POJの問題を解いてみます。

ぱっと見簡単そうですが、入力値が10^9のときも時間内で解くためには、工夫が必要です。包除原理を使って解いてみます。

http://poj.org/problem?id=2407

set<int> primes;

void fact(int n) {
 FOR (i, 2, n) {
     if (n % i == 0) {
         primes.insert(i);
         n /= i;
         i = 1;
     }
 }
 if (n > 1) primes.insert(n);
}

int solve(int n, int bit, VI &vec) {
 int on = 0, mul = 1;

 REP(i, 32)
     if (bit >> i & 1)
         ++on, mul*= vec[i];

 return (on % 2 ? 1 : -1) * n / mul;
}

int main() {
 int n;

 while (cin >> n, n) {
     primes.clear();
     fact(n);

     int ret = 0;

     vector<int>vec(ALL(primes));
     REP(bit, 1<<vec.size()) {
         if (!bit) continue;
         ret += solve(n-1, bit, vec);
     }
     cout << n-1-ret << endl;
 }

 return 0;
}

部分集合をビット演算を使って実装するのがポイントでしょうか～。

2^n - 1パターンの部分集合を考えないといけないので、nが大きい場合は何らかの工夫が必要なようです。数学的アルゴリズムの世界は奥が深い。。

C++でParse

POJで見つけた良問。

難しくはないけど、C/C++で如何に文字列をparseするかが練習できる良い問題です。

あと、無駄なキャスト（upper cast）はしないとか、scanfやcinとgetsの使い分けとか学べます。

取りあえず、問題とソースを見てください。

もっときれいなコードが書ける人は随時募集！（たくさんいそう。。）

問題：

http://poj.org/problem?id=2501

ソース：

char input[256];
int main() {
   int h, m, s, v = 0;
   double x = .0;
   double t0 = .0, t1 = .0;

   while (gets(input)) {
       stringstream ss(input);
       string s1 = "", s2 = "";

       ss >> s1 >> s2;
       sscanf(s1.c_str(), "%02d:%02d:%02d", &h, &m, &s);

       t1 = h + m/60.0 + s/3600.0;
       x += (t1 - t0) * v;
       t0 = t1;

       if (s2 != "")
           sscanf(s2.c_str(), "%d", &v);
       else
           printf("%02d:%02d:%02d %.2lf km\n", h, m, s, x);
   }
   return 0;
}

末尾再帰について

ちょっと夜更かしして末尾再帰について勉強。。

下の解説がかなり分かりやすいです。

In traditional recursion, the typical model is that you perform your recursive calls first, and then you take the return value of the recursive call and calculate the result. In this manner, you don't get the result of your calculation until you have returned from every recursive call.

In tail recursion, you perform your calculations first, and then you execute the recursive call, passing the results of your current step to the next recursive step. This results in the last statement being in the form of "(return (recursive-function params))" (I think that's the syntax for Lisp). Basically, the return value of any given recursive step is the same as the return value of the next recursive call.

The consequence of this is that once you are ready to perform your next recursive step, you don't need the current stack frame any more. This allows for some optimization. In fact, with an appropriately written compiler, you should never have a stack overflow snicker with a tail recursive call. Simply reuse the current stack frame for the next recursive step. I'm pretty sure Lisp does this.

ざっくり日本語で言うと、

1. 普通に再帰を書くと、再帰して呼び出した結果を利用して、値を計算するという流れになる。
2. 末尾再帰では、一番最後に自身を再帰する。最初に計算をして、その値を呼び出し先に伝えるという手法を取る。
3. 末尾再帰を用いると、現在使用しているスタック構造をそのまま再帰呼び出し先で再利用することができ、スタックの節約ができる。

3.は末尾再帰の最適化と呼ばれているそうです。実際は、再帰が繰り返しループに置換されるみたい。末尾再帰の恩恵を受けるには、この最適化機能が使用する言語で提供されているかどうかが肝になるみたい。

とりあえず、言葉の説明はこの程度にしておいて、いつものようにサンプルコードを。

int factorial(int n) {
  if (!n)
      return 1;
  return n * factorial(n-1);
}

上のコードは階乗を求めるソースです。普通は、上のように書くでしょう。。

これは、factorial(n-1)の値を求めて、それにnをかけてreturnするという形になっています。つまり再帰呼び出しの結果を用いて計算を実施したのち、returnという形です。

計算されるイメージはこんな感じ。

factorial(5)

5*factorial(4)

5*(4*factorial(3))

5*(4*(3*factorial(2)))

5*(4*(3*(2*factorial(1))))

5*(4*(3*(2*1)))

これに対して、末尾再帰。

int factorial2(int n, int acc=1) {
   if (!n)
       return acc;
   return factorial2(n-1, n*acc);
}

かるく衝撃を覚えるほどのソースコードです。returnするときには、現在の関数のスタックは必要じゃなくなってますね。。ちょっと動的計画を彷彿とさせるような感じですが。。

この場合の計算イメージはこんな感じです。

factorial2(5,1)

factorial2(4,5)

factorial2(3,20)

factorial2(2,60)

factorial2(1,120)

factorial2(0,120)

末尾再帰の場合は、最後に呼び出された関数の戻り値が答えになります。それぞれの関数は自分の役割を終えたら、後は呼び出し先に任せるよ。ってイメージです。

逆に普通の再帰の場合は、最初に呼ばれた関数の戻り値が答えになります。

それぞれの関数は、自分の仕事を部下に任せるんですが、部下が仕事を終えるのをずっと待っていなくてはなりません。一番上の関数は、部下の部下の部下の・・・・部下の仕事が終って、自分のところまで結果が戻ってくるのを待たないといけないわけですね。。これは大変。（面倒くさそう。。日本のソフト開発の主流であるウォータフォールモデルみたい。。）

残念ながら、私が使用しているC++環境では末尾再帰の最適化はサポートされていないようです。末尾再帰を使って書いた関数（総和を求める関数）でもスタックオーバーを起こしました。。

引用サイト：　http://stackoverflow.com/questions/33923/what-is-tail-recursion

Emacs 入門（２）

前回に引き続き、emacs入門。

今回は、

①ソースコード作成

②コンパイル

③実行

をemacs上でやる方法について整理します。

①ソースコード作成

これは、簡単。

コマンドラインから、「emacs ファイル名」を実行すれば、指定されたファイル名でファイルを編集できます。

保存は、「C-x s」

閉じるは、「C-x C-c」

です。「saveのs」と「closeのc」ですね。

②コンパイル

emacs上でコンパイルができます。

「M-x compile」

と打つと、コンパイルのコマンドを入力できます。makefileがあれば、「make -k」と入力しましょう。

③実行

なんと、実行もemacs上で出来ます。ここまで出来ると単なるエディタではなくて統合開発環境に近い気がします。デバッガも起動できるみたいです。

とりあえず、実行は、「M-x shell」とすればよいです。

すると、コマンドラインが現れます。そこで、いつもやってるように、「./exeファイル名」とすれば実行できます。（emacs画面に表示されるコマンドライン上では、lsやgrepなどのコマンドも使えます。）

プログラミング実行の後、ソースコードに戻りたければ、「バッファ」（※参考サイト参照のこと）を切り替えればよいです。

「C-x b」とすると、ウィンドウに表示されるバッファを変更できます。

参考サイト

バッファとは？　http://cns-guide.sfc.keio.ac.jp/2000/4/2/1.html

Emacs 入門

最近（昨日から）、emacsを始めた。

大学１年生のときさわったことがあったが、当時の私はコンピュータアレルギーだったため、意味が全く分からず挫折。。

とりあえず、これだけ使えれば何とかなりそうなコマンドをまとめてみる。

C-は「ctrl」キーを押しながら

M-は「alt」キーを押しながら

という意味。

①C-x C-c 　　 ファイルを閉じる

②C-x C-s 　 ファイルを保存する

③M-x compile コンパイルを実行する

④C-x 2 画面を上と下に２分割する

⑤C-x 3 画面を左と右に２分割する

⑥C-x 1 分割をやめて、フォーカスのある画面だけを開く

⑦C-x o 画面分割時に他の画面にフォーカスを移す

⑧C-k １行削除する

取りあえず、これだけあれば何とかなるだろう。。コピペとかも必要そうだが。。

あと、.emacsというファイルにて初期設定ができるらしい。

このファイルがなければ、自分のホームディレクトリに作ればいいみたい。

とりあえず、下の２つを設定してみた。

①C/C++の予約語に色を付ける

②C/C++のインデントはタブ４つ分にする

以下、設定ファイル。

(global-font-lock-mode t)

(setq-default indent-tabs-mode nil)

(setq-default tab-width 4)

(add-hook 'c-mode-common-hook

'(lambda ()

(c-set-style "k&r")

(setq indent-tabs-mode t)

(setq c-basic-offset 4)

))

(setq auto-mode-alist

(append

'(("\\.c$" . c-mode))

'(("\\.h$" . c-mode))

'(("\\.cpp$" . c++-mode))

auto-mode-alist))

いろいろな距離概念

日常生活で距離というと、「ユークリッド距離」のことを指す。

他にもいろいろと距離があるので、代表的な3つの距離をまとめてみる。

1. ユークリッド距離
2. マンハッタン距離
3. チェビシェフ距離

まずは、「ユークリッド距離」について。

n次元空間中の2点x, yのユークリッド距離は次のように表される。

これは、ベクトル空間でいう距離や３次元空間における一般的な距離と同じである。

次に、マンハッタン距離。

同じく2点x,yの距離を式で表すとこんなかんじ。

これは、マンハッタンのようなブロックで区画されたまちにおいて、北に○ブロック、西に×ブロック移動した言われた場合に、移動しなければいけないブロックの数を表す。

最後に、チェビシェフ距離。

同じく式で。

これは、あんまり日常生活では使われない。

斜めに移動するという概念を縦横に移動するのと同じと考えた距離概念である。

たとえば、チェスでルークは斜め方向に4移動する場合、縦横二回にわけてそれぞれ4マス移動しなければならない。（マンハッタン距離=4+4=8）

これに対して、クイーンなら一気に斜めに移動できる。（チェビシェフ距離=4）

最後に、それぞれの距離概念において、原点から等距離にある点を結んだ場合にできる線を下図に示す。

赤：チェビシェフ距離

青：ユークリッド距離

緑：マンハッタン距離

である。

ビットDPにチャレンジ

今日は、最近よく解いているビットDPについて。

その名のとおりビット演算を用いた動的計画法。

一番分かりやすい例として循環セールスマン問題を考えます。

セールスマンが、営業所を出発し、n個の都市をまわり、営業所に帰ってくる。

n個の都市を最適な順番でまわったとき、どれくらい時間がかかるか？

以下、n個の都市をa,b,c,・・として考えます。

まず、簡単に思いつく解法は全件探索。

a,b,c,d,・・・

a,b,d,c,・・・

a,c,b,d,・・・

・・・・・

とpermutationで探索すると答えはでます。しかし、階乗オーダーとなるため、現実的ではありません。

ここで、

”ある都市から営業所に戻るまでの最短距離は、それまでにどの都市を訪問したかに依存するが、その都市をどの順番で回ったかには依存しない”

ということに気付けばDPが適用できます。

つまり、これまでに辿った経路が

a,b,c,dだろうが、

b,c,a,dだろうが、

c,a,b,dだろうが、

これから先のdから営業所までの最短距離は同じです。

これは、俗に言う「Principle of Optimality」ってやつですね。これでDPが適用できます。

あとは、どの都市を回ったかを覚えておけばいいのですが、これをビットを使って記憶します。

これで、ビットDPの完成です。

例題を解いてみましょう。

http://poj.org/problem?id=2907

typedef pair<int, int>point;
point beeper[12];
int network[12][12];
int dp[12][1<<12];
int n;

void tsp(int p, int bit) {
   if (bit+1 == 1<<n+1)
       return;

   REP(i, n+1) {
       if (bit >> i & 1) continue;
       int mask = bit | (1 << i);
       if (!i && mask != (1<<n+1)-1) continue;

       if (dp[i][mask] > dp[p][bit] + network[p][i]) {
           dp[i][mask] = dp[p][bit] + network[p][i];
           tsp(i, mask);
       }
   }
}


int main() {
   int sc, xsize, ysize;
   int x, y;

   cin >> sc;
   while (sc--) {
       cin >> xsize >> ysize;
       cin >> x >> y;
       beeper[0] = point(--x, --y);
       cin >> n;

       REP(i, n) {
           cin >> x >> y;
           beeper[i+1] = point(--x, --y);
       }

       REP(i, n+1)REP(j, n+1)
           network[i][j] = ABS(beeper[i].first - beeper[j].first)
               + ABS(beeper[i].second - beeper[j].second);

       REP(i,n+1)REP(j, 1<<n+1)
           dp[i][j] = INF;
       dp[0][0] = 0;
       tsp(0, 0);

       cout << "The shortest path has length " << dp[0][(1<<n+1)-1] << endl;
   }

   return 0;
}

tspの部分は、最後に0に戻ってくるような形にしたくて、ちょっとスパゲッティ気味です。。

よく考えると、以下のコードでOKです。。。

void tsp(int p, int bit) {
    REP(i, n+1) {
        if (bit >> i & 1) continue;
        int mask = bit | (1 << i);
        if (dp[i][mask] > dp[p][bit] + network[p][i]) {
            dp[i][mask] = dp[p][bit] + network[p][i];
            tsp(i, mask);
        }
    }
}

注意）但し、ビットDPでも指数オーダーなので、この方法で解けるのはnが小さい場合のみ。

空間認知能力

今日練習した問題で面白いのがあったので紹介。

問題はこちら。

http://www.topcoder.com/stat?c=problem_statement&pm=9954&rd=13506

簡単に訳するとこんな感じ。

立方体をn*n*n個の小さな立方体に分ける。そのあと、適当に小さな立方体を抜き取る。

そして、この立方体の写真を三方向（x軸、y軸、z軸方向）から取る。

3枚の写真を見て、抜き取られた可能性のある小さな立方体の最小個数を求めよ。

最近買ったルーブックキューブと睨めっこ・・・。

結局解けなかった。。この問題は、テクニックというより、生まれ持った才能による部分が大きい気がする・・・・。そういえば、高校のとき「東大に行きたければ、数学は空間系の問題に強くならないとダメだ。」って先生が言ってたのを思い出した。

まー、話を戻して、友達の回答を見て自分なりに書いたコードはこんな感じ。

int removeCubes(vector<string> A, vector<string> B, vector<string> C) {
   // A x-y
   // B x-z
   // C y-z

   int sz = (int)A.size();

   memset(cube, 1, sizeof(cube));
   REP(i, sz) {
       REP(j, sz) {
           if (A[i][j] == 'N') REP(k, sz) cube[i][j][k] = 0;
           if (B[i][j] == 'N') REP(k, sz) cube[i][k][j] = 0;
           if (C[i][j] == 'N') REP(k, sz) cube[k][i][j] = 0;
       }
   }

   // valid-check
   REP(i, sz) {
       REP(j, sz) {
           if (A[i][j] == 'Y') {
               int t = 0;
               REP(k, sz)
                   t += cube[i][j][k];
               if (!t) return -1;
           }
           if (B[i][j] == 'Y') {
               int t = 0;
               REP(k, sz)
                   t += cube[i][k][j];
               if (!t) return -1;
           }
           if (C[i][j] == 'Y') {
               int t = 0;
               REP(k, sz)
                   t += cube[k][i][j];
               if (!t) return -1;
           }
       }
   }

   int ret = 0;
   REP(i, sz)REP(j,sz)REP(k,sz)
       if (!cube[i][j][k])
           ++ret;

   return ret;
}

なるほど！そういうことか！

ちょっと気になったのが、valid-checkのところ。なんか短くできそうな・・。

Editorialに赤い人のコードが載っていたので、それを参考に書いてみた。

こんな感じ。

int removeCubes(vector<string> A, vector<string> B, vector<string> C) {
    int sz = (int)A.size();
    bool a[sz][sz], b[sz][sz], c[sz][sz];

    memset(a, 0, sizeof(a));
    memset(b, 0, sizeof(b));
    memset(c, 0, sizeof(c));

    int ret = sz*sz*sz;
    REP(i,sz)REP(j,sz)REP(k,sz)
        if (A[i][j] == 'Y' && B[i][k] == 'Y' && C[j][k] == 'Y') {
            --ret;
            a[i][j] = 1;
            b[i][k] = 1;
            c[j][k] = 1;
        }


    REP(i,sz)REP(j,sz) {
        if (A[i][j] == 'Y' && !a[i][j]) return -1;
        if (B[i][j] == 'Y' && !b[i][j]) return -1;
        if (C[i][j] == 'Y' && !c[i][j]) return -1;
    }

    return ret;

}

ほう、逆の発想をすることで、判定のところがシンプルに書けるわけか。。結構このコード理解するのに時間かかりました（笑）

ちょっと、この才能の差を埋めるにはかなり時間がかかりそうだな～～～。。。

誤差にまつわるエトセトラ

今日は、浮動小数点小数の数値計算誤差について書こうと思う。


まず、以下のコードを実行してみましょう。

#define EPS 1e-7

int main() {
    // case 1
    double x = 1e16;
    printf("%lf\n", x + 1);

    // case 2
    double y1 = 123456123456.1234588623046875;
    double y2= 123456123456.1234741210937500;

    printf("%d\n", ABS(y1 - y2) < EPS);
    printf("%d\n", ABS(y1 - y2) < EPS * y1);

    // case 3
    double z1 =  1e-1072;
    double z2 =  -1e-1072;
    printf("%d\n", ABS(z1 - z2) < EPS * z1);

    return 0;
}

結果は以下のようになります。


>10000000000000000.000000
0
>1
>0
>1.000000


まず、case1から見ていきましょう。
10000000000000000.000000に1を足しても
10000000000000001.000000にはなりません。
これは、doubleの有効桁数が2^52 ~10^15だからです。最初の15桁以下は丸められてしまいます。逆に言うと、10^15以下の整数ならdoubleで正確に表すことができます。
この辺の話は昔投稿しているので、こちらをご覧ください。


次に、case2です。
これは、おもしろいです。y1とy2をビット列で表すとそれぞれ以下のようになります。


01000010 00111100 10111110 10001110 11110010 01000000 00011111 10011011
01000010 00111100 10111110 10001110 11110010 01000000 00011111 10011100


実は、この2つの数は、doubleの世界では連続する数値です。それにも関らず、(y1-y2)は1e-5程度です。
これは、大変です。連続する値なので、ほんの少しのエラーで、ある値がもう片方の値になりえます。
数値計算を実施する際は、この２つは、同じ値とみなしてよいでしょう。
しかし、ABS(y1-y2) ＜ eps (eps = 1e-7)

なんてやっちゃうと、この2つの値は、同値とは判定できません。

何がまずかったのでしょう・・。

doubleは浮動小数点なので、小数点の位置は固定ではありません。つまり、絶対的な値を誤差の許容範囲とするには無理があります。

ここでは、

ABS(y1-y2) ＜ eps * y1

とするとよいでしょう。

最後にcase3。これは、やっかいです。case2で

ABS(y1-y2) ＜ eps * y1

を判定条件に用いればよいことを示しましたが、case3では、この判定条件は適切ではありません。

同値とみなすべきz1とz2が同値とみなされません。

符号が異なる2つの微量な数値を比較する場合は、case2の場合は適さないようです。

では、どうすればいいか・・・。

case2とcase3を2つ使って判定しましょう。結構面倒ですが・・。。

競技系プログラミングでは、小数が現れないように分母を他辺に掛けるテクニックが有効のようです。また、比較する値の絶対値が比較的小さいものであれば、絶対評価（case2）だけ十分でしょう。

出典）

http://www.topcoder.com/tc?module=Static&d1=tutorials&d2=integersReals

読んで理解したことを書こうと思ったら、ほぼ和訳みたいになりました（笑）

英語が出来る人は、上の原文を読んだ方が分かりやすいはずです。

C++で関数型プログラミング：accumulate編

マルチパラダイム言語として知られるC++。

競技系プログラミングばかりやっていると、手続き型でばかり書いてしまう。

職場では、ぱっと見オブジェクト指向のようなスパゲッティソースを保守している。。

最近、C++で関数型プログラミングっぽいこともできることを発見。

今日は「畳み込み関数」をC++でやってみます。

畳み込み関数は、Pythonでいうreduce()です。。

では、早速コードと実行結果を。

#include<numeric>

#define INF 999999999

int x[] = {1,2,3,4,5,6,7,8,9,10};
int y[] = {15, 12, 99, 27};
int z[] = {10, 20, 150, 100};
int w[] = {1,3,10,100, -12, 2, 4};

int multiply(int x, int y) {
    return x*y;
}

int gcd(int a, int b) {
    if (!b)
        return a;
    return gcd(b, a%b);
}

int lcd(int a, int b) {
    return a*b/gcd(a,b);
}

int myMin(int a, int b) {
    return (a < b) ?  a : b;
}

int myMax(int a, int b) {
    return (a > b) ? a : b;
}

int main() {
    cout << accumulate(x, x+SIZE(x), 0) << endl;
    cout << accumulate(x, x+SIZE(x), 1, multiply) << endl;
    cout << accumulate(y, y+SIZE(y), *y, gcd) << endl;
    cout << accumulate(z, z+SIZE(z), *z, lcd) << endl;
    cout << accumulate(w, w+SIZE(w), INF, myMin) << endl;
    cout << accumulate(w, w+SIZE(w), -INF, myMax) << endl;

    return 0;
}

[実行結果]

55

3628800

3

300

-12

100

accumulate()のシグネチャーはこんな感じ。

T accumulate(InputIterator first, InputIterator last, T init, BinaryOperation binary_op );

1. input は畳み込み対象の開始位置
2. lastは畳み込み対象の終了位置
3. initは初期値
4. binary_opは引数を２つもつ関数のポインタ

です。第４引数を省略すると、operator+()が畳み込み関数として使用されます。

注意ポイントは、accumulate()を使用する場合は、numericをincludeしないといけないという点です。algorithmではないので間違えないように！！

やばいぞ、C++。C++好き度が上がりました。。

笑

トポロジカルソート

今日は、SRMの過去問をトポロジカルソートで解いてみた。

有向グラフのトポロジーに注目したソートですが、トポロジーって一体・・・、

なんか他にも似たような言葉がたくさんあるような・・・。

ちょっとまとめてみます。

• トポロジー　　 幾何学的な相対位置、ネットワークの接続状況など
• エントロピー　物事の煩雑さ
• オントロジー　物事の概念を意味・関係などで定義したもの

で、今回はトポロジカルソートです。

簡単に言うと、どの接点も自分の子孫より先に来るようにソートすることです。

Makeでコンパイルする順序を決めたり、PERTのスケジュール管理で使われるみたいです。

ではさっそく解いてみます。

問題は、こちら。

http://www.topcoder.com/stat?c=problem_statement&pm=9824&rd=12179

ソースは、これ。

定型的なアルゴリズムがあるようですが、とりあえず、自己流で出力辺の無い接点を取り出して、その接点はグラフ集合からのぞいて・・・を繰り返しています。

bool dn[64];
class CorporationSalary {
public:
 long long totalSalary(vector<string> relations) {
     memset(dn, 0, sizeof(dn));

     VI idx;
     int sz = relations.size();
     while (accumulate(dn, dn+sz, 0) < sz) {
         REP(i, relations.size()) {
             if (dn[i]) continue;
             bool lst = true;
             REP(j, relations[i].size()) {
                 if (dn[j]) continue;
                 if (relations[i][j] == 'Y') {
                     lst = false;
                     break;
                 }
             }
             if (lst) {
                 idx.PB(i);
                 dn[i] = 1;
             }
         }
     }

     long long ret = 0LL;
     long long sl[sz];
     fill(sl, sl+sz, 0LL);

     REP(i, sz) {
         long long m = 0LL;
         REP(j, sz) {
             if (relations[idx[i]][j] == 'Y')
                 m += sl[j];
         }
         sl[idx[i]] = m ? m : 1;
         ret += sl[idx[i]];
     }
     return ret;
 }
};

知ってると便利なSTL(5) fill, fill_n

今日は、便利なSTLについて。

• fill
• fill_n

を紹介します。その名の通り、配列やコンテナをある値で”埋める”ことができる関数です。初期化のときに重宝します。

同じような標準関数にmemset()がありますが、こちらは1byte単位で値を初期化するため、charの初期化には有効ですが、intの初期化にはあまり都合がよくありません。（0で初期化するのであれば問題ありませんが、1や-1で初期化するのは面倒。詳しくは、memset()とmemcpy()を参照。）

これに対して、fill()、fill_n()は、その型に対応した単位で初期化を実行できるためとても便利です。

それでは、実際に動かしてみます。

まずは、fill()から。

void fill ( ForwardIterator first, ForwardIterator last, const T& value );

という形で使用し、[first, last)区間の要素をvalueにセットします。

int main() {
    int x[10];
    vector<string> names(7);

    fill(x, x+10, 100);
    fill(names.begin(), names.end(), "ken-ken");

    REP(i, 10)
        printf("%d\n", x[i]);

    REP(i, 7)
        printf("%s\n", names[i].c_str());

    return 0;
}

次に、fill_n。

void fill_n ( OutputIterator first, Size n, const T& value );

という形で使用し、[first, first+n)までの要素をvalueにセットします。

int main() {
    int x[10];
    vector<int> y(10);
    char z[10];

    fill_n(x, SIZE(x), 10);
    fill_n(y.begin(), 10, -1);
    fill_n(z, 10, 'o');

    REP(i, 10)
        printf("%d\n", x[i]);

    REP(i, 10)
        printf("%d\n", y[i]);

    REP(i, 10)
        printf("%c\n", z[i]);

    return 0;
}

要素数が可変なコンテナ（vectorなど）に対してはfill()を、要素数が固定長の配列の場合は、fill_n()を使うといいかと思います。

初めての・・・・Hard AC！

人生初の、TopCoder Hard(Division 2)自力AC。

直観でこの問題は簡単と思って、1時間もあれば・・って思ったけど、落とし穴があって結局3時間かかってしまった（泣）

でも、初めて何にも頼らず自力で解けたというのは大きな進歩だと思う。

これが、赤色への記念すべき歴史的第一歩となることを願うばかりである。。

ヘタレすぎるコードですが、記念ということで貼っておきます。

問題は、こちら。基本的な動的計画なのですが、・・・、同じ高さのところの処理をどうするかでかなり悩みました。解説や赤色の人のコードを参考にスマートな解法を探ってみます。

Single Round Match 404 Round 1 - Division II, Level Three:

http://www.topcoder.com/stat?c=problem_statement&pm=9755

class Prb {
public:
    int sweet;
    int x;
    int y;
    int len;

    Prb(int s, int _x, int _y, int l) {
        sweet = s;
        x = _x;
        y = _y;
        len = l;
    }
    bool operator<(const Prb &prb) const {
        if (this->y != prb.y)
            return this->y < prb.y;
        if (this->x != prb.x)
            return this->x < prb.x;
        return this->sweet < prb.sweet;
    }
};

int dp[64];
bool dn[64];
vector<Prb> pp;

class GetToTheTop {
public:
    int collectSweets(int K, vector<int> sweets, vector<int> x, vector<int> y, vector<int> stairLength) {
        pp.clear();
        memset(dp, 0, sizeof(dp));

        pp.PB(Prb(0, 1, 0, 11000));
        REP(i, sweets.size())
            pp.PB(Prb(sweets[i], x[i], y[i], stairLength[i]));
        SORT(pp);
        // same y
        FOR (i, 1, pp.size()) {
            if (pp[i].y != pp[i-1].y)
                continue;
            int l1 = pp[i].x, r1 = l1+pp[i].len;
            int l2 = pp[i-1].x, r2 = l2+pp[i-1].len;

            double r;
            if (l1 >r2)
                r = dist(l1, 0, r2, 0);
            else
                r = dist(l2, 0, r1, 0);
            if (r <= K) {
                pp[i].sweet += pp[i-1].sweet;
                pp[i-1].sweet = -1;
            }
        }

        FOR (i, 1, pp.size())
            if (pp[pp.size()-1-i].sweet == -1)
                pp[pp.size()-1-i].sweet = pp[pp.size()-i].sweet;

        int posy = 0;
        REP(i, pp.size()) {
            if (posy != pp[i].y) {
                posy = pp[i].y;
                meanSameY(i, posy, pp.size(), K);
            }
            if (i && !dp[i])
                continue;
            int l1 = pp[i].x;
            int r1 = l1 + pp[i].len;

            REP(j, pp.size()) {
                if (i == j) continue;
                if (pp[i].y >= pp[j].y) continue;

                int l2 = pp[j].x;
                int r2 = l2 + pp[j].len;
                double r;
                if (r1 < l2) r = dist(r1,pp[i].y, l2, pp[j].y);
                else if (l1 > r2) r = dist(l1,pp[i].y, r2, pp[j].y);
                else r = dist(0, pp[i].y, 0, pp[j].y);

                if (r <= K)
                    dp[j] = max(dp[j], dp[i] + pp[j].sweet);
            }
        }

        return *max_element(dp, dp+pp.size());
    }

    double dist(int x1, int y1, int x2, int y2) {
        return sqrt((x1-x2)*(x1-x2) + (y1-y2)*(y1-y2));
    }

    void meanSameY(int pos, int h, int sz, int K) {
        VVI idx;
        VI tmp;
        FOR (i, pos, sz) {
            if (i >= sz  pp[i].y != h) {
                idx.PB(tmp);
                tmp.clear();
                break;
            }
            if (pp[i].x - (pp[i-1].x + pp[i-1].len) > K) {
                idx.PB(tmp);
                tmp.clear();
            }
            tmp.PB(i);
        }
        REP(i, idx.size()) {
            int ret = 0;
            tmp = idx[i];
            REP(j, tmp.size())
                ret = max(ret, dp[tmp[j]]);
            REP(j, tmp.size())
                dp[tmp[j]] = ret;
        }
    }
};

ヒープを作ってみよう！

今日は、ヒープを自分で実装してみます。

昔、「二分木は配列で実装できます。」みたいなことを読んだことがありました。

そのときは、？？だったが、今ならできる！！

では、早速ソースを。

class myHeap {
public:
    myHeap() {
        this->tail = 0;
    }

    void add(int n) {
        hp[tail++] = n;
        int pos = tail-1;
        while (pos) {
            if (hp[(pos-1)/2] < hp[pos])
                swap(hp[(pos-1)/2], hp[pos]);
            pos = (pos-1)/2;
        }
    }

    int pop() {
        if (isEmpty()) {
            cerr << "Heap is empty!" << endl;
            return -1;
        }

        int ret = hp[0];
        hp[0] = hp[--tail];
        int pos = 0;

        while (2 * pos + 2 <= tail) {
            if (hp[pos] > max(hp[2*pos+1], hp[2*pos+2]))
                break;
            if (hp[2*pos+1] > hp[2*pos+2]) {
                swap(hp[pos], hp[2*pos+1]);
                pos = 2*pos+1;
            }
            else {
                swap(hp[pos], hp[2*pos+2]);
                pos = 2*pos+2;
            }
        }
        return ret;
    }

    bool isEmpty() {
        return tail <= 0;
    }

private:
    int hp[1<<10];
    int tail;
};

で、ちゃんと動くか確認してみましょう。

int main() {
    myHeap hp = myHeap();

    REP(i, 100)
        hp.add(rand() % 1000);

    while (!hp.isEmpty())
        printf("%d\n", hp.pop());

    return 0;
}

はい、動きました。

実際に自分で作ってみると、ヒープのオーダーについての理解も深まります。

• 最大値の参照は、O(1)
• push()およびpop()は、O(log n)

なのは明らかですね。

部分和を列挙する

今日は、昨日練習したTopCoderの問題について書きます。

SRM 403 devision2 Hard:

http://www.topcoder.com/stat?c=problem_statement&pm=8571

簡単に要約すると、

4と7のみから構成される数字をlucky numberと呼ぶ。

整数nが与えられたときlucky numberのみを足し合わせてnをつくることができる場合は、そのlucky numberを列挙せよ。ただし、そのような組み合わせが複数個ある場合は、個数が最小のもの、さらにそれが複数個ある場合は、辞書順で最小になるものを求めよ。

詳しくは、上のTopCoderのサイトを参照ください。

ここでポイントは、

1. lucky numberを探索してキャッシュすること
2. nがlucky numberのみで構成できるかどうかDPを使って判定すること
3. 構成出来る場合は、その解を問題の定義どおりに列挙すること

2.までは、出来たんだけど・・・

赤い人のコードを参考に復習しました。

まず、1.ですが再帰で書くとかなりスマートに書けます。

2.は部分和問題と同等の方法で書けます。dp[i]には整数iを構成するのに必要なlucky numberの最小値をキャッシュします。（私はdp[i]にiを構成するために必要なluckyu numberをvectorにぶち込んでキャッシュしてました。これだと、メモリ制限や時間制限でoutになりました。。）

3.は2.で作った配列を元に経路復元みたいなことをやります。その際、なるべく小さい数をたくさん取るように（辞書順で小さくなるように）greedyで復元します。

以下、（赤い人を参考に）私が書いたコードを張っておきます。

int dp[1000001];

class TheSumOfLuckyNumbers {
public:
    void makeLucky(vector<int>&A, int x, int n) {
        if (x)
            A.PB(x);
        if (x * 10 + 4 <= n)
            makeLucky(A, x * 10 + 4, n);
        if (x * 10 + 7 <= n)
            makeLucky(A, x * 10 + 7, n);
    }

    vector<int> sum(int n) {
        VI A;
        makeLucky(A, 0, n);
        SORT(A);

        // DP
        EACH(itr, A) {
            dp[*itr] = 1;
            REP(i, n+1) {
                if (!dp[i])
                    continue;
                if (i + *itr > n)
                    continue;
                if (!dp[i + *itr] || dp[i + *itr] >= dp[i] + 1)
                    dp[i + *itr] = dp[i] + 1;
            }
        }

        VI ret;
        int pos = n;
        // Greedy
        EACH(itr, A) {
            while (dp[pos] - dp[pos - *itr] == 1 && pos) {
                if (!dp[pos - *itr] && pos - *itr)
                    break;
                ret.PB(*itr);
                pos -= *itr;
                if (pos - *itr < 0) break;
            }
            if (!pos)
                break;
        }

        return ret;
    }
};

あと、気付きですが、TopCoderの過去問を解くのはPOJを解くよりも良いかもしれません。

その理由は、

1. 問題のレベル分けがきちんとなされている
2. 赤い人のコードや解説(editorial)が読める
3. すべての問題が種類別にカテゴライズされている。

3.のカテゴライズは今日気付きました。

Problem Archiveのところにcategoryってのがあるので、そこで「DP」とか「String manipulation」とか「recursive」とか選択できるみたいです。

http://www.topcoder.com/tc?module=ProblemArchive

知ってると便利なSTL(4) swap

STLを使いこなすことは、アルゴリズムコンテストを勝ち抜くうえでとても重要である。

まだまだ、しらないSTLがたくさんあるようだ。。

今日のSTLはSTL Algorithmsのswap()。

その名のとおり、変数の値をswapできます。

まずは、簡単な例から。

int main() {
   int a = 10;
   int b = 3;

   swap(a, b);
   printf("%d %d\n", a, b);

   return 0;
}

次に、swapを用いて選択ソートをしてみます。

int main() {
    int x[] = {3,4,5,1,2,10,9,7,8,6};

    REP(i, 10)
        FOR (j, i, 10)
            if (x[i] > x[j])
                swap(x[i], x[j]);
    REP(i, 10)
        printf("%d ", x[i]);

    return 0;
}

これは、かなり便利そう。

これから、重宝しそうです。

知ってると便利なSTL(3) max_element, min_element

便利なSTL紹介第3弾！

今日はSTL Containers ではなくSTL algorithmsの紹介。

配列の最大値、最大解（最大値を与える要素）を出力するプログラムを普通に書いてみます。

#define SIZE(buff) (sizeof(buff)/sizeof(buff[0]))
#define REP(i, n) for(int i=0; i<(int)(n); i++)

int main() {
  int x[] = {1,3,5,7, 10, 0, 4, 20, 4, 6};

  int mx = -1;
  REP(i, SIZE(x)) {
      if (mx == -1 || x[mx] < x[i])
          mx = i;
  }
  printf("max element=%d, max value=%d\n", mx, x[mx]);

  return 0;
}

これが何気に面倒臭かったりします。最大値だけならそうでもない（ループ内でmax使えばいいだけ）ですが、最大値を与えるindexを求める場合は、ちょっとだけソースコードが長くなってしまい面倒です。

そんなときは、STL Algorithmsの出番です。

これでOKです。最大値・最大解と、ついでに最小値・最小解も求めてみましょう。

int main() {
   int x[] = {1,3,5,7, 10, 0, 4, 20, 4, 6};

   // search max
   int *mx = max_element(x, x+SIZE(x));
   printf("max element=%d, max value=%d\n", mx - x, *mx);

   // search min
   int *mn = min_element(x, x+SIZE(x));
   printf("min element=%d, min value=%d\n", mn - x, *mn);

   return 0;
}

かなりシンプルですね。戻り値はiteratorですが、配列の場合はポインタとの引き算もできるみたいです。これを使えば、上の例のように最小値（最大値）を与えるindexも簡単に出ますね！

Pythonで数値計算

Pythonのいいところは、何といってもinteractiveなconsoleモードがあるところ。

ちょっとした数値計算をしたいときに重宝する。

例えば、階乗の計算をしたいとき。。

C/C++だと、64bitまでしか計算できない。。JAVAだとBigIntegerが使えるが書くのが面倒。

Pythonならコマンドラインを立ち上げて、サクサクっとかくことが出来る。

こんな感じ。

>>> def f(x):

... if x == 1:

... return 1

... else:

... return x * f(x-1)

...

>>> f(50)

30414093201713378043612608166064768844377641568960512000000000000L

>>>

さらに、Pythonには強力な数値計算モジュールが存在する。

固有値計算なんかも簡単に出来てしまう。。Page Rankの簡易バージョンを実装しようと考えていたが、Pythonだったら簡単に出来てしまう@_@

でも、固有値・固有ベクトルについてはいつか自分でも勉強しなおさないと、と思っている。。

Pythonで固有値・固有ベクトルを求めるとこんな感じ。。

>>> import numpy

>>> A = [

... [1,2,3],

... [4,5,6],

... [7,8,9]

... ]

>>> l,v = numpy.linalg.eig(A)

>>> l

array([ 1.61168440e+01, -1.11684397e+00, -1.23190993e-15])

>>> v

array([[-0.23197069, -0.78583024, 0.40824829],

[-0.52532209, -0.08675134, -0.81649658],

[-0.8186735 , 0.61232756, 0.40824829]])

>>>

楽すぎる（笑）

プログラミング言語には一長一短あるが、

個人的には、

• 基本はC++
• 多倍長、正規表現の文字列処理を用いる場合はJAVA
• 数値計算、ネットワークプログラミングはPython

という感じだろうか。。

PKU150問突破～～！

PKU150問突破した。

やった～～～。わーーい^^

動的計画にだいぶ慣れてきた。それからグラフの問題が解けるようになった。DijkstraとかWarshall-Floydとか基本的なアルゴリズムを勉強した成果だろう。

あとは、Union-Findとか包除原理とかモジュロ演算とかをもっと勉強したいところ。。

アルゴリズムコンテスト用のコーディングスタイルもだいぶ染みついてきた。

• 変数の名前はなるべく1文字にする
• 変数はグローバル変数として宣言する
• なるべく改行しない
• マクロを駆使する

などなど。。

とりあえず次の目標は、200問。そしてターゲットはKyushu Univ. Topのyuta_ihcarokさん。

あとは、TopCoderの過去問やるのもありかなと思う。こっちは過去問解いても順位が付かないからモチベーションは上がらないけど。。

久々の仕事でコーディング

約1年ぶりの仕事でコーディング。

VC++を使った保守開発が始まった。ちょっとした疑問があったのでここに纏めておく。

.h と .libおよび .dllの違いについて。（以下他サイトからの引用）

A header file contains declaration of something (constants, classes, ...), usually ends with a .h or hpp extension.

A DLL (dynamically linked library) is a binary file (with .dll, .ocx, .a, ... extentions), containing functions, resources, ... and can be linked to your program at run-time. In order to use a DLL you need to include the corresponding header file, which declares things in the DLL, so that your program gets compiled. The DLL file is automatically loaded (by system) at run-time when your executable program invokes functions and asks for resources.

A library file (also called static library, with .lib extension usually) is a binary file which alsocontains functions and variables like DLL, but resolved at compile-time, meaning they need to be provided at link-time.

Dynamic and static libraries are often provided together with header file (but no source file, which contains the implementation) when the provider lets you use his/her functions/services but doesnt give you access to the implementation.

簡単にまとめると、

• libファイルは、コンパイル時にリンクされる。実行時にexeファイル単体で実行できるけど、ファイル容量が重くなる。
• dllファイルは、実行時に動的にリンクする。ファイル容量を小さくできるが、実行時に参照しているdllファイルを配置しなければならない。複数アプリから同様の機能を使用する場合は、dllファイルを用いるとメモリの節約ができる。
• hファイルは、宣言のみ書いたファイル。libファイルやdllファイルを使用する場合は、使用するファイル内に含まれるクラス、関数を宣言しておかなければならない。

てな感じでしょうか。

引用元：

http://bytes.com/topic/c/answers/165103-what-relation-difference-between-header-file-dll-library

Pythonで美人時計の画像を自動取得

昨日に引き続きPythonを用いたネットワークプログラミング。

今日は、Pythonでweb上の画像を自動取得するコードにチャレンジ。

まずは、Pythonでネット上の画像を取得します。

ソースはこちら。

import sys
import os
import urllib2

argvs = sys.argv
img_path = argvs[1]
output = "/home/kenji/python/img/" + os.path.basename(img_path)

print "downloading from", img_path, "to", output, "......"

opener = urllib2.build_opener()
req = urllib2.Request(img_path)
req.add_header('Referer', "http://bijint.com/jp")
file = open(output, 'wb')
file.write(opener.open(req).read())
file.close()

これで画像は取得できました。

次に、上記ソースを使って複数の画像をダウンロードすることを考えてみましょう。

ここでやっかいなのは、美人時計の画像は、現時刻および1分前の画像しかアクセスできないというところ。。cronを使えばいけそうですね。。

Unixのdateコマンドで

imgpath=`date '+http://www.bijint.com/jp/img/clk/%H%M.jpg'`

とすれば、その時間の画像ファイルのパスを生成できます。

あとは、この文字をPythonに渡すシェルスクリプトを作成して、cronで毎分実行してあげれば、24*60枚の美人お姉様の画像をダウンロードすることができます。（※注意1)

Windows OSを使っている人は、cygwinのcronでもOKです。

（注意1）著作権絡みの問題があるかもしれないので実際にやるのはやめましょう。あくまでもこの情報はネットワークプログラミングの勉強のための資料としてご一読ください。

PythonでHTML解析

最近、pagerankに関する論文を読んで、ネットワーク関係のプログラミングをしたくなってきました。今日はpythonを使った簡単なhtml分析を紹介します。（本当に簡単ですいません。。）

pagerankは言わずと知れたgoogleで使われているwebページの重要度を割りあてるための手法。

論文を読むまでは、ページに重み付けて、入力リンクの重み付け和を計算してるだけだろうと、ナメてました。忘れてました。彼らはStanfordの学生だったのです。

上の計算を大規模なネットワークに適用するために、グラフ理論と行列・固有値を駆使して高速に解くという学術的なことをちゃんとやってました。。

簡易版のpagerankなら簡単に実装できそうな気がしたので（気がしただけ）、pythonを使って勉強を始めようかと思います。

とりあえず、今日やったのは、

1. yahooのhtml ソースを取得する
2. 自分のbloggerのページから張っているリンク先のページを列挙する

です。

以下ソース。

まず一つ目。

'''
sample 1
Get html contents and automatically decode its strings
'''
import urllib2

res = urllib2.urlopen("http://yahoo.co.jp")
charset = res.headers.getparam('charset')
html = res.read().decode(charset)
print html

urllib2というモジュールを使うと簡単にhtmlコンテンツを取得できます。あとは、文字コードを取得してそれをデコードすればOK。

次に二つ目。

'''
sample 2
Get html contents and get the link from the page
'''
import urllib2
import re

def getHrefAddress(x):
   x = re.sub(r'^href="|^href=\'', '', x)
   x = re.sub(r'"$|\'$', '', x)
   if re.match(r'http://', x) == None or re.match(url, x):
       x = None
   return x

url = "http://techtipshoge.blogspot.com/"
res = urllib2.urlopen(url)
html = res.read()
links = re.findall(r'href=".+?"|href=\'.+?\'', html)
links = map(getHrefAddress, links)
links = filter(None, links)
for link in links:
   print link

これは、ちょっと汚いです。BeautifulSoupというHTMLパーサライブラリがあるのでそれを使うともっとスマートに書けそうです。

上のソースでは、正規表現でリンクを取得し、自ページ（サブページ含む）へのリンクやjavascriptへのリンクは無視しています。それらしい情報を取り出すことができました。

次は、取得したリンクから再帰的にページを辿っていくようなものを作ってみたいと思います。

最短経路アルゴリズムまとめ

今日、グラフアルゴリズム勉強会をした。

第一弾は、最短経路問題を解くためのアルゴリズムについて。

以下まとめ。

	Bellman-Ford	Dijkstra(matrix)	Dijkstra(list)	Warshall-Floyd
Applicable Problem	単一始点最短路問題	単一始点最短距離問題	単一始点最短距離問題	全頂点対最短経路問題
Input Data Form	隣接リスト	隣接行列	隣接行列	隣接行列
Idea	・すべての距離を毎回更新	・確定した接点に接続する接点の距離のみ更新	・priority_queueを用いて確定接点の探索を高速化 ・暫定距離の更新も実在する枝の数だけで実施	・逐次通過可能な接点集合を増加させながら解く
Calculation Amount	・1回のループではE本の枝を走査 ・負の閉路が含まれないければループは高々V-1回で終わる（最短パスに含まれる頂点数は最高でもV-1であるため） ・O(E * V)	・外側のループ一回につき一つの接点が確定 ・ループ内部では、最小値検索、V個の距離更新を行う ・O(V^2)	・値の更新はO(E)で実施できる ・値の取り出しはO(logV)で実行でき、O(E)回程度実行される ・O(ElogV)	・通過できる接点は1つずつ増やす ・隣接行列の値を毎回すべて更新する ・O(V^3)
Feature	・負の辺が含まれていても解ける ・負の閉路がある場合は更新が無限に続くため解けない。 （※無効グラフの場合は辺も閉路になるので注意） ・上記の性質を利用すれば負の閉路を検出できる	・負の辺が含まれている場合は解けない ・Bellman-Fordの無駄な計算を省略した改良バージョンと解釈できる	・負の辺が含まれている場合は解けない ・EはV^2個に比べて小さい場合が多いため、matrixバージョンの無駄を省略し高速に計算できる	・負の辺が含まれていても解ける ・実装が簡単なので計算量に余裕があれば、単一始点問題でもWarshall-Floydを利用するとよい

ダイクストラくらいしか知らなかったが、今回いろいろなアルゴリズムを知ることができて良かった。

使い分けのポイントを簡単にまとめると、

全頂点間の最小距離が必要なときは、Warshall-Floyd。

単一始点で負の辺がある場合は、Bellman-Ford。（※時間に余裕があればWarshall-FloydでもOK。）

単一始点で、負の辺がなく、計算時間の制限が厳しい場合はDijkstra。

この本はかなりお勧め。アルゴリズマーは是非買うべし！！コンピュータアルゴリズムのバイブルになると言っても過言ではないでしょう。

スマートグリッド、そしてスマートシティへ

恥ずかしながら、『スマートグリッド』と初めて聞いた時はグリッドコンピューティングみたいなものだと思っていた。全然違う。
スマートグリッドとは日本語に訳すと『次世代伝送網』のことで、電力の配送システムを供給・需要の両側から制御し最適化しようという試みである。経済面からも環境面からも有益なソリューションであるため、産官学が協力し研究・計画を実施している。
スマートグリッドをさらに発展させて、『スマートシティ』という概念も現れた。これは伝送網だけではなく、街そのものをスマート化しましょうというアイディアである。
簡単に言うと、“自然に優しく、無駄遣いをしないような賢い街づくりをコンピュータを使ってやってみましょう！”ということだ。具体例を挙げれば、

• 太陽光発電を自宅に取り入れて、
• 自動車は電気自動車に変えて、
• 交通渋滞が起こらないように交通量を制御して、
• 配電はもちろんスマートグリッドで、
• 電気の使用状況は常時モニタリングして
• ・・・

という感じ。最近IBMがテレビCMやってますよね。（参考ページ参照）
日本では、横浜（参考ページ参照）、豊田市、京都府、北九州市がスマートシティのモデル都市として選ばれました。うちの会社もがんばってるみたい！ガンバレー！！！

もともとスマートグリッド構想は、オバマ政権がグリーンニューディール政策の柱として打ち出したもの。さらに起源をたどれば、2003年のアメリカ大停電を経験して、停電を起こしにくい配電網を構築しようというのがそもそもの動機付けだったのかもしれない。それが、エコブームの波に乗り、今日の“スマート○○”ブームを作り上げたのかも。。

とりあえず、この『スマートシティ』、IT関連の会社にとってはとてつもないビジネスチャンスです。市場規模は今後30年累計で3100兆円と言われている@_@
私の会社も、スマートグリッドを『IFRS』と並ぶ今後の二大ビジネスチャンスとして捉えているようです。もちろん、IT会社の他にも自動車会社、電力・ガス会社、電気機器メーカーにとっても大きなビジネスチャンスです。今後、“スマート○○”のパワーで日本の、そして世界の景気が上向くことを期待しましょう！
あーーー、給料あがれー、あがれーー笑。

引用:
http://www.kankyo-business.jp/topix/smartgrid_01.html
http://ja.wikipedia.org/wiki/%E3%82%B9%E3%83%9E%E3%83%BC%E3%83%88%E3%82%B0%E3%83%AA%E3%83%83%E3%83%89
http://blog.goo.ne.jp/thinklive/e/05b191a8ecf7f535934185843f745355

参考ページ：
http://www-06.ibm.com/innovation/jp/smarterplanet/cities/
http://www.jftc.or.jp/shoshaeye/contribute/contrib2010_07e.pdf

ショートコーディング(1) : if文はお嫌い！？

最近、ショートコーディングの本を買った。なかなか面白い。

環境依存のテクニックだったり、コンテスト専用のチートコードの作成方法だったりと、必ずしもすべてが皆に有用かと言われるとそうではないが、学ぶべきものはたくさんある。

今日は、その中で面白いと思ったものを紹介します。

とりあえず、次のコードを見てください。

int main() {
REP(i, 10)
    i % 2 && printf("%d\n", i);

REP(i, 10)
    i % 3  || printf("%d サン！\n", i);

REP(i, 10)
    i ? printf("%lf\n", 1.0 / i) : puts("Zero Devider");

return 0;
}

何をやっているか分かりますか？2つ目のループは3の倍数のときだけアホになってます。

ショートコーダーはif文を基本的に書かないらしいです（笑）

確かにこっちの方が短く書けますね！

以下に対応表をまとめます。

分岐文	対応するショートコーディング文
if(hogehoge) do	hogehoge && do
if(!hogehoge) do	hogehoge || do
if(hogehoge) do1 else do2	hogehoge ? do1 : do2

これで、今日からあなたもショートコーダー！！

ビット演算で遊ぼう

Hi, folks. How's it going?

I made an application that deals with some bit operations. Guess it's downright interesting.

If you're unfamiliar with bit operations, I recommend you try the application below.

どうやら、アメリカからもアクセスがあるみたいなので、ちょっとだけ英語でも書いてみました。

今日はビット演算がテーマです。

世の変態天才プログラマーたちは、ビット演算を巧みに扱います。ビット演算の魔術師と言っても過言ではありません。

私も最近ビット演算のすごさに気付きました。ビットというのは2進数なのでいろいろなことに使えます。例えばバイナリサーチの応用のようなこと(※)もビット演算で出来ます。

（※　「のようなこと」と書いたのは少なくとも私のイメージではという意味です。0か1かを常に2つに1つの可能性を選んでデータを構築していて、探索時はO(log n)で出来るという意味です。詳しくは「ビット演算を用いたべき乗計算の高速化」参照のこと。）

『ビット演算を制するものはプログラムを制する！』

『そうか、ビット演算を使えば、コードの量も1/2になり、実行速度も2倍になる。

つまり、おれがビット演算を使えば、その働きは4倍っていうことか！！！』

桜木花道もびっくり@_@

ちょっと話が脱線しましたが、下がソースです。

ビット演算初心者向けです。かなり基本ですが動きを自分で見て確認するにはいいかと。。

あと、どのような演算を施すと何が起きるのかをしっかり把握できます。一度自分で同じようなものを作ってみるといいでしょう。

void getBin(int n, char bin[]) {
    REP (i, 16)
        bin[15 - i] = (n >> i & 1) + '0';
    bin[16] = '\0';
}

int main() {
    int n, bit;
    char calc, bin[16 + 1];

    cout << "Input an integer." << endl;
    cin >> n;
    getBin(n, bin);
    printf("%s\n", bin);

    // Input one of the manipulations below:
    // q       exit
    // u bit   change bit-th bit to '1'
    // d bit   change bit-th bit to '0'
    // x bit   reverse bit-th bit
    // c       reverse all bits
    // where q, u, d, x and c are characters themselves, bit is an integer and 0-based.
    while (cin >> calc) {
        switch (calc) {
        case 'q':
            cout << "Bye!" << endl;
            return 0;
        case 'u':
            cin >> bit;
            n |= 1 << bit;
            break;
        case 'd':
            cin >> bit;
            n &= ~(1 << bit);
            break;
        case 'x':
            cin >> bit;
            n ^= 1 << bit;
            break;
        case 'c':
            n = ~n;
            break;
        default:
            cerr << "Syntax Error." << endl;
            cerr << "Try again." << endl;
            continue;
        }
        getBin(n, bin);
        printf("%s\n", bin);
    }
    return 0;
}

ポワソン分布

今日、ポアソン分布という分布を知った。

高校の物理でやったような記憶もあるが・・・(ポアソン方程式だった。。しかもやったのは大学。)

まず、以下の問題を考えてみてください。

ある女の人は一日平均2回ナンパされるそうです。（とっても美人です。）

この女の人が、今日3回ナンパされる確率を求めなさい。

この問題、簡単そうで普通には解けません。（少なくとも私には解けませんでした。）

一日にx回ナンパされる回数をP(x)とおくと、

とりあえず、ここまでは分かる。

で、P(3)は？？・・・・・・

そこで登場するのがポアソン分布。ポアソン分布が使えるのは以下の条件（ポアソン過程）を満たす場合のみらしい。

1. 事象はいかなる時点でもランダムに発生しうる。
2. 与えられた時間区間での事象の発生は，それと重複しない他の区間に対して独立である。
3. 微小時間⊿tにおける事象の発生確率は⊿tに比例して小さくなっている。
4. 微小時間⊿tの間に事象が２回以上発生する確率は無視できる。
5. 時間tの間に当該事象が発生する平均発生回数λがおおむね5以下である。

でこれがポアソン分布の式。

数学的な裏付けは分からないが、この式を使用すれば単位時間あたり平均λ回発生する事象が単位時間内にx回発生する確率を求めることができる。

この式考えた人天才すぎます。。。

先ほどのナンパの問題だと、

となります。すごい、すごすぎる、この女。じゃなくてMr. ポアソン。

いろいろなλの値に対するポアソン分布のグラフがwikipedhiaにありましたので、下記の参考サイトからご覧ください。

引用サイト：

http://next1.cc.it-hiroshima.ac.jp/MULTIMEDIA/prob/node14.html

参考サイト：

http://ja.wikipedia.org/wiki/%E3%83%9D%E3%82%A2%E3%82%BD%E3%83%B3%E5%88%86%E5%B8%83

TopCoder惨敗とPKU祝杯！

昨日開催されたTopCoder SRM 483（division 2）.

Easyに手間取り8分くらい使ってしまった。良く見ると、途中で変数の名前が入れ替わってた（笑）焦りは禁物ですね。

Middleは、20分かからないくらいで解けた。これは調子がいい。

初めての本格的なHardの挑戦。なんか行ける気がする。。結局行けなかった・・・・。

そして、チャレンジフェーズ。まさかのMiddle撃沈。ショックを受けている間に、次々と撃沈されるソース達。いったい何が・・

最終的には、部屋の中の上2人以外はMiddleが撃沈されている始末！これは前代未聞の出来事では。。

撃沈されていない人のソースをみると、、、、なるほど1人のときは行と列のうち片方見るだけでいいのか。。確かに。。次回からは、例外っぽいケースはきちんと確認するようにしよう。

そして、PKU。

今日、ついに念願の100問、突破～～～。。

意外と早く達成できて驚いている。次の目標は150問。しかし、そろそろ簡単な問題が無くなってきている。ここからの問題はなかなかの猛者達ばかりだろう。

ワンピースで言うと『グランドライン』到達くらい。。

がんばろう！！targetは東北大学のmatsuさん。。すぐに抜かす。

エラトステネスの脅威

今日は、素数を列挙するアルゴリズムを考えます。

最初は基本的な実装を行い、徐々に高速化していきます。最終的には『エラトステネスの篩』を実装します。

また、それぞれのアルゴリズムのオーダーを見積り、実際に実行速度を計測します。

まずは、基本的な実装から。これ。

void getPrime1(int n) {
    VI prime;

    FOR (i, 2, n) {
        bool isPrime = true;

        FOR (j, 2, i) {
            if (i % j == 0) {
                isPrime = false;
                break;
            }
        }
        if (isPrime)
            prime.PB(i);
    }

    printf("%d\n", prime.size());
}

まー、だれでも思いつくような方法。

少し考えてみると、整数iが素数か否かを判定するのにi-1までの割り算をしなくても、sqrt(i)までの割り算をすればよいことに気付く。これをふまえて、改良するとこうなる。

void getPrime2(int n) {
    VI prime;

    FOR (i, 2, n) {
        bool isPrime = true;

        FOR (j, 2, (int)sqrt(i) + 1) {
            if (i % j == 0) {
                isPrime = false;
                break;
            }
        }
        if (isPrime)
            prime.PB(i);
    }

    printf("%d\n", prime.size());
}

ちょっと変わっただけですが、かなり速くなりそうです。

次に、「エラトステネスの篩」を用いた実装。（アルゴリズムの考え方については、参考ページ参照のこと。）ソースを見てもらえば何をしているか分かると思います。

void getPrime3(int n) {
    bool *isPrime = new bool[n];
    VI prime;

    memset(isPrime, 1, sizeof(bool) * n);
    for (int i = 2; i < sqrt(n) + 1; i++) {
        if (!isPrime[i])
            continue;

        for (int j = i + i; j < n; j += i)
            isPrime[j] = 0;
    }
    FOR(i, 2, n)
        if (isPrime[i])
            prime.PB(i);

    delete [] isPrime;
    printf("%d\n", prime.size());
}

ここまで、3つの方法を記述しました。それぞれの計算量を見積ってみましょう。そして実際にn=1,000,000のときの実行速度を計測してみましょう。

#include <sys/time.h>

inline double gettimeofday_sec() {
    struct timeval t;
    gettimeofday(&t, NULL);
    return (double) t.tv_sec + (double) t.tv_usec * 1e-6;
}

int main() {
    double s = gettimeofday_sec();
    getPrime1(1000000);
    printf("%lf\n", gettimeofday_sec() - s);

    s = gettimeofday_sec();
    getPrime2(1000000);
    printf("%lf\n", gettimeofday_sec() - s);

    s = gettimeofday_sec();
    getPrime3(1000000);
    printf("%lf\n", gettimeofday_sec() - s);

    return 0;
}

結果は以下のとおりになりました。

アルゴリズム	オーダー	実行時間[s]
getPrime1	n^2	157.9
getPrime2	n^1.5	6.943
getPrime3	n	0.016

やっぱりエラトステネス速いですね。。getPrime1, getPrime2のオーダーは見てのとおりですが、getPrime3（エラトステネス）のオーダーは少し分かり難いので少し補足します。

まず、一番外側のループは見てのとおりO(n^0.5)です。次に内側のループですが、n/i回まわります。

iによって値が違うので平均値を考えます。i は 2 から n^0.5 までなので、内側のループの回数は、それぞれのiに対して、n/2, n/3, n/4, ...... n/n^0.5となります。この平均は、およそn/{(2+n^0.5)/2}なので内側のループの回数はO(n^0.5)です。よって2重ループ部分のオーダーは、O(n^0.5) * O (n^0.5) = O(n)です。最後にisPrimeを線形探索する部分もO(n)なので、最終的にエラトステネスのオーダーは、O(n)となります。

参考ページ:

エラトステネスの篩について　http://ja.wikipedia.org/wiki/%E3%82%A8%E3%83%A9%E3%83%88%E3%82%B9%E3%83%86%E3%83%8D%E3%82%B9%E3%81%AE%E7%AF%A9

多倍長演算 --doubleの精度!?

今日は、ついに多倍長演算デビュー。

PKU4問解いてやった。C++で自力で実装しようと思ってたが結局Javaでやることにした。信頼性の高い既成品があれば、当然そちらを使うべき。（あー、いつの間にかビジネスチックな思考になってしまった。。）

Javaは2種類の多倍長演算に対応する型を装備している。

まずは、BigInteger。その名の通り大きな整数。これはC++でいうlong long intの最大値より大きな値を扱いたい場合に使用する。long long intの最大値は、(2^64) / 2 - 1でだいたい8.0 * 10^18。

話はちょっとそれるが、『2^10 ≒ 10^3』と覚えておくと累乗の計算はすぐできるので便利。

2つ目は、BigDecimal。これはdoubleより精度の高い小数を扱うときに用いる。

BigIntegerの概念は分かりやすいが、BigDecimalは・・・。はて？？doubleの精度ってそもそもどのくらいだったっけ？

ということで、今日は、

1. doubleの最大値、最小値
2. doubleの精度
3. 丸め誤差

について、実際にPC上でプログラムを走らせながら検証します。（ちょっと長くなりそうだな～）

その前に、doubleの内部構造についておさらいします。doubleは1ビットの符号部、11ビットの指数部、52ビットの仮数部から成り立っています。

（詳しくは、参考ページをご参照ください。）

簡単に言うと10進数の指数表記（2.999 * 10e+10みたいなやつ）と同様のことを2進数でやってるだけです。

では、準備もできたので、まずは、doubleの最大値、最小値を求めてみましょう。指数部が11ビットなので概算では大体2^1024くらいがdoubleで表すことのできる最大値です。だいたい10^300くらいですねー。最小値は、2^(-1024)なので、10^(-300)程度だと思います。

では実際に計算しましょう。

#include <float.h>

#define EXP_MAX 2046     // 2047 is booked for exception procedure
#define EXP_MIN 1        // 0 is booked for exception procedure
#define EXP_BIAS 1023    // The biased exponent range: [-1022, 1023]
#define SIG_BIT 52       // significand bits

int main() {
   // max
   double dblmax = pow(2.0, EXP_MAX - EXP_BIAS);
   double significand = 1.0;
   REP(i, SIG_BIT)
       significand += pow(2.0, -(i+1));
   dblmax *= significand;
   printf("dblmax = %e, DBL_MAX = %e\n", dblmax, DBL_MAX);

   // min
   double dblmin = pow(2.0, EXP_MIN - EXP_BIAS);
   printf("dblmin = %e, DBL_MIN = %e\n", dblmin, DBL_MIN);

   return 0;
}

んー、実際は、最大値が1.797693e+308、最小値が2.225074e-308でした。ほぼほぼ予想通りです。指数部の値は、0と1024が例外処理のため予約語となっているようです。（参考ページ参照のこと。）やられた。これに気付かず1時間ぐらい悩んでしまった。。。

次に、doubleの精度を考えてみましょう。指数部が数値のスケール、仮数部が数値の精度を表しているので仮数部に注目すればOKです。52ビットなので、およそ10^15程度（有効数字15桁）の精度は持っているのではないでしょうか。（経験的にはそんなに無い気が・・・）実際にPC上で確認してみましょう。

int main() {
   double min_significand = pow(2.0, -SIG_BIT);
   printf("%.16f\n", min_significand);

   double x = 1e15 + 1.0;
   double y = 1e16 + 1.0;
   printf("%.1lf\n", x);
   printf("%.1lf\n", y);
}

んーと、仮数部の最小値は2.0 * 1e-16程度でした。スケールに対して2.0 * e-16程度の大きさの値までしか表現できないという意味なので、有効数字は15桁と考えていいでしょう。

実際に上のソースで検証しましたが、1.0 e+16に1を加算しても無視されてしまいます。

あれ、そうか！これが世間で言われている『情報落ち』ってやつか！！

最後に、「丸め誤差」について。言わずと知れた「10進数を2進数に変換するときにでる誤差」のことです。代表的な例に、0.1を10回足しても1にならないというのがあります。

今回は、0.1をdoubleで表したときにどのような表現になるのか、またその値は10進数ではどの程度なのか調べたいと思います。標準関数frexp()を使えば簡単にできるようですが、なかなか面白そうなテーマなので、今回は自力で実装します。

ポイントは、『仮数部は52ビットあるので全件探索ではNG。何らかの”工夫”が必要』ってことです。まーでもこの場合は工夫と言ってもアレしかないですね（笑）

double getSigVal(LL sigbit) {
   double ret = 1.0;

   REP(i, SIG_BIT) {
       if ((sigbit >> SIG_BIT - 1 - i) & 1)
           ret += pow(2.0, -(i+1));
   }
   return ret;
}

double getDoubleVal(int exp, LL sigbit) {
   return getSigVal(sigbit) * pow(2.0, exp - EXP_BIAS);
}

int main() {
   double x = 0.1;
   LL sig, sigl, sigr;
   int expo = 1;

   while (pow(2.0, expo + 1 - EXP_BIAS) <= x)
       expo++;

   sigl = 0LL;
   sigr = (1LL << 52) - 1LL;
   sig = (sigl+ sigr) / 2;

   REP(i, 200) {
       double val = getDoubleVal(expo, sig);
       if (ABS(val - x) < 1e-18)
           break;
       else if (val > x)
           sigr = sig;
       else
           sigl = sig;
       sig = (sigl + sigr)/2;
   }


   printf("significand = %lld, exponent = %d\n", sig, expo);
   printf("0.1 is expressed something like %.18lf on the PC\n", getDoubleVal(expo, sig));
   return 0;
}

52ビット（10進数で10^15程度）を線形探索していたのでは、時間がありません。Binary Searchで解きましょう。精度は1e-18程度でいいでしょう。ループは200回くらい回せば探索範囲が1e-60倍くらいに縮小するので200回くらい回せばこの場合は十分です。（Binary Searchはバグを生みやすいので、無限ループはやめて必ず有限回の探索で打ち切るようにしましょう。）

参考ページ:

doubleの内部構造

http://www.cc.kyoto-su.ac.jp/~yamada/pB/float.html

指数部の例外処理予約値

http://www.asahi-net.or.jp/~uc3k-ymd/Lesson/Section01/section01_06.html

配列を関数に渡す

C/C++を使用する人なら、値渡しとアドレス渡しがあるのは知っていると思います。

しかし、配列は値渡しをしたつもりでもアドレス渡しになってしまうようです。

また、配列を関数に渡す場合には、どうすればいいか迷うことも多いでしょう。特に多次元配列の場合はそうです。

まず、配列が他の型とは違う挙動をすることを実際にアドレスを表示させることで確認してみましょう。

/*
* int address
*/
void func1(int n) {
   printf("The address of input is %d.\n\n", &n);
}

/*
* buffer address
*/
void func2(int n[]) {
   printf("The address of input is %d.\n\n", n);
}

/*
* STL address
*/
void func3(vector<int> n) {
   printf("The address of input is %d.\n\n", &n);
}

/*
* class address
*/
class Human{
public:
   int age;
   string name;
};
void func4(Human n) {
   printf("The address of input is %d.\n\n", &n);
}

int main() {
   int a;
   printf("The address of input is %d.\n", &a);
   func1(a);

   int b[10];
   printf("The address of input is %d.\n", b);
   func2(b);

   vector<int> c;
   printf("The address of input is %d.\n", &c);
   func3(c);

   Human d;
   printf("The address of input is %d.\n", &d);
   func4(d);


   return 0;
}

配列の場合のみ、呼び出し元と呼び出し先で同じアドレスが出力されます。つまり、配列は値渡ししたつもりでも実際はアドレス渡しになるということです。もちろん関数のプロトタイプ宣言部のint n[]をint *nと書くことも出来ます。こちらの方がアドレスを渡していると明示的に分かるので良いかもしれません。

では、次の問題。多次元配列を関数に渡す時はどうするか。

多次元配列は実際は一次元方向にメモリーが展開されるということを知っておくことが大切です。実際に見てみましょう。

int main() {
  int x[3][3][3];

  REP(i, 3)
      REP(j, 3)
          REP(k, 3)
              printf("%d\n", &x[i][j][k]);

  return 0;
}

どうでしょうか？

メモリが32bit intの環境では4byteずつ増えていると思います。

つまりメモリ上は多次元配列は一次元配列と同じ形をしているということです。

多次元配列を宣言した関数内では、それぞれが何次元で構成されているか分かりますが、関数に渡してしまうと呼び出しさきでは、何次元目がどれくらいの単位で区切られているか分かりません。

よって、関数に渡す際には、呼び出し先の関数のプロトタイプに配列のサイズを明記しなければなりません。

こんな感じ。

void func5(int n[][10]) {
  REP(i, 10)
      REP(j, 10)
          printf("%d%c", n[i][j], (j==9) ? '\n' : ' ');

  return;
}

int main() {
  int buf[10][10];
  memset(buf, 0, sizeof(buf));
  func5(buf);

  return 0;
}

今日の話は基本中の基本って感じですが、やっぱり基本はしっかり押さえておきたいものです。

IFRS対応でITバブル！？

今日は巷を賑す（！？）「IFRS」について書こうと思います。

まず「IFRS」って何ぞや？ということで、簡単に説明します。

International Financial Reporting Standards（ 国際会計基準）の略で、「アイファース」と読まれることが多いです。

読んで字の如く、「企業の財政状態や経営成績を表すための会計基準の世界標準」のことです。

現在、日本の企業は日本独自の会計基準に基づいてレポーティングを実施していますが、近い将来このIFRSが日本の上場企業に強制適用される可能性があります。

2012年を目途に、金融庁によって上場企業へのIFRS強制適用の是非が判断され、強制適用が決定した場合は早ければ2015年よりIFRSに準拠したレポーティングが必須となります。

お試し期間（？）として、2010年3月末より、一定の要件を満たす上場企業の連結財務諸表についてIFRSを任意に適用できるようになりました。既にIFRS取り入れている企業もあるみたいです。。

IFRSの適用されることになれば、決算・財務報告のプロセス及び関連する情報システムの大規模な改革が必要となります。

しかし、会計系システム以外にも影響を受ける部分は多くあります（下表参考）。

影響を受ける業務	変更が必要なシステム
売上計上基準の変更	販売管理システム
棚卸資産の範囲/原価算定方法の変更	購買管理システム/在庫管理システム
固定資産の減価償却の扱いの変更	固定資産管理システム
リース取引の根本的な変更	リース管理システム
未消化有給休暇の費用計上	人事給与システム

いや、これは大変ですね。。(ちょっと難しい経済用語が並んでいるので、その説明については参考サイトを参照下さい。)

上場企業にとっては、お金も時間もかかる一大事ですが、System IntegraterやIT Cosulting Firmにとっては大きなビジネスチャンスとなります。

実際に私の会社でもIFRSは今後もっとも重要視すべきビジネスチャンスとして捉えています。

2010年-2015年にかけてITバブルの再来を期待します。給料あがらないかなー（笑）

引用サイト：

http://www.fsa.go.jp/news/21/sonota/20100423-2.html

http://otona.yomiuri.co.jp/life/kabunushi/kabu100319.htm

http://www.ntts.co.jp/publish/column/ifrs/02/

参考サイト：

棚卸資産とは？http://tax.xrea.jp/auction/tanaoroshi.html

原価とは？http://www.sk-genka.com/kanri/kan02_genka.html

減価償却とは？http://itpro.nikkeibp.co.jp/word/page/10011680/

知ってると便利なSTL(2) pair

intとintのペアを扱いたいとき、ついつい自分でクラスを作ってしまっていた。

何が面倒かというとソート。いちいちoperatorを定義してあげないとsort()が使えない。

最近、pairというSTLが便利だということを知った。

2つの型のペアを格納でき、さらに比較演算子が定義されているため、sort()が使える。かなり便利。

比較演算では、pairの第一要素の比較結果が返される。第一要素が同じ値だった場合は、第二要素の比較結果が返される。

唯一不便そうな、make_pairも下のようなマクロを定義してやれば楽々！！

#define MP make_pair

まー、いつもの如くコードで説明します。誕生日を格納するために、pairを使用しています。そして月日の順でソートをして表示します。

int main() {
vector<pair<int, int> >birthday(20);

REP(i, 20)
 birthday[i] = MP(rand() % 12 + 1, rand() % 30 + 1);

SORT(birthday);

REP(i, 20)
 printf("%02d/%02d is My Birthday.\n", birthday[i].first, birthday[i].second);
return 0;
}

いやー、これはかなり便利だねー。

更に、性と名をpairに格納することも可能。所謂”名前順”ってやつでソートしてみましょう！

int main() {
vector<pair<string, string> >name(10);
name[0] = MP("Tanaka", "Ai");
name[1] = MP("Yamada", "Taro");
name[2] = MP("Yamada", "Hanako");
name[3] = MP("Okada", "Yusuke");
name[4] = MP("Yamaguchi", "Satoshi");
name[5] = MP("Suzuki", "Ichiro");
name[6] = MP("Suzuki", "Jiro");
name[7] = MP("Suzuki", "Saburo");
name[8] = MP("Suzuki", "Shiro");
name[9] = MP("Oda", "Nobuo");

SORT(name);

EACH(itr, name)
  cout << "I'm " << itr->first << " " << itr->second << endl;
return 0;
}

何か気付きましたか？

そうです、鈴木が多いです。

じゃなくて、stringを使っています。const char*ではダメです。

なぜなら、比較演算子が定義されていないからです。（const char*の場合はstrcmp()で比較しますよね！）

TOEIC受験

今日TOEICを受けた。

1年半ぶりだが、そのときよりも良い結果がでそう。。でも正直出来はいまいちだった。。

まず、リスニング。

問題が上から下に設問が書いてあるのに、左から右に設問を読んでしまい、3問ミス笑。。

いやーヘタレすぎる。。

部分的には対策をしたが、きちんと過去問解いておくべきだった。

そして、リーディング。

かなり読むのが速くなってた。でも、最後まで解けなかった。

多分原因は、選択肢をすべて確認していたこと。

選択肢Aが正解だと分かっても、確認のためB, C, Dも吟味していた（笑）

これが無かったら全部解けたかも。。

同じく、過去問を解いておくべきだった。

今回学んだことは、

『何事もリハーサルが必要である。そして、そのリハーサルは可能な限り本番と同じ環境を作り上げて行うこと。』

それから、TOEICのTipsとして以下2点。

1. ちゃんと問題の番号を確認する！！（超基本（笑））
2. これが正解と思ったらすぐにマークする！！

まー、900点取れてなかったら、来年もう一回受けてみようかと思います。

コンマ演算子

世の中には「ショートコーダー」と言われる人達がいる。

コードを短く書くことすべての情熱を注ぐ変態天才たちである。

彼が重宝するテクニックが、「コンマ演算子」。

これを使うと、

1. ループの終了判定が標準入力から入力したものに依存するとき、breakを使わずに簡単に書ける。
2. ループの括弧を省略することができる。（Pythonみたい。）

という利益を享受できる。

すばらしい。。。

百聞は一見に如かずということでソースコードを。

例によって、PKUの問題から。

私の回答はこんな感じ。

int main() {
int L, M, s, e;
while (scanf("%d %d", &L, &M), L++ && M) {
REP(i, M)
    scanf("%d %d", &s, &e),
    L -= e-s+1;
printf("%d\n", L);
}
return 0;
}

注目すべきは、「コンマ演算子」。

まず、whileの中について。ここでコンマを使用することにより、複数の処理を判定条件の中に埋め込むことができます。判定条件として判断されるのは、一番最後の部分（この場合は、L++ && M部分）です。

次に、REP()の中ループについて。ループの内容は2行ですが、括弧を使用していません。

これも「コンマ演算子」によってなせる技。。コンマで繋がれた処理は1行の処理とみなされます。

これは、是非ともマスターしたいテクニック！！！

でもやり過ぎは厳禁。業務でこんなの使いまくってたら殺されますね（笑）