
C	subroutines to sort REAL*16 numbers with "quicksort" and "bubblesort"; 
C       (IARR : old order-nr)

C	sven ivansson   1985





	SUBROUTINE QUICKXD (DARR,IARR,N)

	REAL*16 DARR(*)
	INTEGER IARR(*)

	PARAMETER ( IQDIM = 20000 )
	INTEGER IQ1(iqdim), IQ2(iqdim)
	REAL*16 DCOMPARE

	IF ((1+(N/50)).GT.iqdim) STOP 572            !  can be modified
	IF (N.LE.1) THEN
	  IF (N.EQ.1) IARR(1) = 1
	  RETURN
	END IF

	DO I = 1, N
	  IARR(I) = I
	END DO

	IQUEUE = 1
	IQ1(1) = 1
	IQ2(1) = N
	DO WHILE (IQUEUE.GT.0)
	  ILEFT = IQ1(IQUEUE)
	  IRIGHT = IQ2(IQUEUE)
	  IL = ILEFT
	  IR = IRIGHT
	  ILR = (IL+IR)/2
	  DCOMPARE = DARR(ILR)     ! ta ev ngn mer omsorgsfullt vald komponent
          ICOMP = IARR(ILR)
	  DARR(ILR) = DARR(IL)         ! lucka skapad vid plats "il"
	  IARR(ILR) = IARR(IL)
	  DO WHILE (.TRUE.)
	    DO WHILE (DARR(IR).GE.DCOMPARE .AND. IR.GT.IL)
	      IR = IR-1
	    END DO
	    IF (IL.EQ.IR) GOTO 11
	    DARR(IL) = DARR(IR)
	    IARR(IL) = IARR(IR)
	    IL = IL+1
	    DO WHILE (DARR(IL).LE.DCOMPARE .AND. IL.LT.IR)
	      IL = IL+1
	    END DO
	    IF (IR.EQ.IL) GOTO 11
	    DARR(IR) = DARR(IL)
	    IARR(IR) = IARR(IL)
	    IR = IR-1
	  END DO
11	  ILR = IL    ! samma som "ir"
	  DARR(ILR) = DCOMPARE
	  IARR(ILR) = ICOMP
	  IQUEUE = IQUEUE-1
	  IR = ILR-1
	  IF (IR.GT.ILEFT) THEN
	    LENGTH = 1+IR-ILEFT
	    IF (LENGTH.GT.50) THEN
	      IQUEUE = IQUEUE+1
	      IQ1(IQUEUE) = ILEFT
	      IQ2(IQUEUE) = IR
	    ELSE
	      CALL BUBBLEXD (DARR,IARR,ILEFT,IR)
	    END IF
	  END IF
	  IL = ILR+1
	  IF (IL.LT.IRIGHT) THEN
	    LENGTH = 1+IRIGHT-IL
	    IF (LENGTH.GT.50) THEN
	      IQUEUE = IQUEUE+1
	      IQ1(IQUEUE) = IL
	      IQ2(IQUEUE) = IRIGHT
	    ELSE
	      CALL BUBBLEXD (DARR,IARR,IL,IRIGHT)
	    END IF
	  END IF
	END DO

	RETURN
	END





	SUBROUTINE BUBBLEXD (DARR,IARR,ILEFT,IRIGHT)
	REAL*16 DARR(*)
	INTEGER IARR(*)
	REAL*16 DSLASK

	DO I = (ILEFT+1), IRIGHT
	  DO J = IRIGHT, I, -1
	    IF (DARR(J-1).GT.DARR(J)) THEN
	      DSLASK = DARR(J-1)
	      ISLASK = IARR(J-1)
	      DARR(J-1) = DARR(J)
	      IARR(J-1) = IARR(J)
	      DARR(J) = DSLASK
	      IARR(J) = ISLASK
	    END IF
	  END DO
	END DO

	RETURN
	END
