dfac_mem_stack_aux.F

Go to the documentation of this file.

C
C  This file is part of MUMPS 5.5.1, released
C  on Tue Jul 12 13:17:24 UTC 2022
C
C
C  Copyright 1991-2022 CERFACS, CNRS, ENS Lyon, INP Toulouse, Inria,
C  Mumps Technologies, University of Bordeaux.
C
C  This version of MUMPS is provided to you free of charge. It is
C  released under the CeCILL-C license 
C  (see doc/CeCILL-C_V1-en.txt, doc/CeCILL-C_V1-fr.txt, and
C  https://cecill.info/licences/Licence_CeCILL-C_V1-en.html)
C
      SUBROUTINE dmumps_compact_factors(A, LDA, NPIV, NBROW, KEEP,
     &           SIZEA, IW )
      IMPLICIT NONE
      INTEGER,    INTENT(IN) :: LDA, NPIV, NBROW
      INTEGER(8), INTENT(IN) :: SIZEA
      INTEGER, INTENT(IN) :: IW( NPIV )
      INTEGER :: KEEP(500)
      DOUBLE PRECISION :: A(SIZEA)
      INTEGER(8) :: IOLD, INEW, J8
      INTEGER I , ILAST
      INTEGER NBROW_L_RECTANGLE_TO_MOVE
      INTEGER :: ICOL_BEG, ICOL_END, NBPANELS, NB_TARGET
      INTEGER :: NBCOLS_PANEL, NBROWS_PANEL
      IF ( npiv .EQ. 0 ) GOTO 500
      nb_target = npiv
      IF ( keep(459) .GT. 1 .AND. keep(50) .NE. 0 ) THEN
        CALL mumps_ldltpanel_nbtarget( npiv, nb_target, keep )
      ENDIF
      IF ( keep(50) .EQ.0 .OR. nb_target .EQ. npiv ) THEN
        IF (lda.EQ.npiv) GOTO 500
        IF ( keep(50) .NE. 0 ) THEN
          iold = int(lda  + 1,8)
          inew = int(npiv + 1,8)
          IF (iold .EQ. inew ) THEN
            WRITE(*,*) " Internal error in DMUMPS_COMPACT_FACTORS",
     &      iold, inew, npiv
            CALL mumps_abort()
          ENDIF
          DO i = 1, npiv - 1
            IF ( i .LE. npiv-2 ) THEN
              ilast = i+1
            ELSE
              ilast = i
            ENDIF
            DO j8 = 0_8, int(ilast,8)
              a( inew + j8 ) = a( iold + j8 )
            END DO
            inew = inew + int(npiv,8)
            iold = iold + int(lda,8)
          END DO
          nbrow_l_rectangle_to_move = nbrow
        ELSE 
          inew = 1_8 + int(npiv,8) * int(lda + 1,8)
          iold = 1_8 + int(lda,8) * int(npiv +1,8)
          nbrow_l_rectangle_to_move = nbrow - 1
        ENDIF
      ELSE
        icol_beg = 1               
        nbpanels  = 0
        inew = 1_8
        nbrows_panel = npiv
        DO WHILE ( icol_beg .LE. npiv )
          nbpanels=nbpanels + 1
          icol_end = min(nb_target * nbpanels, npiv)
          IF ( iw( icol_end ) < 0 ) THEN 
            icol_end = icol_end + 1
          ENDIF
          nbcols_panel = icol_end - icol_beg + 1
          iold = int(icol_beg-1,8) * int(lda,8) + int(icol_beg,8)
          DO i =1, nbrows_panel
            IF (iold .NE. inew) THEN
              DO j8=0, min(i+1, nbcols_panel)-1 
              a(inew+j8) = a(iold+j8)
              ENDDO
            ENDIF
            inew = inew + int(nbcols_panel,8)
            iold = iold + int(lda,8)
          ENDDO
          nbrows_panel = nbrows_panel - nbcols_panel
          icol_beg = icol_end + 1
        ENDDO
        iold = 1_8 + int(lda,8)*int(npiv,8)
        nbrow_l_rectangle_to_move = nbrow
      ENDIF
      DO i = 1, nbrow_l_rectangle_to_move
         DO j8 = 0_8, int(npiv - 1,8)
           a( inew + j8 ) = a( iold + j8 )
         END DO
         inew = inew + int(npiv,8)
         iold = iold + int(lda,8)
      ENDDO
 500  RETURN
      END SUBROUTINE dmumps_compact_factors
      SUBROUTINE dmumps_compact_factors_unsym(A, LDA, NPIV, NCONTIG,
     &                                        SIZEA )
      IMPLICIT NONE
      INTEGER,    INTENT(IN)     :: NCONTIG, NPIV, LDA
      INTEGER(8), INTENT(IN)     :: SIZEA
      DOUBLE PRECISION,    INTENT(INOUT)  :: A(SIZEA)
      INTEGER I, J
      INTEGER(8) :: INEW, IOLD
      inew = int(npiv+1,8)
      iold = int(lda+1,8)
      DO i = 2, ncontig
        DO j = 1, npiv
          a(inew)=a(iold)
          inew = inew + 1_8
          iold = iold + 1_8
        ENDDO
        iold = iold + int(lda - npiv,8)
      ENDDO
      RETURN
      END SUBROUTINE dmumps_compact_factors_unsym
      SUBROUTINE dmumps_copy_cb_right_to_left( A, LA, LDA, POSELT,
     &           IPTRLU, NPIV,
     &           NBCOL_STACK, NBROW_STACK,
     &           NBROW_SEND, SIZECB, KEEP, PACKED_CB,
     &           LAST_ALLOWED, NBROW_ALREADY_STACKED )
      IMPLICIT NONE
      INTEGER(8), intent (in) :: POSELT, IPTRLU, LA, SIZECB
      LOGICAL, intent (in) :: PACKED_CB
      DOUBLE PRECISION A(LA)
      INTEGER, intent(in):: LDA, NPIV, NBCOL_STACK, NBROW_STACK,
     &                      NBROW_SEND
      INTEGER, intent(inout) :: NBROW_ALREADY_STACKED
      INTEGER(8), intent(in)    :: LAST_ALLOWED
      INTEGER(8) :: APOS, NPOS
      INTEGER NBROW
      INTEGER(8) :: J
      INTEGER I, KEEP(500)
#if defined(ZERO_TRIANGLE)
      DOUBLE PRECISION ZERO
        parameter( zero = 0.0d0 )
#endif
      nbrow = nbrow_stack + nbrow_send
      IF (nbrow_stack .NE. 0 ) THEN
        npos = iptrlu + sizecb         
        apos = poselt + int(npiv+nbrow,8) 
     &       * int(lda,8) - 1_8 
        IF ( keep(50) .EQ. 0 .OR. .NOT. packed_cb ) THEN
          apos = apos - int(lda,8) * int(nbrow_already_stacked,8)
          npos = npos
     &         - int(nbcol_stack,8) * int(nbrow_already_stacked,8)
        ELSE
          apos = apos - int(lda - 1,8) * int(nbrow_already_stacked,8)
          npos = npos - ( int(nbrow_already_stacked,8) *
     &                    int(nbrow_already_stacked+1,8) ) / 2_8
        ENDIF
        DO i = nbrow - nbrow_already_stacked, nbrow_send+1, -1
          IF (keep(50).EQ.0) THEN
            IF ( npos - int(nbcol_stack,8) + 1_8 .LT.
     &                                  last_allowed ) THEN
              EXIT
            ENDIF
            DO j= 1_8,int(nbcol_stack,8)
              a(npos-j+1_8) = a(apos-j+1_8)
            ENDDO
            npos = npos - int(nbcol_stack,8)
          ELSE
            IF (.NOT. packed_cb) THEN
              IF ( npos - int(nbcol_stack,8) + 1_8 .LT.
     &                                  last_allowed ) THEN
                EXIT
              ENDIF
#if defined(ZERO_TRIANGLE)
              DO j = 1_8, int(nbcol_stack - i,8)
                a(npos - j + 1_8) = zero
              END DO
#endif
              npos = npos + int(- nbcol_stack + i,8)
            ENDIF
            IF ( npos - int(i,8) + 1_8 .LT. last_allowed ) THEN
              EXIT
            ENDIF
            DO j =1_8, int(i,8)
              a(npos-j+1_8) = a(apos-j+1_8)
            ENDDO
            npos = npos - int(i,8)
          ENDIF
          IF (keep(50).EQ.0) THEN
            apos = apos - int(lda,8)
          ELSE
            apos = apos - int(lda + 1,8)
          ENDIF
          nbrow_already_stacked = nbrow_already_stacked + 1
        ENDDO
      END IF
      RETURN
      END SUBROUTINE dmumps_copy_cb_right_to_left
      SUBROUTINE dmumps_copy_cb_left_to_right( A, LA, LDA, POSELT,
     &           IPTRLU, NPIV,
     &           NBCOL_STACK, NBROW_STACK,
     &           NBROW_SEND, SIZECB, KEEP, PACKED_CB)
      IMPLICIT NONE
      INTEGER(8), intent (in) :: POSELT, IPTRLU, LA, SIZECB
      LOGICAL, intent (in) :: PACKED_CB
      DOUBLE PRECISION A(LA)
      INTEGER, intent(in):: LDA, NPIV, NBCOL_STACK, NBROW_STACK,
     &                      nbrow_send
      INTEGER(8) :: APOS, NPOS, APOS_ini, NPOS_ini
      INTEGER I, KEEP(500)
      INTEGER(8) :: J, LDA8
#if defined(ZERO_TRIANGLE)
      DOUBLE PRECISION ZERO
        parameter( zero = 0.0d0 )
#endif
      lda8 = int(lda,8)
      npos_ini = iptrlu + 1_8
      apos_ini = poselt + int(npiv+nbrow_send,8)* lda8 + int(npiv,8)
!$OMP PARALLEL DO PRIVATE(J, NPOS, APOS) IF (NBROW_STACK > KEEP(360))
      DO i = 1, nbrow_stack
        IF (packed_cb) THEN
            npos = npos_ini + int(i-1,8) * int(i,8)/2_8 +
     &             int(i-1,8) * int(nbrow_send,8)
        ELSE
            npos = npos_ini + int(i-1,8) * int(nbcol_stack,8)
        ENDIF
        apos  =  apos_ini + int(i-1,8) * lda8
        IF (keep(50).EQ.0) THEN
          DO j = 1_8, int(nbcol_stack,8)
            a(npos+j-1_8) = a(apos+j-1_8)
          ENDDO
        ELSE
          DO j  = 1_8, int(i + nbrow_send,8)
            a(npos+j-1_8)=a(apos+j-1_8)
          ENDDO
#if defined(ZERO_TRIANGLE)
          IF (.NOT. packed_cb) THEN
            a(npos+int(i+nbrow_send,8):
     &        npos+int(nbcol_stack-1,8))=zero
          ENDIF
#endif
        ENDIF
      ENDDO
!$OMP END PARALLEL DO
      RETURN
      END SUBROUTINE dmumps_copy_cb_left_to_right