ANALYSIS of ALGORITHMS, Bulletin Board

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Quickselect and inversions. (III)




Greetings.

I have found a recurrence that can be used to compute the generating
functions introduced in
  http://pauillac.inria.fr/algo/AofA/mailing_list/msg00432.html.

I found this recurrence by trial and error and it remains to be proved
correct. It does yield the right generating functions for small
problem instances. The recurrence uses the fact that once quickselect
splits the input into two halves, one less and one larger than the
pivot, there are no inversions between elements from different
halves. The half that does not contain the desired element is a
genuinely random permutation. The generating function for inversions
in this case is well known (folklore). The half that does contain the
desired element represents a genuinely random subcase of the
problem. There is a special case when the desired element is the
pivot. In this case the two halves are random permutations.

The recurrence includes a scale factor that I found by trial and error
and can perhaps be simplified.

I am sending an implementation that uses the GiNaC system; the URL is
     http://www.ginac.de/.

Here is the output for n=7.

inversion gfs
1 1
2 1 + z
3 1 + 2z + 2z^2 + z^3
4 1 + 3z + 5z^2 + 6z^3 + 5z^4 + 3z^5 + z^6
5 1 + 4z + 9z^2 + 15z^3 + 20z^4 + 22z^5 + 20z^6 + 15z^7 + 9z^8 + 4z^9 + z^10
6 1 + 5z + 14z^2 + 29z^3 + 49z^4 + 71z^5 + 90z^6 + 101z^7 + 101z^8 + 90z^9 + 71z^10 + 49z^11 + 29z^12 + 14z^13 + 5z^14 + z^15
7 1 + 6z + 20z^2 + 49z^3 + 98z^4 + 169z^5 + 259z^6 + 359z^7 + 455z^8 + 531z^9 + 573z^10 + 573z^11 + 531z^12 + 455z^13 + 359z^14 + 259z^15 + 169z^16 + 98z^17 + 49z^18 + 20z^19 + 6z^20 + z^21
quickselect average gfs
0: 1
1 1
0/1 0 0
0: 2
1: 2
2 4
0/4 0 0
0: 5 + z
1: 6
2: 5 + z
3 16 + 2z
2/18 1/9 0.11111111111111111111
0: 15 + 6z + 2z^2 + z^3
1: 20 + 4z
2: 20 + 4z
3: 15 + 6z + 2z^2 + z^3
4 70 + 20z + 4z^2 + 2z^3
34/96 17/48 0.35416666666666666666
0: 52 + 33z + 15z^2 + 11z^3 + 5z^4 + 3z^5 + z^6
1: 75 + 30z + 10z^2 + 5z^3
2: 86 + 28z + 6z^2
3: 75 + 30z + 10z^2 + 5z^3
4: 52 + 33z + 15z^2 + 11z^3 + 5z^4 + 3z^5 + z^6
5 340 + 154z + 56z^2 + 32z^3 + 10z^4 + 6z^5 + 2z^6
444/600 37/50 0.74
0: 203 + 182z + 109z^2 + 81z^3 + 50z^4 + 40z^5 + 26z^6 + 15z^7 + 9z^8 + 4z^9 + z^10
1: 312 + 198z + 90z^2 + 66z^3 + 30z^4 + 18z^5 + 6z^6
2: 396 + 198z + 76z^2 + 40z^3 + 10z^4
3: 396 + 198z + 76z^2 + 40z^3 + 10z^4
4: 312 + 198z + 90z^2 + 66z^3 + 30z^4 + 18z^5 + 6z^6
5: 203 + 182z + 109z^2 + 81z^3 + 50z^4 + 40z^5 + 26z^6 + 15z^7 + 9z^8 + 4z^9 + z^10
6 1822 + 1156z + 550z^2 + 374z^3 + 180z^4 + 116z^5 + 64z^6 + 30z^7 + 18z^8 + 8z^9 + 2z^10
5508/4320 51/40 1.275
0: 877 + 1034z + 777z^2 + 631z^3 + 434z^4 + 351z^5 + 272z^6 + 206z^7 + 164z^8 + 118z^9 + 78z^10 + 49z^11 + 29z^12 + 14z^13 + 5z^14 + z^15
1: 1421 + 1274z + 763z^2 + 567z^3 + 350z^4 + 280z^5 + 182z^6 + 105z^7 + 63z^8 + 28z^9 + 7z^10
2: 1951 + 1402z + 712z^2 + 477z^3 + 255z^4 + 156z^5 + 72z^6 + 15z^7
3: 2162 + 1466z + 672z^2 + 430z^3 + 210z^4 + 80z^5 + 20z^6
4: 1951 + 1402z + 712z^2 + 477z^3 + 255z^4 + 156z^5 + 72z^6 + 15z^7
5: 1421 + 1274z + 763z^2 + 567z^3 + 350z^4 + 280z^5 + 182z^6 + 105z^7 + 63z^8 + 28z^9 + 7z^10
6: 877 + 1034z + 777z^2 + 631z^3 + 434z^4 + 351z^5 + 272z^6 + 206z^7 + 164z^8 + 118z^9 + 78z^10 + 49z^11 + 29z^12 + 14z^13 + 5z^14 + z^15
7 10660 + 8886z + 5176z^2 + 3780z^3 + 2288z^4 + 1654z^5 + 1072z^6 + 652z^7 + 454z^8 + 292z^9 + 170z^10 + 98z^11 + 58z^12 + 28z^13 + 10z^14 + 2z^15
69264/35280 481/245 1.9632653061224489796


