*************************************************
*** MT1 FORTRAN LIBRARY                       ***
*** Written by Aart J.C. Bik (PhD Student)    ***
*** Copyright (C) 1996                        ***
*************************************************

*************************************************
*** CEILING/FLOOR Functions                   ***
*************************************************

      INTEGER FUNCTION ICEIL__(D, N)
      implicit none
      INTEGER D, N

      ICEIL__ = D / N
      IF (MOD(D,N).NE.0) THEN
        IF ((D.GT.0).EQV.(N.GT.0)) THEN
          ICEIL__ = ICEIL__ + 1
        ENDIF
      ENDIF
      RETURN
      END

      INTEGER FUNCTION IFLOOR__(D, N)
      implicit none
      INTEGER D, N

      IFLOOR__ = D / N 
      IF (MOD(D,N).NE.0) THEN
        IF ((D.LT.0).NEQV.(N.LT.0)) THEN
          IFLOOR__ = IFLOOR__ - 1 
        ENDIF
      ENDIF
      RETURN
      END

*************************************************
*** Initialization Routines                   ***
*************************************************

      SUBROUTINE IINI__(VAL, TMP, IND, LOW, HGH, NP, SZ, LST)
      implicit none
      INTEGER NP, SZ, LST, I, J, K
      INTEGER TMP(SZ), IND(SZ), LOW(NP), HGH(NP)
      INTEGER VAL(SZ), V

      VAL(1) = 0

      DO I = 1, NP
        LOW(I) = 0
      ENDDO

      DO K = 2, LST 
        I      = TMP(K)
        LOW(I) = LOW(I) + 1
      ENDDO

      LOW(1) = LOW(1) + 2
      HGH(1) = LOW(1) - 1
      DO I = 2, NP
        LOW(I) = LOW(I) + LOW(I-1)
        HGH(I) = LOW(I) - 1
      ENDDO

      DO K = 2, LST
        IF (IND(K).GT.0) THEN
          I      = TMP(K)
          LOW(I) = LOW(I) - 1
          J      = LOW(I)

          DO WHILE (J.NE.K)
            V = VAL(J)
            I = IND(J)

            VAL(J) =   VAL(K)
            IND(J) = - IND(K)

            VAL(K) = V
            IND(K) = I

            I      = TMP(J)
            LOW(I) = LOW(I) - 1
            J      = LOW(I)
          ENDDO
        ELSE 
          IND(K) = - IND(K)
        ENDIF
      ENDDO

      DO K = LST+1, SZ
        IND(K) = 0
      ENDDO

      RETURN
      END

      SUBROUTINE SINI__(VAL, TMP, IND, LOW, HGH, NP, SZ, LST)
      implicit none
      INTEGER NP, SZ, LST, I, J, K
      INTEGER TMP(SZ), IND(SZ), LOW(NP), HGH(NP)
      REAL    VAL(SZ), V

      VAL(1) = 0.0

      DO I = 1, NP
        LOW(I) = 0
      ENDDO

      DO K = 2, LST 
        I      = TMP(K)
        LOW(I) = LOW(I) + 1
      ENDDO

      LOW(1) = LOW(1) + 2
      HGH(1) = LOW(1) - 1
      DO I = 2, NP
        LOW(I) = LOW(I) + LOW(I-1)
        HGH(I) = LOW(I) - 1
      ENDDO

      DO K = 2, LST
        IF (IND(K).GT.0) THEN
          I      = TMP(K)
          LOW(I) = LOW(I) - 1
          J      = LOW(I)

          DO WHILE (J.NE.K)
            V = VAL(J)
            I = IND(J)

            VAL(J) =   VAL(K)
            IND(J) = - IND(K)

            VAL(K) = V
            IND(K) = I

            I      = TMP(J)
            LOW(I) = LOW(I) - 1
            J      = LOW(I)
          ENDDO
        ELSE 
          IND(K) = - IND(K)
        ENDIF
      ENDDO

      DO K = LST+1, SZ
        IND(K) = 0
      ENDDO

      RETURN
      END

      SUBROUTINE DINI__(VAL, TMP, IND, LOW, HGH, NP, SZ, LST)
      implicit none
      INTEGER          NP, SZ, LST, I, J, K
      INTEGER          TMP(SZ), IND(SZ), LOW(NP), HGH(NP)
      DOUBLE PRECISION VAL(SZ), V

      VAL(1) = 0.0D0

      DO I = 1, NP
        LOW(I) = 0
      ENDDO

      DO K = 2, LST 
        I      = TMP(K)
        LOW(I) = LOW(I) + 1
      ENDDO

      LOW(1) = LOW(1) + 2
      HGH(1) = LOW(1) - 1
      DO I = 2, NP
        LOW(I) = LOW(I) + LOW(I-1)
        HGH(I) = LOW(I) - 1
      ENDDO

      DO K = 2, LST
        IF (IND(K).GT.0) THEN
          I      = TMP(K)
          LOW(I) = LOW(I) - 1
          J      = LOW(I)

          DO WHILE (J.NE.K)
            V = VAL(J)
            I = IND(J)

            VAL(J) =   VAL(K)
            IND(J) = - IND(K)

            VAL(K) = V
            IND(K) = I

            I      = TMP(J)
            LOW(I) = LOW(I) - 1
            J      = LOW(I)
          ENDDO
        ELSE 
          IND(K) = - IND(K)
        ENDIF
      ENDDO

      DO K = LST+1, SZ
        IND(K) = 0
      ENDDO

      RETURN
      END

      SUBROUTINE CINI__(VAL, TMP, IND, LOW, HGH, NP, SZ, LST)
      implicit none
      INTEGER NP, SZ, LST, I, J, K
      INTEGER TMP(SZ), IND(SZ), LOW(NP), HGH(NP)
      COMPLEX VAL(SZ), V

      VAL(1) = CMPLX(0.0, 0.0) 

      DO I = 1, NP
        LOW(I) = 0
      ENDDO

      DO K = 2, LST 
        I      = TMP(K)
        LOW(I) = LOW(I) + 1
      ENDDO

      LOW(1) = LOW(1) + 2
      HGH(1) = LOW(1) - 1
      DO I = 2, NP
        LOW(I) = LOW(I) + LOW(I-1)
        HGH(I) = LOW(I) - 1
      ENDDO

      DO K = 2, LST
        IF (IND(K).GT.0) THEN
          I      = TMP(K)
          LOW(I) = LOW(I) - 1
          J      = LOW(I)

          DO WHILE (J.NE.K)
            V = VAL(J)
            I = IND(J)

            VAL(J) =   VAL(K)
            IND(J) = - IND(K)

            VAL(K) = V
            IND(K) = I

            I      = TMP(J)
            LOW(I) = LOW(I) - 1
            J      = LOW(I)
          ENDDO
        ELSE 
          IND(K) = - IND(K)
        ENDIF
      ENDDO

      DO K = LST+1, SZ
        IND(K) = 0
      ENDDO

      RETURN
      END

