Actual source code: sortd.c

  1: /*$Id: sortd.c,v 1.29 2001/08/07 21:29:06 bsmith Exp $*/
  2: /*
  3:    This file contains routines for sorting doubles.  Values are sorted in place.
  4:    These are provided because the general sort routines incur a great deal
  5:    of overhead in calling the comparision routines.

  7:    The word "register"  in this code is used to identify data that is not
  8:    aliased.  For some compilers, this can cause the compiler to fail to
  9:    place inner-loop variables into registers.
 10:  */
 11:  #include petsc.h
 12:  #include petscsys.h

 14: #define SWAP(a,b,t) {t=a;a=b;b=t;}
 15: 
 16: /* A simple version of quicksort; taken from Kernighan and Ritchie, page 87 */
 17: static int PetsciDqsort(PetscReal *v,int right)
 18: {
 19:   int       i,last;
 20:   PetscReal vl,tmp;
 21: 
 23:   if (right <= 1) {
 24:     if (right == 1) {
 25:       if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
 26:     }
 27:     return(0);
 28:   }
 29:   SWAP(v[0],v[right/2],tmp);
 30:   vl   = v[0];
 31:   last = 0;
 32:   for (i=1; i<=right; i++) {
 33:     if (v[i] < vl) {last++; SWAP(v[last],v[i],tmp);}
 34:   }
 35:   SWAP(v[0],v[last],tmp);
 36:   PetsciDqsort(v,last-1);
 37:   PetsciDqsort(v+last+1,right-(last+1));
 38:   return(0);
 39: }

 41: /*@
 42:    PetscSortReal - Sorts an array of doubles in place in increasing order.

 44:    Not Collective

 46:    Input Parameters:
 47: +  n  - number of values
 48: -  v  - array of doubles

 50:    Level: intermediate

 52:    Concepts: sorting^doubles

 54: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
 55: @*/
 56: int PetscSortReal(int n,PetscReal v[])
 57: {
 58:   int       j,k;
 59:   PetscReal tmp,vk;

 62:   if (n<8) {
 63:     for (k=0; k<n; k++) {
 64:       vk = v[k];
 65:       for (j=k+1; j<n; j++) {
 66:         if (vk > v[j]) {
 67:           SWAP(v[k],v[j],tmp);
 68:           vk = v[k];
 69:         }
 70:       }
 71:     }
 72:   } else {
 73:     PetsciDqsort(v,n-1);
 74:   }
 75:   return(0);
 76: }