#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <math.h>


const int K = 20;
const int N = K+10;

void randseed();
void packets(int M[N][K]);
void pr(int M[N][K]);
int po(int x, int y);
int independent(int M[N][K], int i);
int recv(int M[N][K]);
int shift(int M[N][K], int i);
double rate(int T, int suc[N+1]);
void newp(int p[K], int q);
void packets_degree(int M[N][K], int d);

main()
{
  randseed();


/*
  int M[N][K];
  packets_degree(M,3);
  pr(M);
  printf("recv=%d\n",recv(M));
*/

  int suc[N+1];
  int T=10000;

  double S;
  S=rate(T,suc);
  printf("S=%f\r\n",S);

  for(int i=K;i<=N;i++)
    printf("%2d, %f\r\n",i-K,(double)suc[i]/T);

}

double rate(int T, int suc[N+1])
{
  int M[N][K];
  for(int i=0;i<=N;i++)
    suc[i]=0;
  for(int t=0;t<T;t++)
  {
    packets_degree(M,3);
    int n=recv(M);
    for(int i=n;i<=N;i++)
      suc[i]++;
    printf("%3d, %2d\n",t,n);
  }
  double S;
  for(int i=0;i<=N;i++)
    S+=(double)suc[i]/T;
  return S;
}

int recv(int M[N][K])
{
  int r=0;
  for(int i=0;i<K;i++)
  {
    while(1)
    {
      if(r==N) return r+1;
      int ind=independent(M,i);
      r++;

      if(ind)
        break;

      shift(M,i);
    }
  }
  return r;
}

int shift(int M[N][K], int i)
{
  for(int k=i;k<N-1;k++)
    for(int j=0;j<K;j++)
      M[k][j]=M[k+1][j];
  for(int j=0;j<K;j++)
    M[N-1][j]=0;
}

int independent(int M[N][K], int i)
{
  int S = po(2,i);
  for(int s=0;s<S;s++)  //all subsets of previous rows
  {
    int ok=0;
    for(int k=0;k<K;k++)
    {
      int a=M[i][k];
      for(int j=0;j<i;j++)
        if(s>>j&1)
          a^=M[j][k];
      if(a)
      {
        ok=1;
        break;
      }
    }
    if(!ok)
      return 0;
  }
  return 1;
}

int po(int x, int y)
{
  int r=1;
  for(int i=0;i<y;i++)
    r*=x;
  return r;
}

void pr(int M[N][K])
{
  printf("M[%d][%d]=\r\n",N,K);
  for(int i=0;i<N;i++)
  {
    printf("%3d:",i);
    for(int j=0;j<K;j++)
      printf(" %s",M[i][j]?"X":"-");
    printf(" %s\r\n",independent(M,i)?"independent":"bad");
  }
}

void packets(int M[N][K])
{
  for(int i=0;i<N;i++)
    for(int j=0;j<K;j++)
      M[i][j]=rand()%2;
}

void newp(int p[K], int q)
{
  for(int i=0;i<K;i++)
    p[i]=0;

  for(int n=0;n<q;n++)
  {
    int i = rand()%(K-n);
    int j = -1;
    for(;i>=0;i--)
    {
      j++;
      for(;p[j];j++)
        ;
    }
    p[j]=1;
  }
}

void packets_degree(int M[N][K], int d)
{
  for(int i=0;i<N;i++)
    newp(M[i],d);
}

void randseed()
{
  int ti=time(0);
  printf("srand=%d\r\n",ti);
  srand(ti);
}