酔歩の応用：種の絶滅のシミュレーション

問題設定

モデル

1. 種の数をNとし、最初はすべて種に同数の個体がいるとする。
2. 各時間ステップで、個体数は増殖率にしたがって個体数を変化させる。ただし、増殖率は平均値のまわりで毎回ランダムに変動するものとする。
3. 個体数が０になった系統はそれ以上変化しない（絶滅）。
4. 所定のステップ数に達するか、 それ以前に全系統の個体数の総和があらかじめ決められた 総個体数に達したら終了する（飽和）。つまり、他の種より先に大増殖できた種が世界を支配する。

プログラム

#include<stdio.h>
#include "rnd.h"
#define NMAX 1000000
#define N 100
main()
{ int i,j,k,nn,ninit,randinit,step;
  int a[N],sum;
  int rseed[RND_P];
  double r[N],rate;

  scanf("%d%d%d%lf%d",&nn,&ninit,&step,&rate,&randinit);

  ranset(randinit, rseed);

  for(i=0; i<nn; ++i){a[i]=ninit;}
  for(j=0; j<step; ++j){
   drnd(r,nn,rseed);
   for(i=0; i<nn; ++i){ a[i] *= (r[i]*2*rate);}
   for(i=0; i<nn; ++i)printf(" %d",a[i]);
    printf("\n");
   
   sum = 0;
   for(i=0; i<nn; ++i) sum += a[i];
   if(sum > NMAX) break;
  }
  }

• defineであらかじめ決められているのは、最大の総個体数NMAXと種の数N
• 実行時に端末から入力するのは、種数nn、最初の個体数ninit、計算するステップ数step、平均増殖率rate、それに乱数の初期値randinit
• 本質的にシミュレーション部分は

  for(i=0; i<nn; ++i){ a[i] *= (r[i]*2*rate);}

  の一行だけである。 ここでは、各種ごとに乱数の値r[i]を用いて個体数を変えている。 r[i]は0から1までの間に均等に分布する乱数（一様乱数）なので、 r[i]*2*rateがそのステップでの増殖率で、これは0からrateの2倍の値までをランダムにとることになる。その平均値がrateである。 個体数は現在の個体数にこの増殖率をかけたものとなる
• 
  

  sum = 0; for(i=0; i<nn; ++i) sum += a[i];

  はすべての種の個体数の総和を求めている。

   for(i=0; i<nn; ++i){ a[i] *= (rate + haba*(r[i]-0.5));}