parini.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE parini
!                   *****************
!
     &(nhp,nhm,indpu,npoin2,nachb,nplan,mesh,nb_neighb,
     & nb_neighb_seg,nelem2,ifapar,modass)
!
!***********************************************************************
! BIEF   V7P1
!***********************************************************************
!
!brief    INITIALISES THE ARRAYS USED IN PARALLEL MODE.
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        13/07/2010
!+        V6P0
!+   Translation of French comments within the FORTRAN sources into
!+   English comments
!
!history  N.DURAND (HRW), S.E.BOURBAN (HRW)
!+        21/08/2010
!+        V6P0
!+   Creation of DOXYGEN tags for automated documentation and
!+   cross-referencing of the FORTRAN sources
!
!history  J-M HERVOUET (EDF R&D, LNHE)
!+        09/05/2014
!+        V7P0
!+   Adding an allocation of BUF_SEND%I and BUF_RECV%I
!+   Adding an allocation of BUF_SENDI8 and BUF_RECVI8
!+   for I4 and I8 integer communications.
!
!history  J-M HERVOUET (EDF R&D, LNHE)
!+        30/07/2015
!+        V7P1
!+   FAC suppressed.
!
!history  R.NHEILI (Univerte de Perpignan, DALI)
!+        24/02/2016
!+        V7P3
!+      ADD ALLOCATION BUF_SEND_ERR AND BUF_RECV_ERR
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| IFAPAR           |-->| IFAPAR(1:3,IELEM)=PROCESSOR NUMBERS BEHIND THE
!|                  |   | 3 ELEMENT EDGES  (NUMBERS FROM 0 TO NCSIZE-1)
!|                  |   | IFAPAR(4:6,IELEM): -LOCAL- ELEMENT NUMBERS
!|                  |   | BEHIND THE 3 EDGES
!| INDPU            |<--| INDEX TABLE : IF 0: NOT AN INTERFACE POINT
!|                  |   |               IF NOT 0: ADDRESS IN THE LIST
!|                  |   |               OF BOUNDARY POINTS.
!| MESH             |-->| MESH STRUCTURE
!| MODASS           |-->| ASSEMBLY MODE 1: NORMAL 2: WITH INTEGERS 3: COMPENSATION
!| NACHB            |-->| IF 'IL' IS THE LOCAL RANK OF A NEIGHBOURING
!|                  |   | SUB-DOMAIN AND 'IP' ONE INTERFACE POINT
!|                  |   | NACHB(IL,IP) WILL BE THE REAL NUMBER OF THIS
!|                  |   | NEIGHBOURING SUB-DOMAIN
!|                  |   | THE LIST IN NACHB IS ORDERED WITH THE
!|                  |   | GLOBAL NUMBERS OF POINTS (HENCE THE POINTS
!|                  |   | WILL BE FOUND IN THE SAME ORDER BY ALL
!|                  |   | PROCESSORS)
!| NB_NEIGHB        |<--| NUMBER OF NEIGHBOURING SUB-DOMAINS (FOR POINTS)
!| NB_NEIGHB_SEG    |<--| NUMBER OF NEIGHBOURING SUB-DOMAINS (FOR EDGES)
!| NB_NEIGHB_PT_SEG |<--| NUMBER OF SEGMENTS SHARED WITH A NEIGHBOUR
!| NELEM2           |-->| NUMBER OF ELEMENTS IN 2D
!| NHM              |<--| NODE NUMBERS OF PROCESSORS WITH SMALLER RANK
!| NHP              |<--| NODE NUMBERS OF PROCESSORS WITH LARGER RANK
!| NPLAN            |-->| NUMBER OF PLANES IN 3D
!| NPOIN2           |-->| NUMBER OF POINTS IN 2D
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief, ex_parini => parini
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN)            :: NPOIN2,NPLAN,NELEM2,MODASS
      INTEGER, INTENT(INOUT)         :: NB_NEIGHB,NB_NEIGHB_SEG
      INTEGER, INTENT(INOUT)         :: NHP(nbmaxdshare,nptir)
      INTEGER, INTENT(INOUT)         :: NHM(nbmaxdshare,nptir)
      INTEGER, INTENT(IN)            :: NACHB(nbmaxnshare,nptir)
      INTEGER, INTENT(IN)            :: IFAPAR(6,nelem2)
      INTEGER, INTENT(INOUT)         :: INDPU(npoin2)
      TYPE(bief_mesh), INTENT(INOUT) :: MESH
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER IKP(nbmaxdshare,2),IKM(nbmaxdshare,2)
      INTEGER I,J,IL,IZH,II,IMAX,IMIN,ILMAX,IELEM,IFACE
      INTEGER ILP,ILM,IPA,IKA,IPB,IKB,NB_PT_MX,DIM1HCOM,CHECKSUM
      LOGICAL NEW
