Ackermann function

アッカーマン関数という関数があります。
バカ正直に計算しようとするとすぐ計算量が膨大になり、計算できなくなるたぐいの関数です。

これをアセンブリでバカ正直に計算させるコードを書いてみました。階乗を計算させるコードのおまけ付きです。

実際に実行させたところ、ack(3, 12) あたりになるとスタックがあふれました。
スタックサイズを調整すればもうちょっといけると思います。

使用した高速化のテクニック

• 末尾再帰を CALL から JMP に置き換え
• TEST+JE でmacro fusion できる命令形式にした
• INC / DEC を ADD / SUB に置き換えて、フラグレジスタの依存関係を解消
• MOV RDX, 1 を MOV EDX, 1 にして命令長を短縮

    .code

ack PROC public

    ; RCX = m, RDX = n
    TEST    RCX, RCX
    JE      m_zero

    TEST    RDX, RDX
    JE      n_zero

    PUSH    RCX
    SUB     RDX, 1
    CALL    ack
    POP     RCX
    SUB     RCX, 1
    MOV     RDX, RAX
    JMP     ack

n_zero:
    MOV     EDX, 1
    SUB     RCX, 1
    JMP     ack

m_zero:
    MOV     RAX, RDX
    ADD     RAX, 1
    RET

    ack ENDP

factorial PROC public
    
    MOV     EAX, 1
    CMP     RCX, 1
    JLE     f_ret

f_loop:
    IMUL    RAX, RCX
    SUB     RCX, 1
    CMP     RCX, 1
    JNLE    f_loop

f_ret:

    RET

    factorial ENDP

    END

一応 C で書いた以下のコードよりは、3割ぐらい速いようです(Visual C++ 2010 で /O2 としてコンパイル)。階乗は微妙に速いだけでほとんど変わりませんでした。CPU は Penryn ベースの Core2 Duo です。

int64_t factorial_c(int64_t n)
{
    int64_t r = 1;
    while(n > 1) r *= n--;
    return r;
}

int64_t ack_c(int64_t m, int64_t n)
{
    return m == 0 ? n + 1 : n == 0 ? ack_c(m - 1, 1) : ack_c(m - 1, ack_c(m, n - 1));
}

Bogosort

世の中にはボゴソート というソートアルゴリズムがあります。ランダムにシャッフルして、順番に並んでいたらソート完了というものです。

魔が差したので実装してみました。

ところで、ボゴソートは完全な乱数を使用した場合にはアルゴリズムが完了することが保証されていますが、疑似乱数には周期があるため同じ状態が繰り返し、完了しない場合があります。メルセンヌツイスタ MT19937 を利用した場合、周期は 2^19937-1 で、およそ 10^6001 です。

一方ボゴソートの平均計算量は O(n * n!) らしいので、138 * 138! = 9.55 * 10^238、つまり 138 要素ぐらいまでならメルセンヌツイスタぐらいの周期を持つ乱数でそろう可能性が高い(50%を越える)ということになります。だいたい計算してみると 2080 * 2080! = 10^6004 ぐらいとのことなので、2079 要素ぐらいまでならメルセンヌツイスタぐらいの周期を持つ乱数でそろう可能性が高い(50%を越える)ということになります。

32bit の線形合同法(普通の rand() 関数)の場合は、最大周期が 2^32-1 なので、わずかに11要素ぐらいまでしかソートできず、それを越えると永久にソートできない可能性が高くなります。

メルセンヌツイスタのような、非常に周期が長いと言われる乱数を利用しても数百千要素までしかソートできないなんて、ボゴソートは恐ろしいですね。

ちなみに 138要素の場合でも、1秒間に1京回計算できる計算機を使っても、完了までに 10^215年以上かかることになり、宇宙の寿命が先に尽きてしまいます。
2000要素とかどんだけかかるのかな...

#include <iostream>
#include <boost/timer.hpp>
#include <boost/random.hpp>
#include <algorithm>
#include <ctime>
#include <boost/foreach.hpp>

class Random
{
public:
  boost::mt19937 gen;
  boost::uniform_int<int> dst;
  boost::variate_generator< boost::mt19937, boost::uniform_int<int> > rand;
  Random( int N ):// call instance:
    gen( static_cast<unsigned long>(std::time(nullptr)) ), dst( 0, N ), rand( gen, dst ) {
  }
  std::ptrdiff_t operator()( std::ptrdiff_t arg ) {
    return static_cast< std::ptrdiff_t >( rand() );
  }
};

template <class InputIter, class Rand>
void bogosort(InputIter first, InputIter last, Rand rand){
  while(!std::is_sorted(first, last))
    std::random_shuffle(first, last, rand);
}

template <class InputIter>
void bogosort(InputIter first, InputIter last){
  while(!std::is_sorted(first, last))
    std::random_shuffle(first, last);
}

int main(int argc, int argv[])
{
  std::vector<int> a;

  while(true){
    int i = 0;
    std::cin >> i;
    if(std::cin.eof() || std::cin.fail()) break;
    a.push_back(i);
  }
  Random rand(a.size());
  boost::timer timer;
  bogosort(a.begin(), a.end(), rand);
  std::cout << "time = " << timer.elapsed() << std::endl;
  BOOST_FOREACH(int i, a){
    std::cout << i << ", ";
  }
  std::cout << std::endl;
  return 0;
}