Actual source code: sortd.c

  1: #define PETSC_DLL
  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:  */
 8:  #include petsc.h
 9:  #include petscsys.h

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

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

 45:    Not Collective

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

 51:    Level: intermediate

 53:    Concepts: sorting^doubles

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

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