!
!-----------------------------------------------------------------------
!
!     INITIALISES THE PROCESSOR NUMBERS FOR 2D MESSAGE-PASSING
!
      DO i=1,nbmaxdshare
        ikp(i,1)=-1
        ikm(i,1)=-1
        ikp(i,2)=0
        ikm(i,2)=0
      ENDDO
!
!     PREPARES COMMUNICATION
!     IN THE FOLLOWING SEQUENCE :
!     1) SENDS     TO   PROCESSORS WITH NUMBER IPID + IL
!     2) RECEIVES  FROM PROCESSORS WITH NUMBER IPID - IL
!     3) SENDS     TO   PROCESSORS WITH NUMBER IPID - IL
!     4) RECEIVES  FROM PROCESSORS WITH NUMBER IPID + IL
!
!     LEVEL IL : SENDS AND RECEIVES
!
!
!     SENDS TO PROCESSORS WITH NUMBER GREATER THAN IPID
!
      imax=ipid
!
      IF (ipid.NE.ncsize-1) THEN
        izh=1
        DO il=ipid+1,ncsize-1
          ii=0
          DO i=1,nptir
            DO j=2,nbmaxnshare
              IF(nachb(j,i).EQ.il) THEN
                IF(izh.GT.nbmaxdshare) THEN
                  WRITE(lu,*) 'PARINI: NBMAXDSHARE TOO SMALL'
                  CALL plante(1)
                  stop
                ENDIF
                ii=ii+1
                nhp(izh,ii)=nachb(1,i)
              ENDIF
            ENDDO ! J
          ENDDO ! I
          IF(ii.NE.0) THEN
            ikp(izh,1)=il
            ikp(izh,2)=ii
            izh=izh+1
            imax=il
          ENDIF
        ENDDO ! IL
      ENDIF
!
!
!     RECEIVES FROM PROCESSORS WITH NUMBER LOWER THAN IPID
!
      imin=ipid
!
      IF (ipid.NE.0) THEN
        izh=1
        DO il=ipid-1,0,-1
          ii=0
          DO i=1,nptir
            DO j=2,nbmaxnshare
              IF(nachb(j,i).EQ.il) THEN
                IF(izh.GT.nbmaxdshare) THEN
                  WRITE(lu,*) 'PARINI: NBMAXDSHARE TOO SMALL'
                  CALL plante(1)
                  stop
                ENDIF
                ii=ii+1
                nhm(izh,ii)=nachb(1,i)
              ENDIF
            ENDDO ! J
          ENDDO ! I
          IF(ii.NE.0) THEN
            ikm(izh,1)=il
            ikm(izh,2)=ii
            izh=izh+1
            imin=il
          ENDIF
        ENDDO ! IL
      ENDIF
!
!**   DETERMINES ILMAX
!
      ilmax=max(imax-ipid,ipid-imin)
!
!-----------------------------------------------------------------------
!
!====   COMPUTES THE NUMBER OF NEIGHBOURS
!
        nb_pt_mx  = 0
        nb_neighb = 0
        ilp = 1
        ilm = 1
!**     PROCESSOR OF HIGHER RANK
        DO il=1,ilmax
          ipa=ikp(ilp,1)
          ika=ikp(ilp,2)
          IF(ipa.EQ.ipid+il.AND.ika.NE.0) THEN
            nb_neighb = nb_neighb + 1
            IF(ika.GT.nb_pt_mx) nb_pt_mx=ika
          ENDIF
          IF(ipa.EQ.ipid+il) ilp=ilp+1
        ENDDO
!**      PROCESSOR OF LOWER RANK
        DO il=1,ilmax
          ipb=ikm(ilm,1)
          ikb=ikm(ilm,2)
          IF(ipb.EQ.ipid-il.AND.ikb.NE.0) THEN
            nb_neighb = nb_neighb + 1
            IF(ikb.GT.nb_pt_mx) nb_pt_mx=ikb
          ENDIF
          IF(ipb.EQ.ipid-il) ilm=ilm+1
        ENDDO