************************************************
*** Insertion Routines                       ***
************************************************

      SUBROUTINE IINS__(VAL, IND, LOW, HGH, E, NP, SZ, LST, AD, F)
      implicit none
      INTEGER E, NP, SZ, LST, AD, F, I, J, L 
      INTEGER LOW(NP), HGH(NP), IND(SZ)
      INTEGER VAL(SZ)

      I = 0
10    IF ((HGH(E).LT.SZ).AND.(IND(HGH(E)+1).EQ.0)) THEN

        HGH(E)  = HGH(E) + 1
        AD      = HGH(E)

        VAL(AD) = 0 
        IND(AD) = F 

        IF (LST.LT.AD) LST = AD 

      ELSEIF ((LOW(E).GT.2).AND.(IND(LOW(E)-1).EQ.0)) THEN

        LOW(E)  = LOW(E) - 1
        AD      = LOW(E)

        VAL(AD) = 0 
        IND(AD) = F 

      ELSEIF ((SZ-LST).GT.(HGH(E)-LOW(E)+1)) THEN

        AD = LST + 1

        DO I = LOW(E), HGH(E)
           VAL(AD) = VAL(I)
           IND(AD) = IND(I)
           IND(I)  = 0
           AD      = AD + 1
        ENDDO

        VAL(AD) = 0 
        IND(AD) = F 

        LOW(E) = LST + 1 
        HGH(E) = AD 
        LST    = AD

      ELSE

        IF (I.EQ.1) THEN
          PRINT *, 'Out of Memory'
          STOP
        ENDIF

        DO I = 1, NP
          IF (HGH(I).GE.LOW(I)) THEN
            J             = IND( LOW(I) )
            IND( LOW(I) ) = -I
            LOW(I)        = J
          ENDIF
        ENDDO

        J = 2
        DO I = 2, LST
          IF (IND(I).GT.0) THEN
            VAL(J) = VAL(I)
            IND(J) = IND(I)
            J      = J + 1
          ELSEIF (IND(I).LT.0) THEN
            L      = - IND(I)
            VAL(J) = VAL(I)
            IND(J) = LOW(L)
            LOW(L) = J
            HGH(L) = J + HGH(L) - I 
            J      = J + 1
          ENDIF
        ENDDO

        DO I = J, LST
          IND(I) = 0
        ENDDO
        LST = J - 1
        I = 1
        GOTO 10

      ENDIF
      RETURN
      END

      SUBROUTINE SINS__(VAL, IND, LOW, HGH, E, NP, SZ, LST, AD, F)
      implicit none
      INTEGER E, NP, SZ, LST, AD, F, I, J,  L
      INTEGER LOW(NP), HGH(NP), IND(SZ)
      REAL    VAL(SZ)

      I = 0
