1e5c89e4eSSatish Balay #define PETSC_DLL 2e5c89e4eSSatish Balay /* 3e5c89e4eSSatish Balay This file contains routines for sorting doubles. Values are sorted in place. 4e5c89e4eSSatish Balay These are provided because the general sort routines incur a great deal 5e5c89e4eSSatish Balay of overhead in calling the comparision routines. 6e5c89e4eSSatish Balay 7e5c89e4eSSatish Balay */ 8e5c89e4eSSatish Balay #include "petsc.h" /*I "petsc.h" I*/ 9e5c89e4eSSatish Balay #include "petscsys.h" /*I "petscsys.h" I*/ 10e5c89e4eSSatish Balay 11e5c89e4eSSatish Balay #define SWAP(a,b,t) {t=a;a=b;b=t;} 12e5c89e4eSSatish Balay 13e5c89e4eSSatish Balay #undef __FUNCT__ 14e5c89e4eSSatish Balay #define __FUNCT__ "PetscSortReal_Private" 15e5c89e4eSSatish Balay /* A simple version of quicksort; taken from Kernighan and Ritchie, page 87 */ 16e5c89e4eSSatish Balay static PetscErrorCode PetscSortReal_Private(PetscReal *v,PetscInt right) 17e5c89e4eSSatish Balay { 18e5c89e4eSSatish Balay PetscInt i,last; 19e5c89e4eSSatish Balay PetscReal vl,tmp; 20e5c89e4eSSatish Balay 21e5c89e4eSSatish Balay PetscFunctionBegin; 22e5c89e4eSSatish Balay if (right <= 1) { 23e5c89e4eSSatish Balay if (right == 1) { 24e5c89e4eSSatish Balay if (v[0] > v[1]) SWAP(v[0],v[1],tmp); 25e5c89e4eSSatish Balay } 26e5c89e4eSSatish Balay PetscFunctionReturn(0); 27e5c89e4eSSatish Balay } 28e5c89e4eSSatish Balay SWAP(v[0],v[right/2],tmp); 29e5c89e4eSSatish Balay vl = v[0]; 30e5c89e4eSSatish Balay last = 0; 31e5c89e4eSSatish Balay for (i=1; i<=right; i++) { 32e5c89e4eSSatish Balay if (v[i] < vl) {last++; SWAP(v[last],v[i],tmp);} 33e5c89e4eSSatish Balay } 34e5c89e4eSSatish Balay SWAP(v[0],v[last],tmp); 35e5c89e4eSSatish Balay PetscSortReal_Private(v,last-1); 36e5c89e4eSSatish Balay PetscSortReal_Private(v+last+1,right-(last+1)); 37e5c89e4eSSatish Balay PetscFunctionReturn(0); 38e5c89e4eSSatish Balay } 39e5c89e4eSSatish Balay 40e5c89e4eSSatish Balay #undef __FUNCT__ 41e5c89e4eSSatish Balay #define __FUNCT__ "PetscSortReal" 42e5c89e4eSSatish Balay /*@ 43e5c89e4eSSatish Balay PetscSortReal - Sorts an array of doubles in place in increasing order. 44e5c89e4eSSatish Balay 45e5c89e4eSSatish Balay Not Collective 46e5c89e4eSSatish Balay 47e5c89e4eSSatish Balay Input Parameters: 48e5c89e4eSSatish Balay + n - number of values 49e5c89e4eSSatish Balay - v - array of doubles 50e5c89e4eSSatish Balay 51e5c89e4eSSatish Balay Level: intermediate 52e5c89e4eSSatish Balay 53e5c89e4eSSatish Balay Concepts: sorting^doubles 54e5c89e4eSSatish Balay 55e5c89e4eSSatish Balay .seealso: PetscSortInt(), PetscSortRealWithPermutation() 56e5c89e4eSSatish Balay @*/ 57e5c89e4eSSatish Balay PetscErrorCode PETSC_DLLEXPORT PetscSortReal(PetscInt n,PetscReal v[]) 58e5c89e4eSSatish Balay { 59e5c89e4eSSatish Balay PetscInt j,k; 60e5c89e4eSSatish Balay PetscReal tmp,vk; 61e5c89e4eSSatish Balay 62e5c89e4eSSatish Balay PetscFunctionBegin; 63e5c89e4eSSatish Balay if (n<8) { 64e5c89e4eSSatish Balay for (k=0; k<n; k++) { 65e5c89e4eSSatish Balay vk = v[k]; 66e5c89e4eSSatish Balay for (j=k+1; j<n; j++) { 67e5c89e4eSSatish Balay if (vk > v[j]) { 68e5c89e4eSSatish Balay SWAP(v[k],v[j],tmp); 69e5c89e4eSSatish Balay vk = v[k]; 70e5c89e4eSSatish Balay } 71e5c89e4eSSatish Balay } 72e5c89e4eSSatish Balay } 73e5c89e4eSSatish Balay } else { 74e5c89e4eSSatish Balay PetscSortReal_Private(v,n-1); 75e5c89e4eSSatish Balay } 76e5c89e4eSSatish Balay PetscFunctionReturn(0); 77e5c89e4eSSatish Balay } 78e5c89e4eSSatish Balay 79*d97c5584SHong Zhang 80*d97c5584SHong Zhang 81*d97c5584SHong Zhang #undef __FUNCT__ 82*d97c5584SHong Zhang #define __FUNCT__ "PetscSortSplit" 83*d97c5584SHong Zhang /*@ 84*d97c5584SHong Zhang PetscSortSplit - Quick-sort split of an array in place. 85*d97c5584SHong Zhang Modified from SPARSEKIT2qsplit() 86*d97c5584SHong Zhang 87*d97c5584SHong Zhang Not Collective 88*d97c5584SHong Zhang 89*d97c5584SHong Zhang Input Parameters: 90*d97c5584SHong Zhang + ncut - splitig index 91*d97c5584SHong Zhang . n - number of values to sort 92*d97c5584SHong Zhang . a - array of values 93*d97c5584SHong Zhang - idx - index for array a 94*d97c5584SHong Zhang 95*d97c5584SHong Zhang Output Parameters: 96*d97c5584SHong Zhang + a - permuted array of values such that its elements satisfy: 97*d97c5584SHong Zhang abs(a[i]) >= abs(a[ncut-1]) for i < ncut and 98*d97c5584SHong Zhang abs(a[i]) <= abs(a[ncut-1]) for i >= ncut 99*d97c5584SHong Zhang - idx - permuted index of array a 100*d97c5584SHong Zhang 101*d97c5584SHong Zhang Level: intermediate 102*d97c5584SHong Zhang 103*d97c5584SHong Zhang Concepts: sorting^doubles 104*d97c5584SHong Zhang 105*d97c5584SHong Zhang .seealso: PetscSortInt(), PetscSortRealWithPermutation() 106*d97c5584SHong Zhang @*/ 107*d97c5584SHong Zhang PetscErrorCode PETSC_DLLEXPORT PetscSortSplit(PetscInt ncut,PetscInt n,PetscScalar a[],PetscInt idx[]) 108*d97c5584SHong Zhang { 109*d97c5584SHong Zhang PetscInt i,mid,last,itmp,j,first; 110*d97c5584SHong Zhang PetscScalar d,tmp; 111*d97c5584SHong Zhang PetscReal abskey; 112*d97c5584SHong Zhang 113*d97c5584SHong Zhang PetscFunctionBegin; 114*d97c5584SHong Zhang first = 0; 115*d97c5584SHong Zhang last = n-1; 116*d97c5584SHong Zhang if (ncut < first || ncut > last) PetscFunctionReturn(0); 117*d97c5584SHong Zhang 118*d97c5584SHong Zhang while (1){ 119*d97c5584SHong Zhang mid = first; 120*d97c5584SHong Zhang abskey = (d = a[mid],PetscAbsScalar(d)); 121*d97c5584SHong Zhang i = last; 122*d97c5584SHong Zhang for (j = first + 1; j <= i; ++j) { 123*d97c5584SHong Zhang if ((d = a[j],PetscAbsScalar(d)) >= abskey) { 124*d97c5584SHong Zhang ++mid; 125*d97c5584SHong Zhang /* interchange */ 126*d97c5584SHong Zhang tmp = a[mid]; itmp = idx[mid]; 127*d97c5584SHong Zhang a[mid] = a[j]; idx[mid] = idx[j]; 128*d97c5584SHong Zhang a[j] = tmp; idx[j] = itmp; 129*d97c5584SHong Zhang } 130*d97c5584SHong Zhang } 131*d97c5584SHong Zhang 132*d97c5584SHong Zhang /* interchange */ 133*d97c5584SHong Zhang tmp = a[mid]; itmp = idx[mid]; 134*d97c5584SHong Zhang a[mid] = a[first]; idx[mid] = idx[first]; 135*d97c5584SHong Zhang a[first] = tmp; idx[first] = itmp; 136*d97c5584SHong Zhang 137*d97c5584SHong Zhang /* test for while loop */ 138*d97c5584SHong Zhang if (mid == ncut) { 139*d97c5584SHong Zhang break; 140*d97c5584SHong Zhang } else if (mid > ncut){ 141*d97c5584SHong Zhang last = mid - 1; 142*d97c5584SHong Zhang } else { 143*d97c5584SHong Zhang first = mid + 1; 144*d97c5584SHong Zhang } 145*d97c5584SHong Zhang } 146*d97c5584SHong Zhang PetscFunctionReturn(0); 147*d97c5584SHong Zhang } 148