Actual source code: sortd.c
1: /*
2: This file contains routines for sorting doubles. Values are sorted in place.
3: These are provided because the general sort routines incur a great deal
4: of overhead in calling the comparision routines.
6: The word "register" in this code is used to identify data that is not
7: aliased. For some compilers, this can cause the compiler to fail to
8: place inner-loop variables into registers.
9: */
10: #include petsc.h
11: #include petscsys.h
13: #define SWAP(a,b,t) {t=a;a=b;b=t;}
14:
17: /* A simple version of quicksort; taken from Kernighan and Ritchie, page 87 */
18: static PetscErrorCode PetsciDqsort(PetscReal *v,PetscInt right)
19: {
20: PetscInt i,last;
21: PetscReal vl,tmp;
22:
24: if (right <= 1) {
25: if (right == 1) {
26: if (v[0] > v[1]) SWAP(v[0],v[1],tmp);
27: }
28: return(0);
29: }
30: SWAP(v[0],v[right/2],tmp);
31: vl = v[0];
32: last = 0;
33: for (i=1; i<=right; i++) {
34: if (v[i] < vl) {last++; SWAP(v[last],v[i],tmp);}
35: }
36: SWAP(v[0],v[last],tmp);
37: PetsciDqsort(v,last-1);
38: PetsciDqsort(v+last+1,right-(last+1));
39: return(0);
40: }
44: /*@
45: PetscSortReal - Sorts an array of doubles in place in increasing order.
47: Not Collective
49: Input Parameters:
50: + n - number of values
51: - v - array of doubles
53: Level: intermediate
55: Concepts: sorting^doubles
57: .seealso: PetscSortInt(), PetscSortRealWithPermutation()
58: @*/
59: PetscErrorCode PetscSortReal(PetscInt n,PetscReal v[])
60: {
61: PetscInt j,k;
62: PetscReal tmp,vk;
65: if (n<8) {
66: for (k=0; k<n; k++) {
67: vk = v[k];
68: for (j=k+1; j<n; j++) {
69: if (vk > v[j]) {
70: SWAP(v[k],v[j],tmp);
71: vk = v[k];
72: }
73: }
74: }
75: } else {
76: PetsciDqsort(v,n-1);
77: }
78: return(0);
79: }