1C
2C  This file is part of MUMPS 5.1.2, released
3C  on Mon Oct  2 07:37:01 UTC 2017
4C
5C
6C  Copyright 1991-2017 CERFACS, CNRS, ENS Lyon, INP Toulouse, Inria,
7C  University of Bordeaux.
8C
9C  This version of MUMPS is provided to you free of charge. It is
10C  released under the CeCILL-C license:
11C  http://www.cecill.info/licences/Licence_CeCILL-C_V1-en.html
12C
13      SUBROUTINE DMUMPS_ROOT_SOLVE( NRHS, DESCA_PAR,
14     &  CNTXT_PAR,LOCAL_M,LOCAL_N,MBLOCK,NBLOCK,
15     &  IPIV,LPIV,MASTER_ROOT,MYID,COMM,
16     &  RHS_SEQ,SIZE_ROOT,A,INFO,MTYPE,LDLT )
17      IMPLICIT NONE
18      INTEGER NRHS, MTYPE
19      INTEGER DESCA_PAR( 9 )
20      INTEGER LOCAL_M, LOCAL_N, MBLOCK, NBLOCK
21      INTEGER CNTXT_PAR, MASTER_ROOT, SIZE_ROOT
22      INTEGER MYID, COMM
23      INTEGER LPIV, IPIV( LPIV )
24      INTEGER INFO(40), LDLT
25      DOUBLE PRECISION RHS_SEQ( SIZE_ROOT *NRHS)
26      DOUBLE PRECISION A( LOCAL_M, LOCAL_N )
27      INTEGER IERR, NPROW, NPCOL, MYROW, MYCOL
28      INTEGER LOCAL_N_RHS
29      DOUBLE PRECISION, ALLOCATABLE, DIMENSION( :,: ) ::RHS_PAR
30      EXTERNAL numroc
31      INTEGER  numroc
32      INTEGER allocok
33      CALL blacs_gridinfo( CNTXT_PAR, NPROW, NPCOL, MYROW, MYCOL )
34      LOCAL_N_RHS = numroc(NRHS, NBLOCK, MYCOL, 0, NPCOL)
35      LOCAL_N_RHS = max(1,LOCAL_N_RHS)
36      ALLOCATE(RHS_PAR(LOCAL_M, LOCAL_N_RHS),stat=allocok)
37      IF (allocok > 0 ) THEN
38        WRITE(*,*) ' Problem during solve of the root.'
39        WRITE(*,*) ' Reduce number of right hand sides.'
40        CALL MUMPS_ABORT()
41      ENDIF
42      CALL DMUMPS_SCATTER_ROOT( MYID, SIZE_ROOT, NRHS, RHS_SEQ,
43     &      LOCAL_M, LOCAL_N_RHS,
44     &      MBLOCK, NBLOCK, RHS_PAR, MASTER_ROOT,
45     &      NPROW, NPCOL, COMM )
46      CALL DMUMPS_SOLVE_2D_BCYCLIC (SIZE_ROOT, NRHS, MTYPE,
47     &     A, DESCA_PAR, LOCAL_M, LOCAL_N, LOCAL_N_RHS,
48     &     IPIV, LPIV, RHS_PAR, LDLT,
49     &     MBLOCK, NBLOCK, CNTXT_PAR,
50     &     IERR)
51      CALL DMUMPS_GATHER_ROOT( MYID, SIZE_ROOT, NRHS,
52     &    RHS_SEQ, LOCAL_M, LOCAL_N_RHS,
53     &    MBLOCK, NBLOCK, RHS_PAR, MASTER_ROOT,
54     &    NPROW, NPCOL, COMM )
55      DEALLOCATE(RHS_PAR)
56      RETURN
57      END SUBROUTINE DMUMPS_ROOT_SOLVE
58      SUBROUTINE DMUMPS_SOLVE_2D_BCYCLIC (SIZE_ROOT, NRHS, MTYPE,
59     &     A, DESCA_PAR, LOCAL_M, LOCAL_N, LOCAL_N_RHS,
60     &     IPIV, LPIV, RHS_PAR, LDLT,
61     &     MBLOCK, NBLOCK, CNTXT_PAR,
62     &     IERR)
63      IMPLICIT NONE
64      INTEGER, intent (in) :: SIZE_ROOT, NRHS, LDLT, LOCAL_M,
65     &                        LOCAL_N, LOCAL_N_RHS,
66     &                        MBLOCK, NBLOCK, CNTXT_PAR, MTYPE
67      INTEGER, intent (in) :: DESCA_PAR( 9 )
68      INTEGER, intent (in) :: LPIV, IPIV( LPIV )
69      DOUBLE PRECISION, intent (in) :: A( LOCAL_M, LOCAL_N )
70      DOUBLE PRECISION, intent (inout) :: RHS_PAR(LOCAL_M, LOCAL_N_RHS)
71      INTEGER, intent (out) :: IERR
72      INTEGER              :: DESCB_PAR( 9 )
73      IERR = 0
74      CALL DESCINIT( DESCB_PAR, SIZE_ROOT,
75     &      NRHS, MBLOCK, NBLOCK, 0, 0,
76     &      CNTXT_PAR, LOCAL_M, IERR )
77            IF (IERR.NE.0) THEN
78              WRITE(*,*) 'After DESCINIT, IERR = ', IERR
79              CALL MUMPS_ABORT()
80            END IF
81      IF ( LDLT .eq. 0 .OR. LDLT .eq. 2 ) THEN
82        IF ( MTYPE .eq. 1 ) THEN
83          CALL pdgetrs('N',SIZE_ROOT,NRHS,A,1,1,DESCA_PAR,IPIV,
84     &      RHS_PAR,1,1,DESCB_PAR,IERR)
85        ELSE
86          CALL pdgetrs('T',SIZE_ROOT,NRHS,A,1,1,DESCA_PAR,IPIV,
87     &      RHS_PAR, 1, 1, DESCB_PAR,IERR)
88        END IF
89      ELSE
90        CALL pdpotrs( 'L', SIZE_ROOT, NRHS, A, 1, 1, DESCA_PAR,
91     &    RHS_PAR, 1, 1, DESCB_PAR, IERR )
92      END IF
93      IF ( IERR .LT. 0 ) THEN
94        WRITE(*,*) ' Problem during solve of the root'
95        CALL MUMPS_ABORT()
96      END IF
97      RETURN
98      END SUBROUTINE DMUMPS_SOLVE_2D_BCYCLIC
99