partitioner.F

Go to the documentation of this file.

!                   **************
                    SUBROUTINE partitioner
!                   **************
     &    (pmethod, nelem, npoin, ndp, nparts, ikles, epart, npart)
!
!
!***********************************************************************
! PARALLEL   V6P2                                   21/08/2010
!***********************************************************************
!
!brief    call to the partionning software
!
!history   R. KOPMANN (BAW)
!+
!+
!+         created
!
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| PMETHOD        |<--| 1: FOR METIS 2: FOR SCOTCH
!| NELEM          |<--| THE NUMBER OF ELEMENTS
!| NDP            |<--| THE NUMBE OF POINT PER ELEMENT
!| NPARTS         |<--| NUMBER OF PARTITIONS
!| IKLES          |<--| CONNECTIVITY TABLE
!| EPART          |-->| PARTITION NUMBER OF AN ELEMENT
!| MYPART         |-->| PARTITION NUMBER OF A POINT
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE iso_c_binding
!
      USE declarations_special
      USE c_binding
      IMPLICIT NONE
#if defined HAVE_SCOTCH
      include 'scotchf.h'
#endif
!
      INTEGER, INTENT(IN)  :: pmethod
      INTEGER, INTENT(IN)  :: nelem
      INTEGER, INTENT(IN)  :: npoin
      INTEGER, INTENT(IN)  :: ndp
      INTEGER, INTENT(IN)  :: ikles(nelem*ndp)
      INTEGER, INTENT(IN)  :: nparts
      INTEGER, INTENT(OUT) :: epart(nelem)
      INTEGER, INTENT(OUT) :: npart(npoin)
!
!-----------------------------------------------------------------------
!
!
!
      INTEGER i,j,k,err,numflag,edgecut,ncommonnodes
      INTEGER, ALLOCATABLE :: eptr(:), eind(:)
!
! SCOTCH
!
#if defined HAVE_SCOTCH
      REAL*8, DIMENSION(SCOTCH_GRAPHDIM) :: scotchgraph
      REAL*8, DIMENSION(SCOTCH_STRATDIM) :: scotchstrat
      TYPE(c_ptr), POINTER :: ptxadj(:), ptadjncy(:)
      INTEGER, ALLOCATABLE :: recvbuf(:)
      INTEGER, POINTER :: xadj2(:), adjncy2(:)
#endif

      IF(pmethod.EQ.1) THEN
        WRITE(lu,*) 'BEGIN PARTITIONING WITH METIS'
!       NEW METIS INTERFACE (>= VERSION 5) :
!
!       EPTR, EIND: THESE ARRAYS SPECIFIES THE ELEMENTS
!       THAT ARE STORED LOCALLY AT EACH PROCESSOR.
!       CF. DISCUSSION IN SECTION 4.3 (?)
!
        ALLOCATE(eptr(nelem+1),stat=err)
        CALL check_allocate(err,'EPTR')
        ALLOCATE(eind(nelem*ndp),stat=err)
        CALL check_allocate(err,'EIND')
!
!
        DO i=1,nelem+1
          eptr(i) = (i-1)*ndp + 1
        ENDDO
!
        k=1
        DO i=1,nelem
          DO j=1,ndp
            eind(k) = ikles((i-1)*ndp+j)
            k = k + 1
          ENDDO
        ENDDO
!
!       END OF CORRECTION JMH 04/07/2012
!
!       SWITCH TO C NUMBERING
        eind = eind -1
        eptr = eptr -1
!
!       METIS REQUIRES THE NUMBER OF COMMON POINT NEEDED BETWEEN 2 ELEMENTS TO MAKE AN EDGE
!       NCOMMONNODES = 2   FOR TRIANGLE OR RECTANGLE
!       NCOMMONNODES = 3   FOR TETRAHEDRE
!       NCOMMONNODES = 4   FOR HEXAHEDRE

!
        IF (ndp==3.OR.ndp==6) THEN
          ncommonnodes = 2 ! FOR TRIANGLE OR PRISM
        ELSE IF(ndp==4) THEN
          ncommonnodes = 3 ! FORT TETRAHEDRON
        ELSE
          WRITE(lu,*) 'METIS: IMPLEMENTED FOR TRIANGLES OR PRISMS ONLY'
          CALL plante(1)
        ENDIF
!
!       WE ONLY USE METIS_PARTMESHDUAL AS ONLY THE FINITE ELEMENTS PARTITION
!       IS RELEVANT HERE.
!
!       IMPORTANT: WE USE FORTRAN-LIKE FIELD ELEMENTS NUMBERING 1...N
!       IN C VERSION, 0...N-1 NUMBERING IS APPLIED!!!
!
        numflag = 1
!
#if defined HAVE_MPI
        CALL mymetis_partmeshdual
     &        (nelem, npoin, eptr, eind,
     &          ncommonnodes, nparts,
     &          edgecut, epart, npart)
#else
        WRITE(lu,*) 'ERROR: TRY TO RUN PARTEL WITH A '//
     &              'SERIAL CONFIGURATION'
        CALL plante(1)
#endif
!
!       DEALLOCATING TEMPORARY ARRAYS FOR METIS
!
!       EPART IS AN ARRAY
        epart = epart+1
        DEALLOCATE(eptr)
        DEALLOCATE(eind)
      ELSE IF(pmethod.EQ.2) THEN
