sfac_root_parallel.F File Reference

Go to the source code of this file.

Functions/Subroutines
subroutine	smumps_facto_root (mpa, myid, master_of_root, root, n, iroot, comm, iw, liw, ifree, a, la, ptrast, ptlust_s, ptrfac, step, info, ldlt, qr, wk, lwk, keep, keep8, dkeep, opeliw, det_exp, det_mant, det_sign)

Function/Subroutine Documentation

smumps_facto_root()

subroutine smumps_facto_root	(	integer, intent(in)	mpa,
		integer	myid,
		integer	master_of_root,
		type ( smumps_root_struc )	root,
		integer	n,
		integer	iroot,
		integer	comm,
		integer, dimension( liw )	iw,
		integer	liw,
		integer	ifree,
		real, dimension( la )	a,
		integer(8)	la,
		integer(8), dimension(keep(28))	ptrast,
		integer, dimension(keep(28))	ptlust_s,
		integer(8), dimension(keep(28))	ptrfac,
		integer, dimension(n)	step,
		integer, dimension( 2 )	info,
		integer	ldlt,
		integer	qr,
		real, dimension( lwk )	wk,
		integer(8)	lwk,
		integer, dimension(500)	keep,
		integer(8), dimension(150)	keep8,
		real, dimension(230)	dkeep,
		double precision, intent(inout)	opeliw,
		integer, intent(inout)	det_exp,
		real, intent(inout)	det_mant,
		integer, intent(inout)	det_sign )

