Numérotation des lignes :

dsortc
C DSORTC    SOURCE    BP208322  15/10/13    21:15:54     8670c-----------------------------------------------------------------------c\BeginDoccc\Name: dsortccc\Description:c  Sorts the complex array in XREAL and XIMAG into the orderc  specified by WHICH and optionally applies the permutation to thec  real array Y. It is assumed that if an element of XIMAG isc  nonzero, then its negative is also an element. In other words,c  both members of a complex conjugate pair are to be sorted and thec  pairs are kept adjacent to each other.cc\Usage:c  call dsortcc     ( WHICH, APPLY, N, XREAL, XIMAG, Y )cc\Argumentsc  WHICH   Character*2.  (Input)c          'LM' -> sort XREAL,XIMAG into increasing order of magnitude.c          'SM' -> sort XREAL,XIMAG into decreasing order of magnitude.c          'LR' -> sort XREAL into increasing order of algebraic.c          'SR' -> sort XREAL into decreasing order of algebraic.c          'LI' -> sort XIMAG into increasing order of magnitude.c          'SI' -> sort XIMAG into decreasing order of magnitude.c          NOTE: If an element of XIMAG is non-zero, then its negativec                is also an element.cc  APPLY   Logical.  (Input)c          APPLY = .TRUE.  -> apply the sorted order to array Y.c          APPLY = .FALSE. -> do not apply the sorted order to array Y.cc  N       Integer.  (INPUT)c          Size of the arrays.cc  XREAL,  Double precision array of length N.  (INPUT/OUTPUT)c  XIMAG   Real and imaginary part of the array to be sorted.cc  Y       REAL*8 array of length N.  (INPUT/OUTPUT)cc\EndDoccc-----------------------------------------------------------------------cc\BeginLibcc\Authorc     Danny Sorensen               Phuong Vuc     Richard Lehoucq              CRPC / Rice Universityc     Dept. of Computational &     Houston, Texasc     Applied Mathematicsc     Rice Universityc     Houston, Texascc\Revision history:c     xx/xx/92: Version ' 2.1'c               Adapted from the sort routine in LANSO.cc\SCCS Information: @(#)c FILE: sortc.F   SID: 2.3   DATE OF SID: 4/20/96   RELEASE: 2cc\EndLibcc-----------------------------------------------------------------------c      subroutine dsortc (which, apply, n, xreal, ximag, y)cc     %------------------%c     | Scalar Arguments |c     %------------------%c      character*2 which      logical    apply      integer    ncc     %-----------------%c     | Array Arguments |c     %-----------------%c      REAL*8     &           xreal(0:n-1), ximag(0:n-1), y(0:n-1)cc     %---------------%c     | Local Scalars |c     %---------------%c      integer    i, igap, j      REAL*8     &           temp, temp1, temp2cc     %--------------------%c     | External Functions |c     %--------------------%c      REAL*8     &           dlapy2      external   dlapy2cc     %-----------------------%c     | Executable Statements |c     %-----------------------%c      igap = n / 2c      if (which .eq. 'LM') thencc        %------------------------------------------------------%c        | Sort XREAL,XIMAG into increasing order of magnitude. |c        %------------------------------------------------------%c   10    continue         if (igap .eq. 0) go to 9000c         do 30 i = igap, n-1            j = i-igap   20       continuec            if (j.lt.0) go to 30c            temp1 = dlapy2(xreal(j),ximag(j))            temp2 = dlapy2(xreal(j+igap),ximag(j+igap))c            if (temp1.gt.temp2) then                temp = xreal(j)                xreal(j) = xreal(j+igap)                xreal(j+igap) = tempc                temp = ximag(j)                ximag(j) = ximag(j+igap)                ximag(j+igap) = tempc                if (apply) then                    temp = y(j)                    y(j) = y(j+igap)                    y(j+igap) = temp                end if            else                go to 30            end if            j = j-igap            go to 20   30    continue         igap = igap / 2         go to 10c      else if (which .eq. 'SM') thencc        %------------------------------------------------------%c        | Sort XREAL,XIMAG into decreasing order of magnitude. |c        %------------------------------------------------------%c   40    continue         if (igap .eq. 0) go to 9000c         do 60 i = igap, n-1            j = i-igap   50       continuec            if (j .lt. 0) go to 60c            temp1 = dlapy2(xreal(j),ximag(j))            temp2 = dlapy2(xreal(j+igap),ximag(j+igap))c            if (temp1.lt.temp2) then               temp = xreal(j)               xreal(j) = xreal(j+igap)               xreal(j+igap) = tempc               temp = ximag(j)               ximag(j) = ximag(j+igap)               ximag(j+igap) = tempc               if (apply) then                  temp = y(j)                  y(j) = y(j+igap)                  y(j+igap) = temp               end if            else               go to 60            endif            j = j-igap            go to 50   60    continue         igap = igap / 2         go to 40c      else if (which .eq. 'LR') thencc        %------------------------------------------------%c        | Sort XREAL into increasing order of algebraic. |c        %------------------------------------------------%c   70    continue         if (igap .eq. 0) go to 9000c         do 90 i = igap, n-1            j = i-igap   80       continuec            if (j.lt.0) go to 90c            if (xreal(j).gt.xreal(j+igap)) then               temp = xreal(j)               xreal(j) = xreal(j+igap)               xreal(j+igap) = tempc               temp = ximag(j)               ximag(j) = ximag(j+igap)               ximag(j+igap) = tempc               if (apply) then                  temp = y(j)                  y(j) = y(j+igap)                  y(j+igap) = temp               end if            else               go to 90            endif            j = j-igap            go to 80   90    continue         igap = igap / 2         go to 70c      else if (which .eq. 'SR') thencc        %------------------------------------------------%c        | Sort XREAL into decreasing order of algebraic. |c        %------------------------------------------------%c  100    continue         if (igap .eq. 0) go to 9000         do 120 i = igap, n-1            j = i-igap  110       continuec            if (j.lt.0) go to 120c            if (xreal(j).lt.xreal(j+igap)) then               temp = xreal(j)               xreal(j) = xreal(j+igap)               xreal(j+igap) = tempc               temp = ximag(j)               ximag(j) = ximag(j+igap)               ximag(j+igap) = tempc               if (apply) then                  temp = y(j)                  y(j) = y(j+igap)                  y(j+igap) = temp               end if            else               go to 120            endif            j = j-igap            go to 110  120    continue         igap = igap / 2         go to 100c      else if (which .eq. 'LI') thencc        %------------------------------------------------%c        | Sort XIMAG into increasing order of magnitude. |c        %------------------------------------------------%c  130    continue         if (igap .eq. 0) go to 9000         do 150 i = igap, n-1            j = i-igap  140       continuec            if (j.lt.0) go to 150c            if (abs(ximag(j)).gt.abs(ximag(j+igap))) then               temp = xreal(j)               xreal(j) = xreal(j+igap)               xreal(j+igap) = tempc               temp = ximag(j)               ximag(j) = ximag(j+igap)               ximag(j+igap) = tempc               if (apply) then                  temp = y(j)                  y(j) = y(j+igap)                  y(j+igap) = temp               end if            else               go to 150            endif            j = j-igap            go to 140  150    continue         igap = igap / 2         go to 130c      else if (which .eq. 'SI') thencc        %------------------------------------------------%c        | Sort XIMAG into decreasing order of magnitude. |c        %------------------------------------------------%c  160    continue         if (igap .eq. 0) go to 9000         do 180 i = igap, n-1            j = i-igap  170       continuec            if (j.lt.0) go to 180c            if (abs(ximag(j)).lt.abs(ximag(j+igap))) then               temp = xreal(j)               xreal(j) = xreal(j+igap)               xreal(j+igap) = tempc               temp = ximag(j)               ximag(j) = ximag(j+igap)               ximag(j+igap) = tempc               if (apply) then                  temp = y(j)                  y(j) = y(j+igap)                  y(j+igap) = temp               end if            else               go to 180            endif            j = j-igap            go to 170  180    continue         igap = igap / 2         go to 160      end ifc 9000 continue      returncc     %---------------%c     | End of dsortc |c     %---------------%c      end  

© Cast3M 2003 - Tous droits réservés.
Mentions légales