eleb3dt.f

Go to the documentation of this file.

!                   ******************
                    SUBROUTINE eleb3dt
!                   ******************
!
     &(ikle3,nbor,nelbor,nelbor2d,iklbor,neleb,nelebx,nulone,
     & nelem2,npoin2,nplan,netage,nptfr,iklbor2d,neleb2d,nelebx2d)
!
!***********************************************************************
! BIEF   V7P0                                   19/03/2014
!***********************************************************************
!
!brief    CASE OF PRISMS SPLIT IN TETRAHEDRONS.
!+                BUILDS THE 3D MESH.
!+
!+            INPUT: 3D MESH ARRAYS FILLED BY A PRELIMINARY CALL
!+                       TO ELEBD.
!+
!+            OUTPUT: ARRAYS COMPLETE IN 3D.
!
!history  J-M HERVOUET(LNH)
!+        23/08/99
!+        V5P3
!+
!
!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(LNHE)
!+        12/08/2011
!+        V6P2
!+   Arguments NELEB and NELEBX added.
!
!history  J-M HERVOUET(EDF LAB, LNHE)
!+        27/01/2014
!+        V7P0
!+   In parallel, boundary elements that correspond to another subdomain
!+   are given a IKLBOR that reduces the element to a point.
!+   The previous version had IKLBOR = 0 which causes a crash when
!+   checking bounds.
!
!history  J-M HERVOUET (EDF LAB, LNHE)
!+        19/03/2014
!+        V7P0
!+   Boundary segments have now their own numbering, independent of
!+   boundary points numbering.
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| IKLBOR         |<--| CONNECTIVITY TABLE FOR BOUNDARY ELEMENTS
!|                |   | HERE DECLARED AS IKLBOR(NELEBX,3)
!|                |   | BUT FILLED AS IKLBOR(NELEB2D,2,NETAGE,3)
!|                |   | 2 IS THE NUMBER OF TETRAHEDRONS FORMING A
!|                |   | VERTICAL RECTANGLE IN THE BORDER.
!|                |   | 3 IS THE NUMBER OF POINTS IN EVERY TRIANGLE
!|                |   | FORMING THIS RECTANGLE
!| IKLBOR2D       |<--| CONNECTIVITY TABLE FOR BOUNDARY ELEMENTS IN 2D
!| IKLE3          |<--| CONNECTIVITY TABLE IN 3D, FOR TETRAHEDRONS
!|                |   | HERE DECLARED AS IKLE3(NELEM2,3,NETAGE,4)
!|                |   | 3 IS THE NUMBER OF TETRAHEDRONS PER PRISM
!|                |   | (TETRAHEDRONS ARE NUMBERED ACCORDINGLY)
!|                |   | 4 IS THE NUMBER OF POINTS IN A TETRAHEDRA.
!| NBOR           |-->| GLOBAL NUMBER OF BOUNDARY POINTS IN 2D
!| NELBOR         |-->| FOR THE KTH BOUNDARY EDGE, GIVES THE CORRESPONDING
!|                |   | ELEMENT (FROM MESH3D).
!| NELBOR2D       |-->| FOR THE KTH BOUNDARY EDGE, GIVES THE CORRESPONDING
!|                |   | ELEMENT (FROM MESH2D).
!| NELEB          |-->| NUMBER OF BOUNDARY ELEMENTS
!| NELEB2D        |-->| NUMBER OF BOUNDARY ELEMENTS OF 2D MESH
!| NELEBX         |-->| MAXIMUM NUMBER OF BOUNDARY ELEMENTS
!|                |   | USED AS FIRST DIMENSION OF IKLBOR
!| NELEBX2D       |-->| MAXIMUM NUMBER OF BOUNDARY ELEMENTS IN 2D MESH
!| NELEM2         |-->| NUMBER OF ELEMENTS IN 2D
!| NETAGE         |-->| NUMBER OF PLANES - 1
!| NPLAN          |-->| NUMBER OF PLANES
!| NPOIN2         |-->| NUMBER OF POINTS IN 2D
!| NPTFR          |-->| NUMBER OF BOUNDARY POINTS
!| NULONE         |-->| GOES WITH ARRAY NELBOR. NELBOR GIVES THE
!|                |   | ADJACENT ELEMENT, NULONE GIVES THE LOCAL
!|                |   | NUMBER OF THE FIRST NODE OF THE BOUNDARY EDGE
!|                |   | I.E. 1, 2 OR 3 FOR TRIANGLES.
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief, ex_eleb3dt => eleb3dt
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN)    :: NELEM2,NPOIN2,NPLAN,NETAGE,NPTFR
      INTEGER, INTENT(IN)    :: NELEBX,NELEB2D,NELEBX2D
      INTEGER, INTENT(INOUT) :: NELEB
      INTEGER, INTENT(INOUT) :: IKLE3(nelem2,3,netage,4)
      INTEGER, INTENT(INOUT) :: IKLBOR(nelebx,3)
      INTEGER, INTENT(IN)    :: IKLBOR2D(nelebx2d,2),NELBOR2D(nelebx2d)
      INTEGER, INTENT(INOUT) :: NULONE(nelebx,3),NELBOR(nelebx)
      INTEGER, INTENT(INOUT) :: NBOR(nptfr*nplan)
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      LOGICAL OK(2)
!
      INTEGER IELEM,IPOIN,T(3),IELEB,IELEB3,IPTFR2
      INTEGER IETAGE,IPTFR,IL1,IL2,IL3,IL4,IG(2,2,3),IL(2,2,3),IPLAN
      INTEGER :: IG1,IG2,IG3,IG4,NUM1(12),NUM2(12),NUM3(12),K,L,M,N