!
!====   ENDS COMPUTATION OF THE NUMBER OF NEIGHBOURS
!
        CALL bief_allvec(2,mesh%NB_NEIGHB_PT,'NBNGPT',
     &                   nb_neighb,1,0,mesh)
        CALL bief_allvec(2,mesh%LIST_SEND   ,'LSSEND',
     &                   nb_neighb,1,0,mesh)
!
! ALIGNMENT ON 16 BYTES
!
        dim1hcom = nb_pt_mx/4
        IF(mod(nb_pt_mx,4).EQ.0) THEN
          dim1hcom = dim1hcom*4
        ELSE
          dim1hcom = dim1hcom*4 + 4
        ENDIF
        CALL bief_allvec(2,mesh%NH_COM,'NH_COM',
     &                   dim1hcom,nb_neighb,0,mesh)
!
!====   COMPUTES THE NUMBER OF INTERFACE POINTS PER NEIGHBOUR
!
        nb_neighb = 0
        ilp = 1
        ilm = 1
        DO il=1,ilmax
          ipa=ikp(ilp,1)
          ika=ikp(ilp,2)
          IF(ipa.EQ.ipid+il.AND.ika.NE.0) THEN
            nb_neighb = nb_neighb + 1
            mesh%NB_NEIGHB_PT%I(nb_neighb) = ika
            mesh%LIST_SEND%I(nb_neighb) = ipa
            DO i=1,ika
              mesh%NH_COM%I(dim1hcom*(nb_neighb-1)+i)=nhp(ilp,i)
            ENDDO
          ENDIF
          IF(ipa.EQ.ipid+il) ilp=ilp+1
        ENDDO
        DO il=1,ilmax
          ipb=ikm(ilm,1)
          ikb=ikm(ilm,2)
          IF(ipb.EQ.ipid-il.AND.ikb.NE.0) THEN
            nb_neighb = nb_neighb + 1
            mesh%NB_NEIGHB_PT%I(nb_neighb) = ikb
            mesh%LIST_SEND%I(nb_neighb) = ipb
            DO i=1,ikb
              mesh%NH_COM%I(dim1hcom*(nb_neighb-1)+i)=nhm(ilm,i)
            ENDDO
          ENDIF
          IF(ipb.EQ.ipid-il) ilm=ilm+1
        ENDDO
!
!==== ENDS COMPUTATION OF THE NUMBER OF INTERFACE POINTS PER NEIGHBOUR
!
!=== POSSIBILITY OF SORTING LIST_SEND AND RECV FOR TORE BG
!
! ALIGNMENT ON 16BYTES BOUNDARIES
!
        nb_pt_mx = nb_pt_mx * nplan
        il = nb_pt_mx/2
        IF(mod(nb_pt_mx,2).EQ.0) THEN
          il = il*2
        ELSE
          il = il*2 + 2
        ENDIF
        CALL bief_allvec(1,mesh%BUF_SEND,'BUSEND',il*3,nb_neighb,0,mesh)
        CALL bief_allvec(1,mesh%BUF_RECV,'BURECV',il*3,nb_neighb,0,mesh)
        IF (modass .EQ.3) THEN
          CALL bief_allvec(1,mesh%BUF_SEND_ERR,'BSSERR',
     &                     il*3,nb_neighb,0,mesh)
          CALL bief_allvec(1,mesh%BUF_RECV_ERR,'BSRERR',
     &                     il*3,nb_neighb,0,mesh)
        ENDIF
!
!       ADDED FOR INTEGER I4 COMMUNICATIONS
!
        ALLOCATE(mesh%BUF_SEND%I(il*3*nb_neighb))
        ALLOCATE(mesh%BUF_RECV%I(il*3*nb_neighb))
!
!       ADDED FOR INTEGER I8 COMMUNICATIONS
!
        IF(modass.EQ.2) THEN
          ALLOCATE(mesh%BUF_SENDI8(il*3*nb_neighb))
          ALLOCATE(mesh%BUF_RECVI8(il*3*nb_neighb))
        ENDIF
