dsol_root_parallel.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_root_solve( NRHS, DESCA_PAR, 
     &  CNTXT_PAR,LOCAL_M,LOCAL_N,MBLOCK,NBLOCK,
     &  IPIV,LPIV,MASTER_ROOT,MYID,COMM,
     &  RHS_SEQ,SIZE_ROOT,A,INFO,MTYPE,LDLT )
      IMPLICIT NONE
      INTEGER NRHS, MTYPE
      INTEGER DESCA_PAR( 9 )
      INTEGER LOCAL_M, LOCAL_N, MBLOCK, NBLOCK
      INTEGER CNTXT_PAR, MASTER_ROOT, SIZE_ROOT
      INTEGER MYID, COMM
      INTEGER LPIV, IPIV( LPIV )
      INTEGER INFO(80), LDLT
      DOUBLE PRECISION RHS_SEQ( SIZE_ROOT *NRHS)
      DOUBLE PRECISION A( LOCAL_M, LOCAL_N )
      INTEGER IERR, NPROW, NPCOL, MYROW, MYCOL
      INTEGER LOCAL_N_RHS
      DOUBLE PRECISION, ALLOCATABLE, DIMENSION( :,: ) ::RHS_PAR
      EXTERNAL numroc
      INTEGER  numroc
      INTEGER allocok
      CALL blacs_gridinfo( cntxt_par, nprow, npcol, myrow, mycol )
      local_n_rhs = numroc(nrhs, nblock, mycol, 0, npcol)
      local_n_rhs = max(1,local_n_rhs)
      ALLOCATE(rhs_par(local_m, local_n_rhs),stat=allocok)
      IF (allocok > 0 ) THEN
        WRITE(*,*) ' Problem during solve of the root.'
        WRITE(*,*) ' Reduce number of right hand sides.'
        CALL mumps_abort()
      ENDIF
      CALL dmumps_scatter_root( myid, size_root, nrhs, rhs_seq,
     &      local_m, local_n_rhs,
     &      mblock, nblock, rhs_par, master_root,
     &      nprow, npcol, comm )
      CALL dmumps_solve_2d_bcyclic (size_root, nrhs, mtype,
     &     a, desca_par, local_m, local_n, local_n_rhs,
     &     ipiv, lpiv, rhs_par, ldlt, 
     &     mblock, nblock, cntxt_par,
     &     ierr)
      CALL dmumps_gather_root( myid, size_root, nrhs,
     &    rhs_seq, local_m, local_n_rhs,
     &    mblock, nblock, rhs_par, master_root,
     &    nprow, npcol, comm )
      DEALLOCATE(rhs_par)
      RETURN
      END SUBROUTINE dmumps_root_solve
      SUBROUTINE dmumps_solve_2d_bcyclic (SIZE_ROOT, NRHS, MTYPE,
     &     A, DESCA_PAR, LOCAL_M, LOCAL_N, LOCAL_N_RHS,
     &     IPIV, LPIV, RHS_PAR, LDLT, 
     &     MBLOCK, NBLOCK, CNTXT_PAR,
     &     IERR)
      IMPLICIT NONE
      INTEGER, intent (in) :: SIZE_ROOT, NRHS, LDLT, LOCAL_M, 
     &                        local_n, local_n_rhs, 
     &                        mblock, nblock, cntxt_par, mtype
      INTEGER, intent (in) :: DESCA_PAR( 9 ) 
      INTEGER, intent (in) :: LPIV, IPIV( LPIV )
      DOUBLE PRECISION, intent (in) :: A( LOCAL_M, LOCAL_N )
      DOUBLE PRECISION, intent (inout) :: RHS_PAR(LOCAL_M, LOCAL_N_RHS)
      INTEGER, intent (out) :: IERR
      INTEGER              :: DESCB_PAR( 9 )
      IERR = 0
      CALL descinit( descb_par, size_root, 
     &      nrhs, mblock, nblock, 0, 0,
     &      cntxt_par, local_m, ierr )
            IF (ierr.NE.0) THEN
              WRITE(*,*) 'After DESCINIT, IERR = ', ierr
              CALL mumps_abort()
            END IF
      IF ( ldlt .eq. 0 .OR. ldlt .eq. 2 ) THEN
        IF ( mtype .eq. 1 ) THEN
          CALL pdgetrs('N',size_root,nrhs,a,1,1,desca_par,ipiv,
     &      rhs_par,1,1,descb_par,ierr)
        ELSE
          CALL pdgetrs('T',size_root,nrhs,a,1,1,desca_par,ipiv,
     &      rhs_par, 1, 1, descb_par,ierr)
        END IF
      ELSE
        CALL pdpotrs( 'L', size_root, nrhs, a, 1, 1, desca_par,
     &    rhs_par, 1, 1, descb_par, ierr )
      END IF
      IF ( ierr .LT. 0 ) THEN
        WRITE(*,*) ' Problem during solve of the root'
        CALL mumps_abort()
      END IF
      RETURN
      END SUBROUTINE dmumps_solve_2d_bcyclic