---

CTBsortQ.icc

Go to the documentation of this file.

/*----------------------------------------------------------------------------*/
/* C Tool Box: Designed and implemented by:                                   */
/*    Walter F.J. Mueller   Gesellschaft fuer Schwerionenforschung (GSI)      */
/*                          Planckstrasse 1, D-64291 Darmstadt, Germany       */
/*                  Email:  W.F.J.Mueller@gsi.de                              */
/*                  WWW:    http://www-kp3.gsi.de/www/kp3/people/mueller.html */
/*------------------------------------------+---------------------------------*/

#include "CTBexceptionBugcheck.hxx"
#include "CTBgrab.hxx"



template <class T, class C>
inline void CTBsortQ<T,C>::operator()(CTBvector<T>& v)
{
  const static CTBint dim_stack = 64;       // stack size
  const static CTBint max_ins   = 7;        // max. size for straight insert

  CTBgrab<T>  grab;
  C           compare;
  T           tmp;
  CTBuint32   i_ran = 4711 + (CTBuint32) &v; // initialize pseudo random number
  CTBint      i_lb  = 0;                    // current partition lower bound
  CTBint      i_ub  = v.Size()-1;           // current partition upper bound
  CTBint      i_lvl = 0;                    // stack pointer
  CTBint      i_stack[dim_stack];

  for (;;) {

    if (i_ub-i_lb < max_ins) {              // if partition small
      for (CTBint i = i_lb+1; i <= i_ub; i++) { // use straight insertion sort
        CTBint j = i-1;
        grab(tmp,v[i]);
        while (j >= i_lb && compare(tmp,v[j])) {
          grab(v[j+1],v[j]);
          j -= 1;
        }
        grab(v[j+1],tmp);
      }
      if (i_lvl == 0) return;               // stack empty -> done
      i_ub = i_stack[--i_lvl];
      i_lb = i_stack[--i_lvl]; 

    } else {                                // if partition large
      CTBuint32 i_s = i_ub - i_lb + 1;      // size of partition
      CTBint    i_p = i_lb + i_ran % i_s;   // index of pivot value
      CTBint    i   = i_lb + 1;
      CTBint    j   = i_ub;

      i_ran = i_ran * 69069 + 1;            // update random number
      grab.Swap(v[i_lb],v[i_p]);            // move pivot to first slot

      for (;;) {
        while (i < j  && compare(v[i],v[i_lb])) i++;
        while (j >= i && compare(v[i_lb],v[j])) j--;
        if (i >= j) break;
        grab.Swap(v[i],v[j]);
        j--;
        i++;
      }
      grab.Swap(v[i_lb],v[j]);              // pivot position is v[j]

      if (i_lvl >= dim_stack) 
        throw CTBexceptionBugcheck("stack overflow","CTBsortQ");

      if (i_ub-j >= j-i_lb) {               // stack largest, process smallest
        i_stack[i_lvl++] = j+1;             // stack upper
        i_stack[i_lvl++] = i_ub;
        i_ub = j-1;                         // process lower
      } else {
        i_stack[i_lvl++] = i_lb;            // stack lower
        i_stack[i_lvl++] = j-1;
        i_lb = j+1;                         // process upper
      }
    }
    
  }
  
  return;
}

// Note:
//  1. This implementation was inspired by various sources. Note that the
//     code given in Numerical Recipes in C, Chapter 8.4 has a bug !!

//------------------------------------------+-----------------------------------

template <class T, class C>
inline void CTBsortQ<T,C>::operator()(CTBvector<CTBint>& ind,
                                      const CTBvector<T>& v)
{
  const static CTBint dim_stack = 64;       // stack size
  const static CTBint max_ins   = 7;        // max. size for straight insert

  CTBgrab<CTBint> grab;
  C               compare;
  CTBint      tmp;
  CTBuint32   i_ran = 4711 + (CTBuint32) &v; // initialize pseudo random number
  CTBint      i_lb  = 0;                    // current partition lower bound
  CTBint      i_ub  = v.Size()-1;           // current partition upper bound
  CTBint      i_lvl = 0;                    // stack pointer
  CTBint      i_stack[dim_stack];
  
  ind.Resize(v.Size());
  for (CTBint i = 0; i < v.Size(); i++) ind[i] = i;
  
  for (;;) {

    if (i_ub-i_lb < max_ins) {              // if partition small
      for (CTBint i = i_lb+1; i <= i_ub; i++) { // use straight insertion sort
        CTBint j = i-1;
        tmp = ind[i];
        while (j >= i_lb && compare(v[tmp],v[ind[j]])) {
          ind[j+1] = ind[j];
          j -= 1;
        }
        ind[j+1] = tmp;
      }
      if (i_lvl == 0) return;               // stack empty -> done
      i_ub = i_stack[--i_lvl];
      i_lb = i_stack[--i_lvl]; 

    } else {                                // if partition large
      CTBuint32 i_s = i_ub - i_lb + 1;      // size of partition
      CTBint    i_p = i_lb + i_ran % i_s;   // index of pivot value
      CTBint    i   = i_lb + 1;
      CTBint    j   = i_ub;

      i_ran = i_ran * 69069 + 1;            // update random number
      grab.Swap(ind[i_lb],ind[i_p]);        // move pivot to first slot

      for (;;) {
        while (i < j  && compare(v[ind[i]],v[ind[i_lb]])) i++;
        while (j >= i && compare(v[ind[i_lb]],v[ind[j]])) j--;
        if (i >= j) break;
        grab.Swap(ind[i],ind[j]);
        j--;
        i++;
      }
      grab.Swap(ind[i_lb],ind[j]);          // pivot position is v[j]

      if (i_lvl >= dim_stack) 
        throw CTBexceptionBugcheck("stack overflow","CTBsortQ");

      if (i_ub-j >= j-i_lb) {               // stack largest, process smallest
        i_stack[i_lvl++] = j+1;             // stack upper
        i_stack[i_lvl++] = i_ub;
        i_ub = j-1;                         // process lower
      } else {
        i_stack[i_lvl++] = i_lb;            // stack lower
        i_stack[i_lvl++] = j-1;
        i_lb = j+1;                         // process upper
      }
    }
    
  }
  
  return;
}

---