#if defined HAVE_SCOTCH
!
!       SCOTCH PARTITIONER
!
        WRITE(lu,*) 'BEGIN PARTITIONING WITH SCOTCH'
        ALLOCATE (eptr(nelem+1),stat=err)
        CALL check_allocate(err,"EPTR")
        ALLOCATE (eind(nelem*ndp),stat=err)
        CALL check_allocate(err,"EIND")
!
!       bUILDING EPTR AND EIND SAME AS BEFORE.
!
        DO i=1,nelem+1
          eptr(i) = (i-1)*ndp + 1
        ENDDO
!
        k=1
        DO i=1,nelem
          DO j=1,ndp
            eind(k) = ikles((i-1)*ndp+j)
            k = k + 1
          ENDDO
        ENDDO
!       THE NUMBER OF COMMON POINT NEEDED BETWEEN 2 ELEMENTS TO MAKE AN
!       EDGE
!       NCOMMONNODES = 2 FOR TRIANGLE OR RECTANGLE
!       NCOMMONNODES = 3 FOR TETRAHEDRE
!       NCOMMONNODES = 4 FOR HEXAHEDRE
        IF (ndp==3.OR.ndp==6) THEN
          ncommonnodes = 2 ! FOR TRIANGLE OR RECTANGLE
        ELSE
          WRITE(lu,*)'SCOTCH: IMPLEMENTED FOR TRIANGLES OR PRISMS ONLY'
          CALL plante(1)
        ENDIF

!
!       mETIS_MESHTODUAL USES POINTER OF POINTER (PTXADJ, PTADJNCY) SO IN ORDER TO USE THAT
!       IN FORTRAN WE USE THE ISO_BINDING MODULE
!
        ALLOCATE(ptxadj(1))
        ALLOCATE(ptadjncy(1))
!
!       BUILD THE DUAL GRAPH BY CALLING THE METIS SUBROUTINE
!
        numflag = 1
        CALL metis_meshtodual
     &       (nelem, npoin, eptr, eind, ncommonnodes, numflag,
     &        ptxadj, ptadjncy)

!
!       TRANSFORM THE C POINTER OF POINTER BACK TO A FORTRAN ARRAY
!
        ALLOCATE(recvbuf(1))
        recvbuf(1) = nelem+1
        CALL c_f_pointer(ptxadj(1),xadj2,recvbuf)
!
        recvbuf(1) = xadj2(nelem+1)-1
        CALL c_f_pointer(ptadjncy(1),adjncy2,recvbuf)
        DEALLOCATE(recvbuf)
!
!       INITIALIZE THE GRAPH AND THE STRATEGY USED BY SCOTCH FOR THE
!       PARTITIONING WE USE THE DEFAULT STRATEGY
!
        CALL scotchfstratinit(scotchstrat, err)
        IF (err.NE.0) THEN
          WRITE(lu,*) 'SCOTCH ERROR: CANNOT INITIALIZE STRAT'
          CALL plante(1)
        ENDIF
!
        CALL scotchfgraphinit(scotchgraph, err)
        IF (err.NE.0) THEN
          WRITE(lu,*) 'SCOTCH ERROR: CANNOT INITIALIZE GRAPH'
          CALL plante(1)
        ENDIF
!
        CALL scotchfgraphbuild ( scotchgraph,         ! GRAFDAT
     &                           numflag,             ! BASEVAL
     &                           nelem,              ! VERTNBR
     &                           xadj2(1:nelem),     ! VERTTAB
     &                           xadj2(2:nelem+1),   ! VENDTAB
     &                           xadj2(1:nelem),     ! VELOTAB, VERTEX WEIGHTS (nULL)
     &                           xadj2(1:nelem),     ! VLBLTAB, VERTEX LABELS (nULL)
     &                           xadj2(nelem+1)-1,   ! EDGENBR
     &                           adjncy2,             ! EDGETAB
     &                           adjncy2,             ! EDLOTAB (nULL)
     &                           err )
!
        IF (err.NE.0) THEN
          WRITE(lu,*) 'SCOTCH ERROR: CANNOT BUILD GRAPH'
          CALL plante(1)
        ENDIF
!
        CALL scotchfgraphcheck(scotchgraph,err)
        IF (err.NE.0) THEN
          WRITE(lu,*) 'SCOTCH ERROR: GRAPH NOT CONSISTANT'
          CALL plante(1)
        ENDIF
!
        ! RUN PARTITIONING
        CALL scotchfgraphpart ( scotchgraph,
     &                          nparts,
     &                          scotchstrat,
     &                          epart,
     &                          err )
!
        IF (err.NE.0) THEN
          WRITE(lu,*) 'SCOTCH ERROR: CANNOT PARTITION GRAPH'
          CALL plante(1)
        END IF
!!!! CHANGING EPART NUMBERING TO 1-NPART
        epart = epart + 1

        CALL scotchfgraphexit (scotchgraph)
        CALL scotchfstratexit (scotchstrat)

!
!        THOSE DEALLOCATION CAN GENERATE ERROR WITH SOME COMPILER
!        DEALLOCATE(XADJ2)
!        DEALLOCATE(ADJNCY2)
        DEALLOCATE(ptxadj)
        DEALLOCATE(ptadjncy)
        DEALLOCATE(eind)
        DEALLOCATE(eptr)
#else
        WRITE(lu,*) "TRYING TO USE SCOTCH TO PARTIONNE WHEN SCOTCH",
     &              " IS NOT INSTALLED"
        CALL plante(1)
#endif
        ENDIF
      END SUBROUTINE partitioner