src/gmxlib/gmx_lapack/dlasrt.c

   1 #include "gmx_lapack.h"
   2
   3 void
   4 F77_FUNC(dlasrt,DLASRT)(const char *id,
   5         int *n,
   6         double *d__,
   7         int *info)
   8 {
   9     int i__1, i__2;
  10
  11     int i__, j;
  12     double d1, d2, d3;
  13     int dir;
  14     double tmp;
  15     int endd;
  16     int stack[64];
  17     double dmnmx;
  18     int start;
  19     int stkpnt;
  20
  21     --d__;
  22
  23     *info = 0;
  24     dir = -1;
  25     if (*id=='D' || *id=='d')
  26         dir = 0;
  27     else if (*id=='I' || *id=='i')
  28         dir = 1;
  29
  30     if (dir == -1) {
  31         *info = -1;
  32     } else if (*n < 0) {
  33         *info = -2;
  34     }
  35     if (*info != 0) {
  36         i__1 = -(*info);
  37         return;
  38     }
  39     if (*n <= 1) {
  40         return;
  41     }
  42
  43     stkpnt = 1;
  44     stack[0] = 1;
  45     stack[1] = *n;
  46 L10:
  47     start = stack[(stkpnt << 1) - 2];
  48     endd = stack[(stkpnt << 1) - 1];
  49     --stkpnt;
  50     if (endd - start <= 20 && endd - start > 0) {
  51
  52
  53         if (dir == 0) {
  54
  55             i__1 = endd;
  56             for (i__ = start + 1; i__ <= i__1; ++i__) {
  57                 i__2 = start + 1;
  58                 for (j = i__; j >= i__2; --j) {
  59                     if (d__[j] > d__[j - 1]) {
  60                         dmnmx = d__[j];
  61                         d__[j] = d__[j - 1];
  62                         d__[j - 1] = dmnmx;
  63                     } else {
  64                         goto L30;
  65                     }
  66                 }
  67 L30:
  68                 ;
  69             }
  70
  71         } else {
  72
  73             i__1 = endd;
  74             for (i__ = start + 1; i__ <= i__1; ++i__) {
  75                 i__2 = start + 1;
  76                 for (j = i__; j >= i__2; --j) {
  77                     if (d__[j] < d__[j - 1]) {
  78                         dmnmx = d__[j];
  79                         d__[j] = d__[j - 1];
  80                         d__[j - 1] = dmnmx;
  81                     } else {
  82                         goto L50;
  83                     }
  84                 }
  85 L50:
  86                 ;
  87             }
  88
  89         }
  90
  91     } else if (endd - start > 20) {
  92
  93         d1 = d__[start];
  94         d2 = d__[endd];
  95         i__ = (start + endd) / 2;
  96         d3 = d__[i__];
  97         if (d1 < d2) {
  98             if (d3 < d1) {
  99                 dmnmx = d1;
 100             } else if (d3 < d2) {
 101                 dmnmx = d3;
 102             } else {
 103                 dmnmx = d2;
 104             }
 105         } else {
 106             if (d3 < d2) {
 107                 dmnmx = d2;
 108             } else if (d3 < d1) {
 109                 dmnmx = d3;
 110             } else {
 111                 dmnmx = d1;
 112             }
 113         }
 114
 115         if (dir == 0) {
 116
 117             i__ = start - 1;
 118             j = endd + 1;
 119 L60:
 120 L70:
 121             --j;
 122             if (d__[j] < dmnmx) {
 123                 goto L70;
 124             }
 125 L80:
 126             ++i__;
 127             if (d__[i__] > dmnmx) {
 128                 goto L80;
 129             }
 130             if (i__ < j) {
 131                 tmp = d__[i__];
 132                 d__[i__] = d__[j];
 133                 d__[j] = tmp;
 134                 goto L60;
 135             }
 136             if (j - start > endd - j - 1) {
 137                 ++stkpnt;
 138                 stack[(stkpnt << 1) - 2] = start;
 139                 stack[(stkpnt << 1) - 1] = j;
 140                 ++stkpnt;
 141                 stack[(stkpnt << 1) - 2] = j + 1;
 142                 stack[(stkpnt << 1) - 1] = endd;
 143             } else {
 144                 ++stkpnt;
 145                 stack[(stkpnt << 1) - 2] = j + 1;
 146                 stack[(stkpnt << 1) - 1] = endd;
 147                 ++stkpnt;
 148                 stack[(stkpnt << 1) - 2] = start;
 149                 stack[(stkpnt << 1) - 1] = j;
 150             }
 151         } else {
 152
 153             i__ = start - 1;
 154             j = endd + 1;
 155 L90:
 156 L100:
 157             --j;
 158             if (d__[j] > dmnmx) {
 159                 goto L100;
 160             }
 161 L110:
 162             ++i__;
 163             if (d__[i__] < dmnmx) {
 164                 goto L110;
 165             }
 166             if (i__ < j) {
 167                 tmp = d__[i__];
 168                 d__[i__] = d__[j];
 169                 d__[j] = tmp;
 170                 goto L90;
 171             }
 172             if (j - start > endd - j - 1) {
 173                 ++stkpnt;
 174                 stack[(stkpnt << 1) - 2] = start;
 175                 stack[(stkpnt << 1) - 1] = j;
 176                 ++stkpnt;
 177                 stack[(stkpnt << 1) - 2] = j + 1;
 178                 stack[(stkpnt << 1) - 1] = endd;
 179             } else {
 180                 ++stkpnt;
 181                 stack[(stkpnt << 1) - 2] = j + 1;
 182                 stack[(stkpnt << 1) - 1] = endd;
 183                 ++stkpnt;
 184                 stack[(stkpnt << 1) - 2] = start;
 185                 stack[(stkpnt << 1) - 1] = j;
 186             }
 187         }
 188     }
 189     if (stkpnt > 0) {
 190         goto L10;
 191     }
 192     return;
 193
 194 }