!
      parameter( num1 = (/
     &    1 , 2 , 4 , 1 , 3 , 2 , 2 , 3 , 4 , 3 , 1 , 4 /) )
      parameter( num2 = (/
     &    2 , 4 , 1 , 3 , 2 , 1 , 3 , 4 , 2 , 1 , 4 , 3 /) )
      parameter( num3 = (/
     &    4 , 1 , 2 , 2 , 1 , 3 , 4 , 2 , 3 , 4 , 3 , 1 /) )
!
!***********************************************************************
!
! CONNECTIVITY TABLES FOR BOUNDARY FACES --> IKLBOR , NBOR3 ,
! CORRESPONDENCE BETWEEN LOCAL BOUNDARY NUMBERS AND 3D LOCAL NUMBERS --> NULONE
!
!     COMPLETING NBOR
!
      DO iptfr = 1,nptfr
        ipoin = nbor(iptfr)
        DO iplan = 2,nplan
          nbor(iptfr +(iplan-1)*nptfr)=ipoin+(iplan-1)*npoin2
        ENDDO
      ENDDO
!
!     LATERAL BOUNDARIES :
!     FOR EACH RECTANGULAR FACE SPLIT IN TWO TRIANGLES
!     THE LOWER TRIANGLE IS NUMBER 1, THE HIGHER IS NUMBER 2
!
      DO ieleb = 1,neleb2d
!
        iptfr =iklbor2d(ieleb,1)
        iptfr2=iklbor2d(ieleb,2)
!
!       TRIANGLE TOUCHING THE BOUNDARY IN 2D
        ielem = nelbor2d(ieleb)
!
        DO ietage = 1,netage
!
!         3D BOUNDARY NUMBERING OF THE 4 POINTS OF THE RECTANGULAR FACE
!
          il1 = iptfr  + (ietage-1)*nptfr
          il2 = iptfr2 + (ietage-1)*nptfr
          il3 = il2 + nptfr
          il4 = il1 + nptfr
!
!         3D GLOBAL NUMBERING OF THE 4 POINTS OF THE RECTANGULAR FACE
!
          ig1 = nbor(iptfr)  + (ietage-1)*npoin2
          ig2 = nbor(iptfr2) + (ietage-1)*npoin2
          ig3 = ig2 + npoin2
          ig4 = ig1 + npoin2
!
!         NUMBERS OF THE 3 TETRAHEDRONS POSSIBLY TOUCHING THE FACE
!
          t(1) = (ietage-1)*3*nelem2+ielem
          t(2) = t(1) + nelem2
          t(3) = t(2) + nelem2
!
!         LOOKS FOR THE LOWER TRIANGLE (CAN BE 1-2-4 OR 1-2-3)
!
!         2 POSSIBLE FORMS OF THE LOWER TRIANGLE (GLOBAL AND BOUNDARY)
          ig(1,1,1)=ig1
          ig(1,1,2)=ig2
          ig(1,1,3)=ig4
          ig(1,2,1)=ig1
          ig(1,2,2)=ig2
          ig(1,2,3)=ig3
          il(1,1,1)=il1
          il(1,1,2)=il2
          il(1,1,3)=il4
          il(1,2,1)=il1
          il(1,2,2)=il2
          il(1,2,3)=il3
!         2 POSSIBLE FORMS OF THE HIGHER TRIANGLE (GLOBAL AND BOUNDARY)
          ig(2,1,1)=ig1
          ig(2,1,2)=ig3
          ig(2,1,3)=ig4
          ig(2,2,1)=ig2
          ig(2,2,2)=ig3
          ig(2,2,3)=ig4
          il(2,1,1)=il1
          il(2,1,2)=il3
          il(2,1,3)=il4
          il(2,2,1)=il2
          il(2,2,2)=il3
          il(2,2,3)=il4
!
          ok(1)=.false.
          ok(2)=.false.
!
!         K=1 LOWER TRIANGLE   K=2 HIGHER TRIANGLE
          DO k=1,2
!           2 POSSIBLE SPLITS
            DO l=1,2
!             12 WAYS FOR A TETRAHEDRON OF PRESENTING ITS FACES
              DO m=1,12
!               3 POSSIBLE TETRAHEDRONS
                DO n=1,3
                  IF(ig(k,l,1).EQ.ikle3(ielem,n,ietage,num1(m)).AND.
     &               ig(k,l,2).EQ.ikle3(ielem,n,ietage,num2(m)).AND.
     &               ig(k,l,3).EQ.ikle3(ielem,n,ietage,num3(m))) THEN
!                   STORAGE LIKE IKLBOR(NELEB2D,2,NETAGE,3)
                    ieleb3=(2*ietage+k-3)*neleb2d+ieleb
                    iklbor(ieleb3,1) = il(k,l,1)
                    iklbor(ieleb3,2) = il(k,l,2)
                    iklbor(ieleb3,3) = il(k,l,3)
                    nelbor(ieleb3)   = t(n)
                    nulone(ieleb3,1) = num1(m)
                    nulone(ieleb3,2) = num2(m)
                    nulone(ieleb3,3) = num3(m)
                    ok(k) = .true.
                  ENDIF
                ENDDO
              ENDDO
            ENDDO
          ENDDO
          IF(.NOT.ok(1).OR..NOT.ok(2)) THEN
            WRITE(lu,*) 'PB IN ELEB3DT IELEM=',ielem,' IPTFR=',iptfr
            CALL plante(1)
            stop
          ENDIF
!
        ENDDO
!
      ENDDO
!
!-----------------------------------------------------------------------
!
!     NELEB IS THE 3D VALUE
!
      neleb=2*neleb2d*netage
!
!-----------------------------------------------------------------------
!
      RETURN
      END