10    IF ((HGH(E).LT.SZ).AND.(IND(HGH(E)+1).EQ.0)) THEN

        HGH(E)  = HGH(E) + 1
        AD      = HGH(E)

        VAL(AD) = 0.0
        IND(AD) = F

        IF (LST.LT.AD) LST = AD 

      ELSEIF ((LOW(E).GT.2).AND.(IND(LOW(E)-1).EQ.0)) THEN

        LOW(E)  = LOW(E) - 1
        AD      = LOW(E)

        VAL(AD) = 0.0 
        IND(AD) = F 

      ELSEIF ((SZ-LST).GT.(HGH(E)-LOW(E)+1)) THEN

        AD = LST + 1

        DO I = LOW(E), HGH(E)
           VAL(AD) = VAL(I)
           IND(AD) = IND(I)
           IND(I)  = 0
           AD      = AD + 1
        ENDDO

        VAL(AD) = 0.0
        IND(AD) = F 

        LOW(E) = LST + 1 
        HGH(E) = AD 
        LST    = AD

      ELSE

        IF (I.EQ.1) THEN
          PRINT *, 'Out of Memory'
          STOP
        ELSE
          PRINT *, 'Left Compression, Insertion at', E, F
        ENDIF

        DO I = 1, NP
          IF (HGH(I).GE.LOW(I)) THEN
            J             = IND( LOW(I) )
            IND( LOW(I) ) = -I
            LOW(I)        = J
          ENDIF
        ENDDO

        J = 2
        DO I = 2, LST
          IF (IND(I).GT.0) THEN
            VAL(J) = VAL(I)
            IND(J) = IND(I)
            J      = J + 1
          ELSEIF (IND(I).LT.0) THEN
            L      = - IND(I)
            VAL(J) = VAL(I)
            IND(J) = LOW(L)
            LOW(L) = J
            HGH(L) = J + HGH(L) - I 
            J      = J + 1
          ENDIF
        ENDDO

        DO I = J, LST
          IND(I) = 0
          VAL(I) = 0.0
        ENDDO
        LST = J - 1
        I = 1
        GOTO 10

      ENDIF
      RETURN
      END

      SUBROUTINE DINS__(VAL, IND, LOW, HGH, E, NP, SZ, LST, AD, F)
      implicit none
      INTEGER          E, NP, SZ, LST, AD, F, I, J, L 
      INTEGER          LOW(NP), HGH(NP), IND(SZ)
      DOUBLE PRECISION VAL(SZ)

      I = 0
