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 #undef __FUNCT__ 80 #define __FUNCT__ "PetscSortSplit" 81 /*@ 82 PetscSortSplit - Quick-sort split of an array of PetscScalars in place. 83 84 Not Collective 85 86 Input Parameters: 87 + ncut - splitig index 88 . n - number of values to sort 89 . a - array of values 90 - idx - index for array a 91 92 Output Parameters: 93 + a - permuted array of values such that its elements satisfy: 94 abs(a[i]) >= abs(a[ncut-1]) for i < ncut and 95 abs(a[i]) <= abs(a[ncut-1]) for i >= ncut 96 - idx - permuted index of array a 97 98 Level: intermediate 99 100 Concepts: sorting^doubles 101 102 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 103 @*/ 104 PetscErrorCode PETSC_DLLEXPORT PetscSortSplit(PetscInt ncut,PetscInt n,PetscScalar a[],PetscInt idx[]) 105 { 106 PetscInt i,mid,last,itmp,j,first; 107 PetscScalar d,tmp; 108 PetscReal abskey; 109 110 PetscFunctionBegin; 111 first = 0; 112 last = n-1; 113 if (ncut < first || ncut > last) PetscFunctionReturn(0); 114 115 while (1){ 116 mid = first; 117 abskey = (d = a[mid],PetscAbsScalar(d)); 118 i = last; 119 for (j = first + 1; j <= i; ++j) { 120 if ((d = a[j],PetscAbsScalar(d)) >= abskey) { 121 ++mid; 122 /* interchange */ 123 tmp = a[mid]; itmp = idx[mid]; 124 a[mid] = a[j]; idx[mid] = idx[j]; 125 a[j] = tmp; idx[j] = itmp; 126 } 127 } 128 129 /* interchange */ 130 tmp = a[mid]; itmp = idx[mid]; 131 a[mid] = a[first]; idx[mid] = idx[first]; 132 a[first] = tmp; idx[first] = itmp; 133 134 /* test for while loop */ 135 if (mid == ncut) { 136 break; 137 } else if (mid > ncut){ 138 last = mid - 1; 139 } else { 140 first = mid + 1; 141 } 142 } 143 PetscFunctionReturn(0); 144 } 145 146 #undef __FUNCT__ 147 #define __FUNCT__ "PetscSortSplitReal" 148 /*@ 149 PetscSortSplitReal - Quick-sort split of an array of PetscReals in place. 150 151 Not Collective 152 153 Input Parameters: 154 + ncut - splitig index 155 . n - number of values to sort 156 . a - array of values in PetscReal 157 - idx - index for array a 158 159 Output Parameters: 160 + a - permuted array of real values such that its elements satisfy: 161 abs(a[i]) >= abs(a[ncut-1]) for i < ncut and 162 abs(a[i]) <= abs(a[ncut-1]) for i >= ncut 163 - idx - permuted index of array a 164 165 Level: intermediate 166 167 Concepts: sorting^doubles 168 169 .seealso: PetscSortInt(), PetscSortRealWithPermutation() 170 @*/ 171 PetscErrorCode PETSC_DLLEXPORT PetscSortSplitReal(PetscInt ncut,PetscInt n,PetscReal a[],PetscInt idx[]) 172 { 173 PetscInt i,mid,last,itmp,j,first; 174 PetscReal d,tmp; 175 PetscReal abskey; 176 177 PetscFunctionBegin; 178 first = 0; 179 last = n-1; 180 if (ncut < first || ncut > last) PetscFunctionReturn(0); 181 182 while (1){ 183 mid = first; 184 abskey = (d = a[mid],PetscAbsScalar(d)); 185 i = last; 186 for (j = first + 1; j <= i; ++j) { 187 if ((d = a[j],PetscAbsScalar(d)) >= abskey) { 188 ++mid; 189 /* interchange */ 190 tmp = a[mid]; itmp = idx[mid]; 191 a[mid] = a[j]; idx[mid] = idx[j]; 192 a[j] = tmp; idx[j] = itmp; 193 } 194 } 195 196 /* interchange */ 197 tmp = a[mid]; itmp = idx[mid]; 198 a[mid] = a[first]; idx[mid] = idx[first]; 199 a[first] = tmp; idx[first] = itmp; 200 201 /* test for while loop */ 202 if (mid == ncut) { 203 break; 204 } else if (mid > ncut){ 205 last = mid - 1; 206 } else { 207 first = mid + 1; 208 } 209 } 210 PetscFunctionReturn(0); 211 } 212 213