Definition at line 14 of file sfac_root_parallel.F.

        USE smumps_lr_stats, ONLY: upd_flop_root
        USE smumps_struc_def, ONLY : smumps_root_struc
         IMPLICIT NONE
      include 'mpif.h'
      TYPE ( SMUMPS_ROOT_STRUC ) :: root
      INTEGER, INTENT(IN) :: MPA
      INTEGER N, IROOT, COMM, LIW, MYID, IFREE, MASTER_OF_ROOT
      INTEGER(8) :: LA
      INTEGER(8) :: LWK
      REAL WK( LWK )
      INTEGER KEEP(500)
      REAL    DKEEP(230)
      INTEGER(8) KEEP8(150)
      INTEGER(8) :: PTRAST(KEEP(28))
      INTEGER(8) :: PTRFAC(KEEP(28))
      INTEGER PTLUST_S(KEEP(28)), STEP(N), IW( LIW )
      INTEGER INFO( 2 ), LDLT, QR
      REAL A( LA )
      DOUBLE PRECISION, intent(inout) :: OPELIW
      INTEGER, INTENT(INOUT) :: DET_SIGN, DET_EXP
      REAL, INTENT(INOUT) :: DET_MANT
      INTEGER IOLDPS
      INTEGER(8) :: IAPOS
      DOUBLE PRECISION :: FLOPS_ROOT
      INTEGER(8) :: ENTRIES_ROOT
      INTEGER LOCAL_M, LOCAL_N, LPIV, IERR, allocok
      INTEGER FWD_LOCAL_N_RHS, FWD_MTYPE 
      include 'mumps_headers.h'
      EXTERNAL numroc
      INTEGER numroc
        IF ( .NOT. root%yes ) RETURN
        IF ( keep(60) .NE. 0 ) THEN
          IF ((ldlt == 1 .OR. ldlt == 2) .AND. keep(60) == 3 ) THEN
            CALL smumps_symmetrize( wk, root%MBLOCK,
     &      root%MYROW, root%MYCOL, root%NPROW, root%NPCOL,
     &      root%SCHUR_POINTER(1),
     &      root%SCHUR_LLD, root%SCHUR_NLOC,
     &      root%TOT_ROOT_SIZE, myid, comm )
          ENDIF
        RETURN
        ENDIF
        IF (mpa.GT.0) THEN
         IF (myid.EQ.master_of_root) THEN
           CALL  mumps_get_flops_cost
     &     (root%TOT_ROOT_SIZE, root%TOT_ROOT_SIZE, root%TOT_ROOT_SIZE,
     &      ldlt, 3, flops_root)
           WRITE(mpa,'(A, A, 1PD10.3)')
     &     " ... Start processing the root node with ScaLAPACK, ", 
     &     " remaining flops                = ", flops_root
         ENDIF
        ENDIF
        ioldps  = ptlust_s(step(iroot))+keep(ixsz)
        iapos   = ptrast(step(iroot))
        local_m = iw( ioldps + 2 )
        local_n = iw( ioldps + 1 )
        iapos = ptrfac(iw( ioldps + 4 ))
        IF ( ldlt.EQ.0 .OR. ldlt.EQ.2 .OR. qr.ne.0 ) THEN
         lpiv = local_m + root%MBLOCK
        ELSE
         lpiv = 1
        END IF
        IF (associated( root%IPIV )) DEALLOCATE(root%IPIV)
        root%LPIV = lpiv
        ALLOCATE( root%IPIV( lpiv ), stat = allocok )
        IF ( allocok .GT. 0 ) THEN
          info(1) = -13
          info(2) = lpiv
          WRITE(*,*) myid,': problem allocating IPIV(',lpiv,') in root'
          CALL mumps_abort()
        END IF
        CALL descinit( root%DESCRIPTOR(1), root%TOT_ROOT_SIZE,
     &      root%TOT_ROOT_SIZE, root%MBLOCK, root%NBLOCK,
     &      0, 0, root%CNTXT_BLACS, local_m, ierr )
        IF ( ldlt.EQ.2 ) THEN
           IF(root%MBLOCK.NE.root%NBLOCK) THEN
              WRITE(*,*) ' Error: symmetrization only works for'
              WRITE(*,*) ' square block sizes, MBLOCK/NBLOCK=',
     &        root%MBLOCK, root%NBLOCK
              CALL mumps_abort()
            END IF
            IF ( lwk .LT. min(
     &           int(root%MBLOCK,8) * int(root%NBLOCK,8),
     &           int(root%TOT_ROOT_SIZE,8)* int(root%TOT_ROOT_SIZE,8 )
     &           )) THEN
               WRITE(*,*) 'Not enough workspace for symmetrization.'
               CALL mumps_abort()
            END IF
            CALL smumps_symmetrize( wk, root%MBLOCK,
     &      root%MYROW, root%MYCOL, root%NPROW, root%NPCOL,
     &      a( iapos ), local_m, local_n,
     &      root%TOT_ROOT_SIZE, myid, comm )
        END IF
        IF (ldlt.EQ.0.OR.ldlt.EQ.2) THEN
          CALL psgetrf( root%TOT_ROOT_SIZE, root%TOT_ROOT_SIZE,
     &      a( iapos ),
     &      1, 1, root%DESCRIPTOR(1), root%IPIV(1), ierr )
          IF ( ierr .GT. 0 ) THEN
              info(1)=-10
              info(2)=ierr-1
          END IF
        ELSE
          CALL pspotrf('L',root%TOT_ROOT_SIZE,a(iapos),
     &      1,1,root%DESCRIPTOR(1),ierr)
            IF ( ierr .GT. 0 ) THEN
              info(1)=-40
              info(2)=ierr-1
            END IF
        END IF
        IF (ierr .GT. 0) THEN
          CALL mumps_update_flops_root( opeliw, ldlt,
     &                          root%TOT_ROOT_SIZE, info(2),
     &                          root%NPROW, root%NPCOL, myid )
          IF (keep(486) .GT. 0) THEN
            CALL upd_flop_root( ldlt,
     &                          root%TOT_ROOT_SIZE, info(2),
     &                          root%NPROW, root%NPCOL, myid )
          ENDIF
        ELSE
          CALL mumps_update_flops_root( opeliw, ldlt,
     &                          root%TOT_ROOT_SIZE, root%TOT_ROOT_SIZE,
     &                          root%NPROW, root%NPCOL, myid )
          IF (keep(486) .GT. 0) THEN
            CALL upd_flop_root( ldlt,
     &                          root%TOT_ROOT_SIZE, root%TOT_ROOT_SIZE,
     &                          root%NPROW, root%NPCOL, myid )
          ENDIF
        ENDIF
        IF ( ldlt .EQ. 0 ) THEN
          entries_root = int(root%TOT_ROOT_SIZE,8)
     &                 * int(root%TOT_ROOT_SIZE,8)
        ELSE
          entries_root = int(root%TOT_ROOT_SIZE,8)
     &                 * int(root%TOT_ROOT_SIZE+1,8)/2_8
        ENDIF
        keep8(10)=keep8(10) + entries_root /
     &                        int(root%NPROW * root%NPCOL,8)
        IF (myid .eq. master_of_root) THEN
          keep8(10)=keep8(10) +
     &    mod(entries_root, int(root%NPROW*root%NPCOL,8))
        ENDIF
        CALL smumps_par_root_minmax_piv_upd (
     &         root%MBLOCK, root%IPIV(1),root%MYROW,
     &         root%MYCOL, root%NPROW, root%NPCOL, a(iapos), local_m,
     &         local_n, root%TOT_ROOT_SIZE, myid, dkeep, keep, ldlt)
        IF (keep(258).NE.0) THEN
           IF (root%MBLOCK.NE.root%NBLOCK) THEN
            write(*,*) "Internal error in SMUMPS_FACTO_ROOT:",
     &      "Block size different for rows and columns",
     &      root%MBLOCK, root%NBLOCK
            CALL mumps_abort()
          ENDIF
          CALL smumps_getdeter2d(root%MBLOCK, root%IPIV(1),root%MYROW,
     &         root%MYCOL, root%NPROW, root%NPCOL, a(iapos), local_m,
     &         local_n, root%TOT_ROOT_SIZE, myid, det_mant, det_exp,
     &         ldlt)
        ENDIF
        IF (keep(252) .NE. 0) THEN
          fwd_local_n_rhs = numroc(keep(253), root%NBLOCK,
     &                      root%MYCOL, 0, root%NPCOL)
          fwd_local_n_rhs = max(1,fwd_local_n_rhs)
          fwd_mtype       = 1 
          CALL smumps_solve_2d_bcyclic(
     &         root%TOT_ROOT_SIZE,
     &         keep(253),          
     &         fwd_mtype,
     &         a(iapos),
     &         root%DESCRIPTOR(1),
     &         local_m, local_n, fwd_local_n_rhs,
     &         root%IPIV(1), lpiv,
     &         root%RHS_ROOT(1,1), ldlt,
     &         root%MBLOCK, root%NBLOCK,
     &         root%CNTXT_BLACS, ierr)
        ENDIF
        RETURN