← Back to branch summary

~ubuntu-branches/ubuntu/hoary/scilab/hoary

« back to all changes in this revision

Viewing changes to routines/arpack/zsortc.f

  • Committer: Bazaar Package Importer
  • Author(s): Torsten Werner
  • Date: 2005-01-09 22:58:21 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20050109225821-473xr8vhgugxxx5j
Tags: 3.0-12
http://bugs.debian.org/255869
changed configure.in to build scilab's own malloc.o, closes: #255869

Show diffs side-by-side

added added

removed removed

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