Regards,

Marko Riedel


----------------------------------------------------------------------


#include <stdio.h>
#include <ginac/ginac.h>

using namespace std;
using namespace GiNaC;

static symbol z("z");

#define MAXN 256

static ex invgftable[MAXN];

void invgf(int k, ex *result){
  ex factor=0, term;
  int j;

  if(invgftable[k]==-1){
    for(j=0; j<k; j++){
      factor+=pow(z, j);
    }

    if(k>1){
      invgf(k-1, &term);
    }
    else{
      term=1;
    }
    term*=factor; 
    invgftable[k]=term.expand();
  }

  *result=invgftable[k];
  return;
}

static ex qselgftable[MAXN][MAXN];

void qselgfl(int n, int l, ex *result){
  ex term, prod, gf, f=factorial(n-1);
  int p;

  if(qselgftable[n][l]==-1){
    if(n==1 || n==2){
      qselgftable[n][l]=n;
    }
    else{
      term=0;
      for(p=0; p<n; p++){
	prod=1;
	if(p<l){
	  if(p>0){
	    invgf(p, &gf);
	    prod*=gf;
	  }
	  if(n-1-p>0){
	    qselgfl(n-1-p, l-p-1, &gf);
	    prod*=gf;
	  }
	}
	else if(p>l){
	  if(p>0){
	    qselgfl(p, l, &gf);
	    prod*=gf;
	  }
	  if(n-1-p>0){
	    invgf(n-1-p, &gf);
	    prod*=gf;
	  }
	}
	else{
	  if(p>0){
	    invgf(p, &gf);
	    prod*=gf;
	  }
	  if(n-1-p>0){
	    invgf(n-1-p, &gf);
	    prod*=gf;
	  }
	}
	term+=prod*f/prod.subs(z==1);
      }

      qselgftable[n][l]=term.expand();
    }
  }

  *result=qselgftable[n][l];
  return;
}

void qselgf(int n, ex *result){
  ex sum, term;
  int l;

  sum=0;
  for(l=0; l<n; l++){
    qselgfl(n, l, &term);
    sum+=term;
  }

  *result=sum;
  return;
}

void printgf(ex gf){
  int i, maxdeg=gf.degree(z);
  ex c;

  for(i=gf.ldegree(z); i<=maxdeg; i++){
    c=gf.coeff(z, i);

    if(i==0){
      cout << c;
    }
    else if(i==1){
      if(c!=1){
	cout << c << "z";
      }
      else{
	cout << "z";
      }
    }
    else{
      if(c!=1){
	cout << c << "z^" << i;
      }
      else{
	cout << "z^" << i;
      }
    }

    if(i<maxdeg){
      cout << " + ";
    }
  }
  cout << endl;

  return;
}

int main(int argc, char **argv){
  int i, j, l, n;
  ex gf, u, v, w;

  if(argc!=2){
    printf("usage: %s <n>\n", argv[0]);
    exit(1);
  }

  sscanf(argv[1], "%d", &n);
  if(n<1 || n>=MAXN){
    printf("range: 1 to %d\n", MAXN-1);
    exit(2);
  }

  for(i=0; i<MAXN; i++){
    invgftable[i]=-1;
  }

  cout << "inversion gfs" << endl;
  for(i=1; i<=n; i++){
    invgf(i, &gf);
    cout << i << " ";
    printgf(gf);
  }
  
  for(i=0; i<MAXN; i++){
    for(j=0; j<MAXN; j++){
      qselgftable[i][j]=-1;
    }
  }

  cout << "quickselect average gfs" << endl;
  for(i=1; i<=n; i++){
    for(l=0; l<i; l++){
      qselgfl(i, l, &gf);
      cout << l << ": ";
      printgf(gf);
    }
   
    qselgf(i, &gf);
    cout << i << " ";
    printgf(gf);

    u=gf.diff(z).subs(z==1);
    v=gf.subs(z==1);
    w=u/v;
    cout << u << "/" << v << " " << w << " " 
	 << evalf(u/v) << endl;
  }
}

Date Prev | Date Next | Date Index | Thread Index