10    IF ((HGH(E).LT.SZ).AND.(IND(HGH(E)+1).EQ.0)) THEN

        HGH(E)  = HGH(E) + 1
        AD      = HGH(E)

        VAL(AD) = 0.0D0 
        IND(AD) = F 

        IF (LST.LT.AD) LST = AD 

      ELSEIF ((LOW(E).GT.2).AND.(IND(LOW(E)-1).EQ.0)) THEN

        LOW(E)  = LOW(E) - 1
        AD      = LOW(E)

        VAL(AD) = 0.0D0
        IND(AD) = F 

      ELSEIF ((SZ-LST).GT.(HGH(E)-LOW(E)+1)) THEN

        AD = LST + 1

        DO I = LOW(E), HGH(E)
           VAL(AD) = VAL(I)
           IND(AD) = IND(I)
           IND(I)  = 0
           AD      = AD + 1
        ENDDO

        VAL(AD) = 0.0D0
        IND(AD) = F 

        LOW(E) = LST + 1 
        HGH(E) = AD 
        LST    = AD

      ELSE

        IF (I.EQ.1) THEN
          PRINT *, 'Out of Memory'
          STOP
        ENDIF

        DO I = 1, NP
          IF (HGH(I).GE.LOW(I)) THEN
            J             = IND( LOW(I) )
            IND( LOW(I) ) = -I
            LOW(I)        = J
          ENDIF
        ENDDO

        J = 2
        DO I = 2, LST
          IF (IND(I).GT.0) THEN
            VAL(J) = VAL(I)
            IND(J) = IND(I)
            J      = J + 1
          ELSEIF (IND(I).LT.0) THEN
            L      = - IND(I)
            VAL(J) = VAL(I)
            IND(J) = LOW(L)
            LOW(L) = J
            HGH(L) = J + HGH(L) - I 
            J      = J + 1
          ENDIF
        ENDDO

        DO I = J, LST
          IND(I) = 0
        ENDDO
        LST = J - 1
        I = 1
        GOTO 10

      ENDIF
      RETURN
      END

      SUBROUTINE CINS__(VAL, IND, LOW, HGH, E, NP, SZ, LST, AD, F)
      implicit none
      INTEGER E, NP, SZ, LST, AD, F, I, J, L 
      INTEGER LOW(NP), HGH(NP), IND(SZ)
      COMPLEX VAL(SZ)

      I = 0
10    IF ((HGH(E).LT.SZ).AND.(IND(HGH(E)+1).EQ.0)) THEN

        HGH(E)  = HGH(E) + 1
        AD      = HGH(E)

        VAL(AD) = CMPLX(0.0, 0.0) 
        IND(AD) = F 

        IF (LST.LT.AD) LST = AD 

      ELSEIF ((LOW(E).GT.2).AND.(IND(LOW(E)-1).EQ.0)) THEN

        LOW(E)  = LOW(E) - 1
        AD      = LOW(E)

        VAL(AD) = CMPLX(0.0, 0.0) 
        IND(AD) = F 

      ELSEIF ((SZ-LST).GT.(HGH(E)-LOW(E)+1)) THEN

        AD = LST + 1

        DO I = LOW(E), HGH(E)
           VAL(AD) = VAL(I)
           IND(AD) = IND(I)
           IND(I)  = 0
           AD      = AD + 1
        ENDDO

        VAL(AD) = CMPLX(0.0, 0.0) 
        IND(AD) = F 

        LOW(E) = LST + 1 
        HGH(E) = AD 
        LST    = AD

      ELSE

        IF (I.EQ.1) THEN
          PRINT *, 'Out of Memory'
          STOP
        ENDIF

        DO I = 1, NP
          IF (HGH(I).GE.LOW(I)) THEN
            J             = IND( LOW(I) )
            IND( LOW(I) ) = -I
            LOW(I)        = J
          ENDIF
        ENDDO

        J = 2
        DO I = 2, LST
          IF (IND(I).GT.0) THEN
            VAL(J) = VAL(I)
            IND(J) = IND(I)
            J      = J + 1
          ELSEIF (IND(I).LT.0) THEN
            L      = - IND(I)
            VAL(J) = VAL(I)
            IND(J) = LOW(L)
            LOW(L) = J
            HGH(L) = J + HGH(L) - I 
            J      = J + 1
          ENDIF
        ENDDO

        DO I = J, LST
          IND(I) = 0
        ENDDO
        LST = J - 1
        I = 1
        GOTO 10

      ENDIF
      RETURN
      END

