← Back to branch summary

~ubuntu-branches/ubuntu/saucy/python-scipy/saucy

« back to all changes in this revision

Viewing changes to Lib/sandbox/arpack/ARPACK/SRC/dsortc.f

  • Committer: Bazaar Package Importer
  • Author(s): Ondrej Certik
  • Date: 2008-06-16 22:58:01 UTC
  • mfrom: (2.1.24 intrepid)
  • Revision ID: james.westby@ubuntu.com-20080616225801-irdhrpcwiocfbcmt
Tags: 0.6.0-12
http://bugs.debian.org/489149
* The description updated to match the current SciPy (Closes: #489149).
* Standards-Version bumped to 3.8.0 (no action needed)
* Build-Depends: netcdf-dev changed to libnetcdf-dev

Show diffs side-by-side

added added

removed removed

Lines of Context:
Lib/sandbox/arpack/ARPACK/SRC/dsortc.f
1
 
c-----------------------------------------------------------------------
2
 
c\BeginDoc
3
 
c
4
 
c\Name: dsortc
5
 
c
6
 
c\Description:
7
 
c  Sorts the complex array in XREAL and XIMAG into the order 
8
 
c  specified by WHICH and optionally applies the permutation to the
9
 
c  real array Y. It is assumed that if an element of XIMAG is
10
 
c  nonzero, then its negative is also an element. In other words,
11
 
c  both members of a complex conjugate pair are to be sorted and the
12
 
c  pairs are kept adjacent to each other.
13
 
c
14
 
c\Usage:
15
 
c  call dsortc
16
 
c     ( WHICH, APPLY, N, XREAL, XIMAG, Y )
17
 
c
18
 
c\Arguments
19
 
c  WHICH   Character*2.  (Input)
20
 
c          'LM' -> sort XREAL,XIMAG into increasing order of magnitude.
21
 
c          'SM' -> sort XREAL,XIMAG into decreasing order of magnitude.
22
 
c          'LR' -> sort XREAL into increasing order of algebraic.
23
 
c          'SR' -> sort XREAL into decreasing order of algebraic.
24
 
c          'LI' -> sort XIMAG into increasing order of magnitude.
25
 
c          'SI' -> sort XIMAG into decreasing order of magnitude.
26
 
c          NOTE: If an element of XIMAG is non-zero, then its negative
27
 
c                is also an element.
28
 
c
29
 
c  APPLY   Logical.  (Input)
30
 
c          APPLY = .TRUE.  -> apply the sorted order to array Y.
31
 
c          APPLY = .FALSE. -> do not apply the sorted order to array Y.
32
 
c
33
 
c  N       Integer.  (INPUT)
34
 
c          Size of the arrays.
35
 
c
36
 
c  XREAL,  Double precision array of length N.  (INPUT/OUTPUT)
37
 
c  XIMAG   Real and imaginary part of the array to be sorted.
38
 
c
39
 
c  Y       Double precision array of length N.  (INPUT/OUTPUT)
40
 
c
41
 
c\EndDoc
42
 
c
43
 
c-----------------------------------------------------------------------
44
 
c
45
 
c\BeginLib
46
 
c
47
 
c\Author
48
 
c     Danny Sorensen               Phuong Vu
49
 
c     Richard Lehoucq              CRPC / Rice University
50
 
c     Dept. of Computational &     Houston, Texas
51
 
c     Applied Mathematics
52
 
c     Rice University           
53
 
c     Houston, Texas            
54
 
c
55
 
c\Revision history:
56
 
c     xx/xx/92: Version ' 2.1'
57
 
c               Adapted from the sort routine in LANSO.
58
 
c
59
 
c\SCCS Information: @(#) 
60
 
c FILE: sortc.F   SID: 2.3   DATE OF SID: 4/20/96   RELEASE: 2
61
 
c
62
 
c\EndLib
63
 
c
64
 
c-----------------------------------------------------------------------
65
 
c
66
 
      subroutine dsortc (which, apply, n, xreal, ximag, y)
67
 
c
68
 
c     %------------------%
69
 
c     | Scalar Arguments |
70
 
c     %------------------%
71
 
c
72
 
      character*2 which
73
 
      logical    apply
74
 
      integer    n
75
 
c
76
 
c     %-----------------%
77
 
c     | Array Arguments |
78
 
c     %-----------------%
79
 
c
80
 
      Double precision     
81
 
     &           xreal(0:n-1), ximag(0:n-1), y(0:n-1)
82
 
c
83
 
c     %---------------%
84
 
c     | Local Scalars |
85
 
c     %---------------%
86
 
c
87
 
      integer    i, igap, j
88
 
      Double precision     
89
 
     &           temp, temp1, temp2
90
 
c
91
 
c     %--------------------%
92
 
c     | External Functions |
93
 
c     %--------------------%
94
 
c
95
 
      Double precision     
96
 
     &           dlapy2
97
 
      external   dlapy2
98
 
c
99
 
c     %-----------------------%
100
 
c     | Executable Statements |
101
 
c     %-----------------------%
102
 
c
103
 
      igap = n / 2
104
 
c 
105
 
      if (which .eq. 'LM') then
106
 
c
107
 
c        %------------------------------------------------------%
108
 
c        | Sort XREAL,XIMAG into increasing order of magnitude. |
109
 
c        %------------------------------------------------------%
110
 
c
111
 
   10    continue
112
 
         if (igap .eq. 0) go to 9000
113
 
c
114
 
         do 30 i = igap, n-1
115
 
            j = i-igap
116
 
   20       continue
117
 
c
118
 
            if (j.lt.0) go to 30
119
 
c
120
 
            temp1 = dlapy2(xreal(j),ximag(j))
121
 
            temp2 = dlapy2(xreal(j+igap),ximag(j+igap))
122
 
c
123
 
            if (temp1.gt.temp2) then
124
 
                temp = xreal(j)
125
 
                xreal(j) = xreal(j+igap)
126
 
                xreal(j+igap) = temp
127
 
c
128
 
                temp = ximag(j)
129
 
                ximag(j) = ximag(j+igap)
130
 
                ximag(j+igap) = temp
131
 
c
132
 
                if (apply) then
133
 
                    temp = y(j)
134
 
                    y(j) = y(j+igap)
135
 
                    y(j+igap) = temp
136
 
                end if
137
 
            else
138
 
                go to 30
139
 
            end if
140
 
            j = j-igap
141
 
            go to 20
142
 
   30    continue
143
 
         igap = igap / 2
144
 
         go to 10
145
 
c
146
 
      else if (which .eq. 'SM') then
147
 
c
148
 
c        %------------------------------------------------------%
149
 
c        | Sort XREAL,XIMAG into decreasing order of magnitude. |
150
 
c        %------------------------------------------------------%
151
 
c
152
 
   40    continue
153
 
         if (igap .eq. 0) go to 9000
154
 
c
155
 
         do 60 i = igap, n-1
156
 
            j = i-igap
157
 
   50       continue
158
 
c
159
 
            if (j .lt. 0) go to 60
160
 
c
161
 
            temp1 = dlapy2(xreal(j),ximag(j))
162
 
            temp2 = dlapy2(xreal(j+igap),ximag(j+igap))
163
 
c
164
 
            if (temp1.lt.temp2) then
165
 
               temp = xreal(j)
166
 
               xreal(j) = xreal(j+igap)
167
 
               xreal(j+igap) = temp
168
 
c
169
 
               temp = ximag(j)
170
 
               ximag(j) = ximag(j+igap)
171
 
               ximag(j+igap) = temp
172
 
c 
173
 
               if (apply) then
174
 
                  temp = y(j)
175
 
                  y(j) = y(j+igap)
176
 
                  y(j+igap) = temp
177
 
               end if
178
 
            else
179
 
               go to 60
180
 
            endif
181
 
            j = j-igap
182
 
            go to 50
183
 
   60    continue
184
 
         igap = igap / 2
185
 
         go to 40
186
 
c 
187
 
      else if (which .eq. 'LR') then
188
 
c
189
 
c        %------------------------------------------------%
190
 
c        | Sort XREAL into increasing order of algebraic. |
191
 
c        %------------------------------------------------%
192
 
c
193
 
   70    continue
194
 
         if (igap .eq. 0) go to 9000
195
 
c
196
 
         do 90 i = igap, n-1
197
 
            j = i-igap
198
 
   80       continue
199
 
c
200
 
            if (j.lt.0) go to 90
201
 
c
202
 
            if (xreal(j).gt.xreal(j+igap)) then
203
 
               temp = xreal(j)
204
 
               xreal(j) = xreal(j+igap)
205
 
               xreal(j+igap) = temp
206
 
c
207
 
               temp = ximag(j)
208
 
               ximag(j) = ximag(j+igap)
209
 
               ximag(j+igap) = temp
210
 
c 
211
 
               if (apply) then
212
 
                  temp = y(j)
213
 
                  y(j) = y(j+igap)
214
 
                  y(j+igap) = temp
215
 
               end if
216
 
            else
217
 
               go to 90
218
 
            endif
219
 
            j = j-igap
220
 
            go to 80
221
 
   90    continue
222
 
         igap = igap / 2
223
 
         go to 70
224
 
c 
225
 
      else if (which .eq. 'SR') then
226
 
c
227
 
c        %------------------------------------------------%
228
 
c        | Sort XREAL into decreasing order of algebraic. |
229
 
c        %------------------------------------------------%
230
 
c
231
 
  100    continue
232
 
         if (igap .eq. 0) go to 9000
233
 
         do 120 i = igap, n-1
234
 
            j = i-igap
235
 
  110       continue
236
 
c
237
 
            if (j.lt.0) go to 120
238
 
c
239
 
            if (xreal(j).lt.xreal(j+igap)) then
240
 
               temp = xreal(j)
241
 
               xreal(j) = xreal(j+igap)
242
 
               xreal(j+igap) = temp
243
 
c
244
 
               temp = ximag(j)
245
 
               ximag(j) = ximag(j+igap)
246
 
               ximag(j+igap) = temp
247
 
c 
248
 
               if (apply) then
249
 
                  temp = y(j)
250
 
                  y(j) = y(j+igap)
251
 
                  y(j+igap) = temp
252
 
               end if
253
 
            else
254
 
               go to 120
255
 
            endif
256
 
            j = j-igap
257
 
            go to 110
258
 
  120    continue
259
 
         igap = igap / 2
260
 
         go to 100
261
 
c 
262
 
      else if (which .eq. 'LI') then
263
 
c
264
 
c        %------------------------------------------------%
265
 
c        | Sort XIMAG into increasing order of magnitude. |
266
 
c        %------------------------------------------------%
267
 
c
268
 
  130    continue
269
 
         if (igap .eq. 0) go to 9000
270
 
         do 150 i = igap, n-1
271
 
            j = i-igap
272
 
  140       continue
273
 
c
274
 
            if (j.lt.0) go to 150
275
 
c
276
 
            if (abs(ximag(j)).gt.abs(ximag(j+igap))) then
277
 
               temp = xreal(j)
278
 
               xreal(j) = xreal(j+igap)
279
 
               xreal(j+igap) = temp
280
 
c
281
 
               temp = ximag(j)
282
 
               ximag(j) = ximag(j+igap)
283
 
               ximag(j+igap) = temp
284
 
c 
285
 
               if (apply) then
286
 
                  temp = y(j)
287
 
                  y(j) = y(j+igap)
288
 
                  y(j+igap) = temp
289
 
               end if
290
 
            else
291
 
               go to 150
292
 
            endif
293
 
            j = j-igap
294
 
            go to 140
295
 
  150    continue
296
 
         igap = igap / 2
297
 
         go to 130
298
 
c 
299
 
      else if (which .eq. 'SI') then
300
 
c
301
 
c        %------------------------------------------------%
302
 
c        | Sort XIMAG into decreasing order of magnitude. |
303
 
c        %------------------------------------------------%
304
 
c
305
 
  160    continue
306
 
         if (igap .eq. 0) go to 9000
307
 
         do 180 i = igap, n-1
308
 
            j = i-igap
309
 
  170       continue
310
 
c
311
 
            if (j.lt.0) go to 180
312
 
c
313
 
            if (abs(ximag(j)).lt.abs(ximag(j+igap))) then
314
 
               temp = xreal(j)
315
 
               xreal(j) = xreal(j+igap)
316
 
               xreal(j+igap) = temp
317
 
c
318
 
               temp = ximag(j)
319
 
               ximag(j) = ximag(j+igap)
320
 
               ximag(j+igap) = temp
321
 
c 
322
 
               if (apply) then
323
 
                  temp = y(j)
324
 
                  y(j) = y(j+igap)
325
 
                  y(j+igap) = temp
326
 
               end if
327
 
            else
328
 
               go to 180
329
 
            endif
330
 
            j = j-igap
331
 
            go to 170
332
 
  180    continue
333
 
         igap = igap / 2
334
 
         go to 160
335
 
      end if
336
 
c 
337
 
 9000 continue
338
 
      return
339
 
c
340
 
c     %---------------%
341
 
c     | End of dsortc |
342
 
c     %---------------%
343
 
c
344
 
      end