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. 6 7 */ 8 #include "petsc.h" /*I "petsc.h" I*/ 9 #include "petscsys.h" /*I "petscsys.h" I*/ 10 11 #define SWAP(a,b,t) {t=a;a=b;b=t;} 12 13 #undef __FUNCT__ 14 #define __FUNCT__ "PetscSortReal_Private" 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 21 PetscFunctionBegin; 22 if (right <= 1) { 23 if (right == 1) { 24 if (v[0] > v[1]) SWAP(v[0],v[1],tmp); 25 } 26 PetscFunctionReturn(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 PetscFunctionReturn(0); 38 } 39 40 #undef __FUNCT__ 41 #define __FUNCT__ "PetscSortReal" 42 /*@ 43 PetscSortReal - Sorts an array of doubles in place in increasing order. 44 45 Not Collective 46 47 Input Parameters: 48 + n - number of values 49 - v - array of doubles 50 51 Level: intermediate 52 53 Concepts: sorting^doubles 54 55 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 56 @*/ 57 PetscErrorCode PETSC_DLLEXPORT PetscSortReal(PetscInt n,PetscReal v[]) 58 { 59 PetscInt j,k; 60 PetscReal tmp,vk; 61 62 PetscFunctionBegin; 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 PetscFunctionReturn(0); 77 } 78 79 80 81 #undef __FUNCT__ 82 #define __FUNCT__ "PetscSortSplit" 83 /*@ 84 PetscSortSplit - Quick-sort split of an array in place. 85 Modified from SPARSEKIT2qsplit() 86 87 Not Collective 88 89 Input Parameters: 90 + ncut - splitig index 91 . n - number of values to sort 92 . a - array of values 93 - idx - index for array a 94 95 Output Parameters: 96 + a - permuted array of values such that its elements satisfy: 97 abs(a[i]) >= abs(a[ncut-1]) for i < ncut and 98 abs(a[i]) <= abs(a[ncut-1]) for i >= ncut 99 - idx - permuted index of array a 100 101 Level: intermediate 102 103 Concepts: sorting^doubles 104 105 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 106 @*/ 107 PetscErrorCode PETSC_DLLEXPORT PetscSortSplit(PetscInt ncut,PetscInt n,PetscScalar a[],PetscInt idx[]) 108 { 109 PetscInt i,mid,last,itmp,j,first; 110 PetscScalar d,tmp; 111 PetscReal abskey; 112 113 PetscFunctionBegin; 114 first = 0; 115 last = n-1; 116 if (ncut < first || ncut > last) PetscFunctionReturn(0); 117 118 while (1){ 119 mid = first; 120 abskey = (d = a[mid],PetscAbsScalar(d)); 121 i = last; 122 for (j = first + 1; j <= i; ++j) { 123 if ((d = a[j],PetscAbsScalar(d)) >= abskey) { 124 ++mid; 125 /* interchange */ 126 tmp = a[mid]; itmp = idx[mid]; 127 a[mid] = a[j]; idx[mid] = idx[j]; 128 a[j] = tmp; idx[j] = itmp; 129 } 130 } 131 132 /* interchange */ 133 tmp = a[mid]; itmp = idx[mid]; 134 a[mid] = a[first]; idx[mid] = idx[first]; 135 a[first] = tmp; idx[first] = itmp; 136 137 /* test for while loop */ 138 if (mid == ncut) { 139 break; 140 } else if (mid > ncut){ 141 last = mid - 1; 142 } else { 143 first = mid + 1; 144 } 145 } 146 PetscFunctionReturn(0); 147 } 148