*************************************************
*** Expansion Routines                        ***
*************************************************

      SUBROUTINE ISCT__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER LOW, HGH, I
      INTEGER IND(*)
      INTEGER VAL(*), P(*)
      LOGICAL S(*)

      DO I = LOW, HGH
        P( IND(I) ) = VAL(I)
        S( IND(I) ) = .TRUE. 
      ENDDO
      RETURN
      END

      SUBROUTINE SSCT__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER LOW, HGH, I
      INTEGER IND(*)
      REAL    VAL(*), P(*)
      LOGICAL S(*)

      DO I = LOW, HGH
        P( IND(I) ) = VAL(I)
        S( IND(I) ) = .TRUE. 
      ENDDO
      RETURN
      END

      SUBROUTINE DSCT__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER          LOW, HGH, I
      INTEGER          IND(*)
      DOUBLE PRECISION VAL(*), P(*)
      LOGICAL          S(*)

      DO I = LOW, HGH
        P( IND(I) ) = VAL(I)
        S( IND(I) ) = .TRUE. 
      ENDDO
      RETURN
      END

      SUBROUTINE CSCT__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER LOW, HGH, I
      INTEGER IND(*)
      COMPLEX VAL(*), P(*)
      LOGICAL S(*)

      DO I = LOW, HGH
        P( IND(I) ) = VAL(I)
        S( IND(I) ) = .TRUE. 
      ENDDO
      RETURN
      END

*************************************************
*** Compression Routines                      ***
*************************************************

      SUBROUTINE IGTH__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER LOW, HGH, I
      INTEGER IND(*)
      INTEGER VAL(*), P(*)
      LOGICAL S(*)

      DO I = LOW, HGH
        VAL(I)      = P( IND(I) ) 
        P( IND(I) ) = 0
        S( IND(I) ) = .FALSE. 
      ENDDO
      RETURN
      END

      SUBROUTINE SGTH__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER LOW, HGH, I
      INTEGER IND(*)
      REAL    VAL(*), P(*)
      LOGICAL S(*)

      DO I = LOW, HGH
        VAL(I)      = P( IND(I) ) 
        P( IND(I) ) = 0.0
        S( IND(I) ) = .FALSE. 
      ENDDO
      RETURN
      END

      SUBROUTINE DGTH__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER          LOW, HGH, I
      INTEGER          IND(*)
      DOUBLE PRECISION VAL(*), P(*)
      LOGICAL          S(*)

      DO I = LOW, HGH
        VAL(I)      = P( IND(I) ) 
        P( IND(I) ) = 0.0D0
        S( IND(I) ) = .FALSE. 
      ENDDO
      RETURN
      END

      SUBROUTINE CGTH__(VAL, IND, LOW, HGH, P, S)
      implicit none
      INTEGER LOW, HGH, I
      INTEGER IND(*)
      COMPLEX VAL(*), P(*)
      LOGICAL S(*)

      DO I = LOW, HGH
        VAL(I)      = P( IND(I) ) 
        P( IND(I) ) = CMPLX(0.0, 0.0) 
        S( IND(I) ) = .FALSE. 
      ENDDO
      RETURN
      END

*************************************************
*** Lookup Function                           ***
*************************************************

      INTEGER FUNCTION LKP__(IND, LOW, HGH, F)
      implicit none
      INTEGER LOW, HGH, F 
      INTEGER IND(*)

      DO LKP__ = LOW, HGH
        IF (IND(LKP__).EQ.F) GOTO 10 
      ENDDO
      LKP__ = 1 
10    RETURN 
      END

*************************************************
*** Output Routines                    S Only ***
*************************************************

      SUBROUTINE SSHW__(VAL, IND, LOW, HGH, NP, SZ, LST)
      implicit none
      INTEGER NP, SZ, LST, I, J
      INTEGER IND(SZ), LOW(NP), HGH(NP)
      REAL    VAL(SZ)

      PRINT *, '*** SHOW *** Last Used Element ', LST
      DO I = 1, NP
        PRINT *, 'Access Pattern (', I, ')'
        PRINT 1, (VAL(J), J=LOW(I),HGH(I))
        PRINT 2, (IND(J), J=LOW(I),HGH(I))
      ENDDO

1     FORMAT (30(F5.1))
2     FORMAT (30(I5))

      RETURN
      END

      SUBROUTINE SDMP__(VAL, IND, LOW, HGH, NP, SZ, LST)
      implicit none
      INTEGER NP, SZ, LST, I
      INTEGER IND(SZ), LOW(NP), HGH(NP)
      REAL    VAL(SZ)

      PRINT *, '*** DUMP *** Last Used Element ', LST

      PRINT 2, (LOW(I), I=1, NP)
      PRINT 2, (HGH(I), I=1, NP)

      PRINT *, '***'

      PRINT 1, (VAL(I), I=1, SZ)
      PRINT 2, (IND(I), I=1, SZ)

      PRINT *, '***'

1     FORMAT (30(F5.1))
2     FORMAT (30(I5))

      RETURN 
      END