!
!-----------------------------------------------------------------------
!
!     FOR SEGMENTS
!
!     WE ASSUME HERE THAT NB_NEIGHB.GE.NB_NEIGHB_SEG
!
!     NOTE: NH_COM_SEG IS FILLED WITH 4*IELEM+IFACE
!           THIS IS TO RETRIEVE IELEM AND IFACE ONCE ELTSEG IS KNOWN
!           THE FINAL VALUE OF NH_COM_SEG IS ELTSEG(IELEM,IFACE)
!
      CALL bief_allvec(2,mesh%NB_NEIGHB_PT_SEG,'NBNGSG',
     &                 nb_neighb,1,0,mesh)
      CALL bief_allvec(2,mesh%LIST_SEND_SEG   ,'LSSESG',
     &                 nb_neighb,1,0,mesh)
      CALL bief_allvec(2,mesh%NH_COM_SEG      ,'NH_CSG',
     &                 dim1hcom,nb_neighb,0,mesh)
!
      nb_neighb_seg=0
!
!     INITIALISES NH_COM_SEG (SEE COMP_NH_COM_SEG)
!
      DO i=1,dim1hcom*nb_neighb
        mesh%NH_COM_SEG%I(i)=-999999
      ENDDO
!
      DO ielem=1,nelem2
!
!       LOOKS FOR A FACE WITH THE OTHER SIDE IN ANOTHER SUB-DOMAIN
!
!       ELEMENTS WITHOUT ANY INTERFACE SEGMENT HAVE 3 ZEROS
        checksum=ifapar(1,ielem)**2+
     &           ifapar(2,ielem)**2+
     &           ifapar(3,ielem)**2
!
        IF(checksum.NE.0) THEN
        DO iface=1,3
!
          ilm=ifapar(iface,ielem)
          IF(ilm.GE.0.AND.ilm.NE.ipid) THEN
!           NEW INTERFACE SEGMENT FOUND
            IF(nb_neighb_seg.EQ.0) THEN
!             THE FIRST ONE
              nb_neighb_seg=1
              mesh%NB_NEIGHB_PT_SEG%I(1)=1
              mesh%LIST_SEND_SEG%I(1)=ilm
              mesh%NH_COM_SEG%I(1)=4*ielem+iface
            ELSE
!             FROM THE SECOND ON
!             IS IT A NEW PROCESSOR
              new=.true.
              DO il=1,nb_neighb_seg
                IF(ilm.EQ.mesh%LIST_SEND_SEG%I(il)) THEN
!                 NEW SEGMENT, OLD PROCESSOR
                  mesh%NB_NEIGHB_PT_SEG%I(il)=
     &            mesh%NB_NEIGHB_PT_SEG%I(il)+1
                  i=mesh%NB_NEIGHB_PT_SEG%I(il)
                  mesh%NH_COM_SEG%I(dim1hcom*(il-1)+i)=4*ielem+iface
                  new=.false.
                  EXIT
                ENDIF
              ENDDO
              IF(new) THEN
!               NEW SEGMENT, NEW PROCESSOR
                nb_neighb_seg=nb_neighb_seg+1
                mesh%NB_NEIGHB_PT_SEG%I(nb_neighb_seg)=1
                mesh%LIST_SEND_SEG%I(nb_neighb_seg)=ilm
                mesh%NH_COM_SEG%I(dim1hcom*(nb_neighb_seg-1)+1)=
     &                            4*ielem+iface
              ENDIF
            ENDIF
          ENDIF
!
        ENDDO
        ENDIF
!
      ENDDO
!
      IF(nb_neighb_seg.GT.nb_neighb) THEN
        WRITE(lu,*) 'IN PARINI NB_NEIGHB    =',nb_neighb
        WRITE(lu,*) '          NB_NEIGHB_SEG=',nb_neighb_seg
        CALL plante(1)
        stop
      ENDIF
!
!-----------------------------------------------------------------------
!
!     INDEX TABLE FOR BUFFER IN COMMUNICATION
!
      DO i=1,npoin2
        indpu(i)=0
      ENDDO
!
!  COEFFICIENTS FOR THE SCALAR PRODUCT:
!
      IF(nptir.GT.0) THEN
        DO i=1,nptir
          indpu(nachb(1,i))=i
        ENDDO
      ENDIF
!
!-----------------------------------------------------------------------
!
      RETURN
      END