Hi,

Sorry - another heapsort problem... I would really appreciate any suggestions. I am trying to sort an array of 80,000 integers and heapsort almost works except in a 3 or 4 places. I have shortened the array now to:
ra[0]=13697;
ra[1]=114269;
ra[2]=114269;
ra[3]=98514;
ra[4]=114167;
ra[5]=105169;
ra[6]=112051;
ra[7]=102445;
ra[8]=102965;
ra[9]=105585;
ra[10]=106043;
ra[11]=16995;
ra[12]=113773;
ra[13]=113863;
ra[14]=100632;
ra[15]=52743;
ra[16]=52723;

and hpsort from NRinC sorts it to
ra[0]=13697
ra[1]=1.02295e-43
ra[2]=16995
ra[3]=52723
ra[4]=52743
ra[5]=98514
ra[6]=100632
ra[7]=102445
ra[8]=102965
ra[9]=105169
ra[10]=105585
ra[11]=106043
ra[12]=112051
ra[13]=113773
ra[14]=113863
ra[15]=114167
ra[16]=114269

Is this because of the repeated 114269?

Thank you!

Hi,

hpsort.c works perfectly on the given array. Make sure you call it with the elements in a vector from nrutil, not a built-in C array. Remember that NR vectors are 1-based in C.

Saul Teukolsky

Thank you very much! You were right, I was using C-built arrays. The routine now works fine. I think the final bit that was throwing it was that I included an index array as done for the quicksort in NR. (This didn't work, but does now.)
Thanks for the pointers.

void hpsort(unsigned long n, int ra[], unsigned long indx[]) {
unsigned long i,ir,j,l,indxt;
int rra,x,y;
for(j=1;j<=n;j++){indx[j]=j;}
if (n < 2) return;
l=(n >> 1)+1;
ir=n;
for (;;) {
if (l > 1) {
y=--l;
rra=ra[y];
indxt=indx[y];
}
else {
rra=ra[ir];
indxt=indx[ir];
ra[ir]=ra[1];
indx[ir]=indx[1];
if (--ir == 1) {
ra[1]=rra;
indx[1]=indxt; break; }
}
i=l; j=l+l;
while (j <= ir) {
if (j < ir && ra[j] < ra[j+1]) j++;
if (rra < ra[j]) {
ra[i]=ra[j];
indx[i]=indx[j];
i=j;
j <<= 1; }
else j=ir+1;
}
ra[i]=rra;
indx[i]=indxt;
} }