1c-----------------------------------------------------------------------
2c\BeginDoc
3c
4c\Name: dsesrt
5c
6c\Description:
7c  Sort the array X in the order specified by WHICH and optionally
8c  apply the permutation to the columns of the matrix A.
9c
10c\Usage:
11c  call dsesrt
12c     ( WHICH, APPLY, N, X, NA, A, LDA)
13c
14c\Arguments
15c  WHICH   Character*2.  (Input)
16c          'LM' -> X is sorted into increasing order of magnitude.
17c          'SM' -> X is sorted into decreasing order of magnitude.
18c          'LA' -> X is sorted into increasing order of algebraic.
19c          'SA' -> X is sorted into decreasing order of algebraic.
20c
21c  APPLY   Logical.  (Input)
22c          APPLY = .TRUE.  -> apply the sorted order to A.
23c          APPLY = .FALSE. -> do not apply the sorted order to A.
24c
25c  N       Integer.  (INPUT)
26c          Dimension of the array X.
27c
28c  X      Double precision array of length N.  (INPUT/OUTPUT)
29c          The array to be sorted.
30c
31c  NA      Integer.  (INPUT)
32c          Number of rows of the matrix A.
33c
34c  A      Double precision array of length NA by N.  (INPUT/OUTPUT)
35c
36c  LDA     Integer.  (INPUT)
37c          Leading dimension of A.
38c
39c\EndDoc
40c
41c-----------------------------------------------------------------------
42c
43c\BeginLib
44c
45c\Routines
46c     dswap  Level 1 BLAS that swaps the contents of two vectors.
47c
48c\Authors
49c     Danny Sorensen               Phuong Vu
50c     Richard Lehoucq              CRPC / Rice University
51c     Dept. of Computational &     Houston, Texas
52c     Applied Mathematics
53c     Rice University
54c     Houston, Texas
55c
56c\Revision history:
57c     12/15/93: Version ' 2.1'.
58c               Adapted from the sort routine in LANSO and
59c               the ARPACK code dsortr
60c
61c\SCCS Information: @(#)
62c FILE: sesrt.F   SID: 2.3   DATE OF SID: 4/19/96   RELEASE: 2
63c
64c\EndLib
65c
66c-----------------------------------------------------------------------
67c
68      subroutine dsesrt (which, apply, n, x, na, a, lda)
69c
70c     %------------------%
71c     | Scalar Arguments |
72c     %------------------%
73c
74      character*2 which
75      logical    apply
76      integer    lda, n, na
77c
78c     %-----------------%
79c     | Array Arguments |
80c     %-----------------%
81c
82      Double precision
83     &           x(0:n-1), a(lda, 0:n-1)
84c
85c     %---------------%
86c     | Local Scalars |
87c     %---------------%
88c
89      integer    i, igap, j
90      Double precision
91     &           temp
92c
93c     %----------------------%
94c     | External Subroutines |
95c     %----------------------%
96c
97      external   dswap
98c
99c     %-----------------------%
100c     | Executable Statements |
101c     %-----------------------%
102c
103      igap = n / 2
104c
105      if (which .eq. 'SA') then
106c
107c        X is sorted into decreasing order of algebraic.
108c
109   10    continue
110         if (igap .eq. 0) go to 9000
111         do 30 i = igap, n-1
112            j = i-igap
113   20       continue
114c
115            if (j.lt.0) go to 30
116c
117            if (x(j).lt.x(j+igap)) then
118               temp = x(j)
119               x(j) = x(j+igap)
120               x(j+igap) = temp
121               if (apply) call dswap( na, a(1, j), 1, a(1,j+igap), 1)
122            else
123               go to 30
124            endif
125            j = j-igap
126            go to 20
127   30    continue
128         igap = igap / 2
129         go to 10
130c
131      else if (which .eq. 'SM') then
132c
133c        X is sorted into decreasing order of magnitude.
134c
135   40    continue
136         if (igap .eq. 0) go to 9000
137         do 60 i = igap, n-1
138            j = i-igap
139   50       continue
140c
141            if (j.lt.0) go to 60
142c
143            if (abs(x(j)).lt.abs(x(j+igap))) then
144               temp = x(j)
145               x(j) = x(j+igap)
146               x(j+igap) = temp
147               if (apply) call dswap( na, a(1, j), 1, a(1,j+igap), 1)
148            else
149               go to 60
150            endif
151            j = j-igap
152            go to 50
153   60    continue
154         igap = igap / 2
155         go to 40
156c
157      else if (which .eq. 'LA') then
158c
159c        X is sorted into increasing order of algebraic.
160c
161   70    continue
162         if (igap .eq. 0) go to 9000
163         do 90 i = igap, n-1
164            j = i-igap
165   80       continue
166c
167            if (j.lt.0) go to 90
168c
169            if (x(j).gt.x(j+igap)) then
170               temp = x(j)
171               x(j) = x(j+igap)
172               x(j+igap) = temp
173               if (apply) call dswap( na, a(1, j), 1, a(1,j+igap), 1)
174            else
175               go to 90
176            endif
177            j = j-igap
178            go to 80
179   90    continue
180         igap = igap / 2
181         go to 70
182c
183      else if (which .eq. 'LM') then
184c
185c        X is sorted into increasing order of magnitude.
186c
187  100    continue
188         if (igap .eq. 0) go to 9000
189         do 120 i = igap, n-1
190            j = i-igap
191  110       continue
192c
193            if (j.lt.0) go to 120
194c
195            if (abs(x(j)).gt.abs(x(j+igap))) then
196               temp = x(j)
197               x(j) = x(j+igap)
198               x(j+igap) = temp
199               if (apply) call dswap( na, a(1, j), 1, a(1,j+igap), 1)
200            else
201               go to 120
202            endif
203            j = j-igap
204            go to 110
205  120    continue
206         igap = igap / 2
207         go to 100
208      end if
209c
210 9000 continue
211      return
212c
213c     %---------------%
214c     | End of dsesrt |
215c     %---------------%
216c
217      end
218