1 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 <petscsys.h> /*I "petscsys.h" I*/ 9 10 #define SWAP(a,b,t) {t=a;a=b;b=t;} 11 12 #undef __FUNCT__ 13 #define __FUNCT__ "PetscSortReal_Private" 14 /* A simple version of quicksort; taken from Kernighan and Ritchie, page 87 */ 15 static PetscErrorCode PetscSortReal_Private(PetscReal *v,PetscInt right) 16 { 17 PetscInt i,last; 18 PetscReal vl,tmp; 19 20 PetscFunctionBegin; 21 if (right <= 1) { 22 if (right == 1) { 23 if (v[0] > v[1]) SWAP(v[0],v[1],tmp); 24 } 25 PetscFunctionReturn(0); 26 } 27 SWAP(v[0],v[right/2],tmp); 28 vl = v[0]; 29 last = 0; 30 for (i=1; i<=right; i++) { 31 if (v[i] < vl) {last++; SWAP(v[last],v[i],tmp);} 32 } 33 SWAP(v[0],v[last],tmp); 34 PetscSortReal_Private(v,last-1); 35 PetscSortReal_Private(v+last+1,right-(last+1)); 36 PetscFunctionReturn(0); 37 } 38 39 #undef __FUNCT__ 40 #define __FUNCT__ "PetscSortReal" 41 /*@ 42 PetscSortReal - Sorts an array of doubles in place in increasing order. 43 44 Not Collective 45 46 Input Parameters: 47 + n - number of values 48 - v - array of doubles 49 50 Level: intermediate 51 52 Concepts: sorting^doubles 53 54 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 55 @*/ 56 PetscErrorCode PetscSortReal(PetscInt n,PetscReal v[]) 57 { 58 PetscInt j,k; 59 PetscReal tmp,vk; 60 61 PetscFunctionBegin; 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 PetscSortReal_Private(v,n-1); 73 PetscFunctionReturn(0); 74 } 75 76 #undef __FUNCT__ 77 #define __FUNCT__ "PetscSortSplit" 78 /*@ 79 PetscSortSplit - Quick-sort split of an array of PetscScalars in place. 80 81 Not Collective 82 83 Input Parameters: 84 + ncut - splitig index 85 . n - number of values to sort 86 . a - array of values 87 - idx - index for array a 88 89 Output Parameters: 90 + a - permuted array of values such that its elements satisfy: 91 abs(a[i]) >= abs(a[ncut-1]) for i < ncut and 92 abs(a[i]) <= abs(a[ncut-1]) for i >= ncut 93 - idx - permuted index of array a 94 95 Level: intermediate 96 97 Concepts: sorting^doubles 98 99 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 100 @*/ 101 PetscErrorCode PetscSortSplit(PetscInt ncut,PetscInt n,PetscScalar a[],PetscInt idx[]) 102 { 103 PetscInt i,mid,last,itmp,j,first; 104 PetscScalar d,tmp; 105 PetscReal abskey; 106 107 PetscFunctionBegin; 108 first = 0; 109 last = n-1; 110 if (ncut < first || ncut > last) PetscFunctionReturn(0); 111 112 while (1) { 113 mid = first; 114 abskey = (d = a[mid],PetscAbsScalar(d)); 115 i = last; 116 for (j = first + 1; j <= i; ++j) { 117 if ((d = a[j],PetscAbsScalar(d)) >= abskey) { 118 ++mid; 119 /* interchange */ 120 tmp = a[mid]; itmp = idx[mid]; 121 a[mid] = a[j]; idx[mid] = idx[j]; 122 a[j] = tmp; idx[j] = itmp; 123 } 124 } 125 126 /* interchange */ 127 tmp = a[mid]; itmp = idx[mid]; 128 a[mid] = a[first]; idx[mid] = idx[first]; 129 a[first] = tmp; idx[first] = itmp; 130 131 /* test for while loop */ 132 if (mid == ncut) break; 133 else if (mid > ncut) last = mid - 1; 134 else first = mid + 1; 135 } 136 PetscFunctionReturn(0); 137 } 138 139 #undef __FUNCT__ 140 #define __FUNCT__ "PetscSortSplitReal" 141 /*@ 142 PetscSortSplitReal - Quick-sort split of an array of PetscReals in place. 143 144 Not Collective 145 146 Input Parameters: 147 + ncut - splitig index 148 . n - number of values to sort 149 . a - array of values in PetscReal 150 - idx - index for array a 151 152 Output Parameters: 153 + a - permuted array of real values such that its elements satisfy: 154 abs(a[i]) >= abs(a[ncut-1]) for i < ncut and 155 abs(a[i]) <= abs(a[ncut-1]) for i >= ncut 156 - idx - permuted index of array a 157 158 Level: intermediate 159 160 Concepts: sorting^doubles 161 162 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 163 @*/ 164 PetscErrorCode PetscSortSplitReal(PetscInt ncut,PetscInt n,PetscReal a[],PetscInt idx[]) 165 { 166 PetscInt i,mid,last,itmp,j,first; 167 PetscReal d,tmp; 168 PetscReal abskey; 169 170 PetscFunctionBegin; 171 first = 0; 172 last = n-1; 173 if (ncut < first || ncut > last) PetscFunctionReturn(0); 174 175 while (1) { 176 mid = first; 177 abskey = (d = a[mid],PetscAbsReal(d)); 178 i = last; 179 for (j = first + 1; j <= i; ++j) { 180 if ((d = a[j],PetscAbsReal(d)) >= abskey) { 181 ++mid; 182 /* interchange */ 183 tmp = a[mid]; itmp = idx[mid]; 184 a[mid] = a[j]; idx[mid] = idx[j]; 185 a[j] = tmp; idx[j] = itmp; 186 } 187 } 188 189 /* interchange */ 190 tmp = a[mid]; itmp = idx[mid]; 191 a[mid] = a[first]; idx[mid] = idx[first]; 192 a[first] = tmp; idx[first] = itmp; 193 194 /* test for while loop */ 195 if (mid == ncut) break; 196 else if (mid > ncut) last = mid - 1; 197 else first = mid + 1; 198 } 199 PetscFunctionReturn(0); 200 } 201 202