inbief.f

Go to the documentation of this file.

!                   *****************
                    SUBROUTINE inbief
!                   *****************
!
     &(lihbor,klog,it1,it2,it3,lvmac,ielmx,
     & lambd0,spheri,mesh,t1,t2,optass,produc,equa,mesh2d)
!
!***********************************************************************
! BIEF   V7P3
!***********************************************************************
!
!brief    PREPARES THE DATA STRUCTURE FOR BIEF.
!+                THE INTEGER AND REAL ARRAYS DESCRIBING THE MESH
!+                ARE BUILT AND STORED IN MESH.
!
!history  J-M HERVOUET (LNHE) ; REGINA NEBAUER; LAM MINH PHUONG; EMILE RAZAFINDRAKOTO
!+        05/02/2010
!+        V6P0
!+   COMP_SEG UPDATED FOR QUADRATIC TRIANGLES
!
!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) ; REGINA NEBAUER; LAM MINH PHUONG; EMILE RAZAFINDRAKOTO
!+        07/12/2011
!+        V6P2
!+   Checking elements for building GLOSEG. New optional argument MESH2D.
!+   for prisms split into tetrahedrons (call of STOSEG51)
!
!history  J-M HERVOUET (LNHE)
!+        10/09/2012
!+        V6P2
!+   Different call to VOISIN31, and call added for IELM=51
!
!history  J-M HERVOUET (LNHE)
!+        06/03/2013
!+        V6P3
!+   XEL, YEL now built in 3D. SURFAC based on XEL and YEL. Mercator
!+   projection treated in 3D.
!
!history R. ATA (EDF R&D - LNHE)
!+        21/05/2013
!+        V6P3
!+   add centre_mass_seg, new infcel, hloc.
!
!history  J-M HERVOUET (EDF LAB, LNHE)
!+        28/03/2014
!+        V7P0
!+   Now written to enable different numbering of boundary points and
!+   boundary segments.
!
!history  J-M HERVOUET (EDF LAB, LNHE)
!+        09/05/2014
!+        V7P0
!+   Allocation of new vectors of I8 integers MESH%WI8 and MESH%TI8
!+   for finite element assembly with integers.
!history  R.NHEILI (Univerte de Perpignan, DALI)
!+        24/02/2016
!+        V7P3
!+        ALLOCATE BUF_SEND_ERR AND BUD_RECV_ERR
!
!history  J-M HERVOUET (EDF LAB, LNHE)
!+        30/07/2015
!+        V7P1
!+   Call of make_eltcar modified to compute IFAC. FAC suppressed.
!
!history  M.S.TURNBULL (HRW)
!+        24/09/2015
!+        V7P1
!+   Correction to the computation of the both Y and XDEL in spherical
!+      coordinates, accounting for the fact that the Earth is round.
!
!history  J-M HERVOUET (Jubilado)
!+        10/09/2017
!+        V7P3
!+   MESH%IKLE%DIM1 replaced with MESH%IKLE%MAXDIM1 in calls to VOISIN,
!+   VOISIN31 and ELEBD. NELEM replaced with MESH%ELTSEG%DIM1 in the
!+   call to COMP_NH_COM_SEG.
!
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!| EQUA           |-->| IDENTIFICATION OF PROGRAM OR EQUATIONS SOLVED
!| IELMX          |-->| THE MORE COMPLEX ELEMENT USED (FOR MEMORY)
!| IT1            |<->| INTEGER WORK ARRAY IN A BIEF_OBJ STRUCTURE
!| IT2            |<->| INTEGER WORK ARRAY IN A BIEF_OBJ STRUCTURE
!| IT3            |<->| INTEGER WORK ARRAY IN A BIEF_OBJ STRUCTURE
!| KLOG           |-->| CONVENTION FOR SOLID BOUNDARY
!| LAMBD0         |-->| LATITUDE OF ORIGIN POINT (SPHERICAL COORDINATES)
!| LIHBOR         |-->| TYPES OF BOUNDARY CONDITIONS ON DEPTH
!| LVMAC          |-->| VECTOR LENGTH (IF VECTOR MACHINE)
!| MESH           |-->| MESH STRUCTURE
!| MESH2D         |-->| UNDERLYING 2D MESH (FOR PRISMS AND PRISMS SPLIT
!|                |   | INTO TETRAHEDRONS)
!| OPTASS         |-->| OPTION FOR MATRIX STORAGE.
!| PRODUC         |-->| OPTION FOR MATRIX x VECTOR PRODUCT.
!| SPHERI         |-->| LOGICAL, IF YES : SPHERICAL COORDINATES.
!| T1             |<->| WORK BIEF_OBJ STRUCTURE
!| T2             |<->| WORK BIEF_OBJ STRUCTURE
!~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
!
      USE bief, ex_inbief => inbief
      USE declarations_telemac, ONLY : modass, check_mesh
!
      USE declarations_special
      IMPLICIT NONE
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER, INTENT(IN)            :: IELMX,OPTASS,PRODUC,KLOG,LVMAC
      INTEGER, INTENT(IN)            :: LIHBOR(*)
      DOUBLE PRECISION, INTENT(IN)   :: LAMBD0
      LOGICAL, INTENT(IN)            :: SPHERI
      CHARACTER(LEN=20)              :: EQUA
      TYPE(bief_mesh), INTENT(INOUT) :: MESH
      TYPE(bief_obj), INTENT(INOUT)  :: T1,T2,IT1,IT2,IT3
      TYPE(bief_mesh), INTENT(INOUT), OPTIONAL :: MESH2D
!
!+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
!
      INTEGER I,IELEM,NELEM,NELMAX,NPTFR,NPOIN,IELM,IPLAN,I3D
      INTEGER MXPTVS,NPLAN
      INTEGER LV,NDP,IDP,I1,I2,I3,NPOIN2
      INTEGER NPTFR2,NELEM2,NELMAX2,NELEB2,NELEB
      INTEGER NELCOU,IELEM3D
!
      DOUBLE PRECISION X2,X3,Y2,Y3
!
!-----------------------------------------------------------------------
!     FOR CALL TO VOISIN31
      INTEGER IKLESTR(1,3)
!
!     DEPLOYMENT OF THE DATA STRUCTURE
!
      nelem = mesh%NELEM
      nelmax= mesh%NELMAX
      npoin = mesh%NPOIN
      ielm  = mesh%X%ELM
      ndp   = bief_nbpel(ielm,mesh)
      nptfr = mesh%NPTFR
      neleb = mesh%NELEB
!
!     WITH PRISMS, DIFFERENT FROM 2D VALUES, OTHERWISE
!
      IF(ielm.EQ.41.OR.ielm.EQ.51) THEN
        npoin2  =bief_nbpts(11,mesh)
        nelem2  =bief_nbpts(10,mesh)
        nelmax2 =bief_nbmpts(10,mesh)
        nptfr2  =bief_nbpts(1,mesh)
        nplan   =npoin/npoin2
      ELSEIF(ielm.EQ.11.OR.ielm.EQ.31) THEN
        npoin2  =npoin
        nelem2  =nelem
        nelmax2 =nelmax
        nptfr2  =nptfr
        neleb2  =neleb
        nplan   =1
      ELSE
        WRITE(lu,*) 'UNEXPECTED ELEMENT IN INBIEF:',ielm
        CALL plante(1)
        stop
      ENDIF
!
!     FINITE ELEMENT ASSEMBLY WITH I8 INTEGERS
!
      IF(modass.EQ.2) THEN
        ALLOCATE(mesh%WI8(nelmax*ndp))
        ALLOCATE(mesh%TI8(npoin))
      ENDIF
!
!     PARALLEL MODE : INITIALISES THE ARRAYS NHP,NHM
!                        INDPU,FAC, ETC.
!
!
      IF(ncsize.GT.1) THEN
!
        CALL parini(mesh%NHP%I,mesh%NHM%I,mesh%INDPU%I,
     &              npoin2,mesh%NACHB%I,nplan,mesh,
     &              mesh%NB_NEIGHB,mesh%NB_NEIGHB_SEG,
     &              nelem2,mesh%IFAPAR%I,modass)
!
      ELSE
!       THESE STUCTURES ARE ALLOCATED IN PARINI
        CALL bief_allvec(2,mesh%NB_NEIGHB_PT,'NBNGPT',0,1,0,mesh)
        CALL bief_allvec(2,mesh%LIST_SEND   ,'LSSEND',0,1,0,mesh)
        CALL bief_allvec(2,mesh%NH_COM      ,'NH_COM',0,1,0,mesh)
        CALL bief_allvec(2,mesh%NB_NEIGHB_PT_SEG,'NBNGSG',0,1,0,mesh)
        CALL bief_allvec(2,mesh%LIST_SEND_SEG,'LSSESG',0,1,0,mesh)
        CALL bief_allvec(2,mesh%NH_COM_SEG  ,'NH_CSG',0,1,0,mesh)
        CALL bief_allvec(1,mesh%BUF_SEND    ,'BUSEND',0,1,0,mesh)
        CALL bief_allvec(1,mesh%BUF_RECV    ,'BURECV',0,1,0,mesh)
        IF (modass .EQ.3) THEN
          CALL bief_allvec(1,mesh%BUF_SEND_ERR ,'BUSEND_ERR',0,1,0,mesh)
          CALL bief_allvec(1,mesh%BUF_RECV_ERR ,'BURECV_ERR',0,1,0,mesh)
        ENDIF
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!     COMPUTES THE NEIGHBOURS OF THE BOUNDARY FACES (TRIANGULAR MESH)
!
!     NOTE: SEE CPIKLE2 AND CPIKLE3 IN 3D. IKLE CAN HERE BE 3D BECAUSE
!           THE BEGINNING OF IKLE IN 3D IS THE SAME AS THAT IN 2D (THE
!           FIRST 3 POINTS OF THE PRISMS OR TETRAHEDRONS CORRESPOND
!           TO THE 3 POINTS OF THE BOTTOM TRIANGLES)
!
!
      IF(ielm.EQ.11.OR.ielm.EQ.41.OR.ielm.EQ.51) THEN
        CALL voisin(mesh%IFABOR%I,nelem2,nelmax2,ielm,mesh%IKLE%I,
     &              mesh%IKLE%MAXDIM1,
     &              npoin2,mesh%NACHB%I,mesh%NBOR%I,nptfr2,it1%I,it2%I)
!
      ELSEIF(ielm.NE.31) THEN
        WRITE(lu,*) 'UNEXPECTED ELEMENT IN INBIEF:',ielm
        CALL plante(1)
        stop
      ENDIF
!
!-----------------------------------------------------------------------
!
      IF(ielm.EQ.11.OR.ielm.EQ.41.OR.ielm.EQ.51) THEN
!
!       CASES WITH A FIRST CALL IN 2D
!
        mxptvs = mesh%MXPTVS
!       HERE IFABOR FOR IELM=51 MUST STILL BE 2D
!       SO VOISIN31 CALLED LATER
!
!       NOTE: IN 3D IKLBOR BUILT HERE WITH NELEB 2D AND NELEBX 3D.
!
        CALL elebd(mesh%NELBOR%I,mesh%NULONE%I,mesh%KP1BOR%I,
     &             mesh%IFABOR%I,mesh%NBOR%I,mesh%IKLE%I,
     &             mesh%IKLE%MAXDIM1,
     &             mesh%IKLBOR%I,nelem2,nelmax2,npoin2,nptfr2,ielm,
     &             lihbor,klog,mesh%ISEG%I,
     &             it1%I,it2%I,it3%I,mesh%NELEBX,mesh%NELEB)
!
      ENDIF
!
!     3D CASES
!
      IF(ielm.EQ.31) THEN
!
!       BUILDING ARRAYS FOR TETRAHEDRONS
!
        CALL voisin31(mesh%IFABOR%I,nelem,nelmax,ielm,mesh%IKLE%I,
     &                mesh%IKLE%MAXDIM1,npoin,mesh%NBOR%I,nptfr,
     &                lihbor,klog,mesh%INDPU%I,iklestr,neleb2)
!
        CALL elebd31(mesh%NELBOR%I,mesh%NULONE%I,mesh%IKLBOR%I,
     &               mesh%IFABOR%I,mesh%NBOR%I,mesh%IKLE%I,
     &               nelem,neleb,nelmax,npoin,nptfr,ielm)
!
      ELSEIF(ielm.EQ.41) THEN
!
!       COMPLETES ARRAYS FOR PRISMS
!
        CALL eleb3d(mesh%IKLE%I,mesh%NBOR%I,
     &              mesh%NELBOR%I,mesh%IKLBOR%I,
     &              mesh%NELEB,mesh%NELEBX,
     &              mesh%NULONE%I,nelmax2,npoin2,nplan,nplan-1,nptfr2)
!
      ELSEIF(ielm.EQ.51) THEN
!
!       COMPLETES ARRAYS FOR PRISMS SPLIT INTO TETRAHEDRONS
!
        IF(PRESENT(mesh2d)) THEN
!         NOTE THE USE OF MESH2D FOR NELBOR AND NULONE
!         THIS IS FOR CALLING STOSEG
          CALL eleb3dt(mesh%IKLE%I,mesh%NBOR%I,mesh%NELBOR%I,
     &                 mesh2d%NELBOR%I,mesh%IKLBOR%I,
     &                 mesh%NELEB,mesh%NELEBX,mesh%NULONE%I,nelmax2,
     &                 npoin2,nplan,nplan-1,nptfr2,
     &                 mesh2d%IKLBOR%I,mesh2d%NELEB,mesh2d%NELEBX)
        ELSE
          WRITE(lu,*) 'ARGUMENT MESH2D SHOULD BE ADDED TO INBIEF'
          WRITE(lu,*) 'FOR A CALL WITH IELM=51'
          CALL plante(1)
          stop
        ENDIF
!
      ELSEIF(ielm.NE.11) THEN
!
        WRITE(lu,*) 'INBIEF UNEXPECTED ELEMENT: ',ielm
        CALL plante(1)
        stop
!
      ENDIF
!
!-----------------------------------------------------------------------
!
! LOOKS FOR VECTORISATION POSSIBILITIES
!
      IF(ielm.EQ.11) THEN
!
      IF(lvmac.NE.1) THEN
        WRITE(lu,201) lvmac
201     FORMAT(1x,'INBIEF (BIEF): VECTOR MACHINE',/,1x,
     &  'WITH VECTOR LENGTH :',1i6,
     &  ' (ACCORDING TO YOUR DATA OR IN THE DICTIONNARY OF KEY-WORDS)')
        CALL veclen(lv,ndp,mesh%IKLE%I,nelem,nelmax,npoin,t1%R)
        IF(lv.LT.lvmac) THEN
          WRITE(lu,301) lv
301       FORMAT(1x,'THIS LENGTH IS REDUCED TO ',1i4,' BY THE NUMBERING
     &OF THE ELEMENTS (SEE STBTEL DOCUMENTATION)')
        ENDIF
      ELSE
        lv = 1
        WRITE(lu,401)
401     FORMAT(1x,'INBIEF (BIEF): NOT A VECTOR MACHINE',
     &                                      ' (ACCORDING TO YOUR DATA)')
      ENDIF
!
      mesh%LV = lv
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!     MERCATOR PROJECTION (TRIANGLES AND PRISMS ONLY)
!
      IF(spheri.AND.ielm.NE.11.AND.ielm.NE.41) THEN
        WRITE(lu,399)
399     FORMAT(1x,'INBIEF (BIEF) : ELEMENT NOT IMPLEMENTED WITH',/,1x,
     &            'MERCATOR PROJECTION:',1i3)
        CALL plante(1)
        stop
      ENDIF
!
      IF(spheri) THEN
!
        CALL cpstvc(mesh%X,t1)
        CALL cpstvc(mesh%Y,t2)
        CALL latitu(t2%R,mesh%COSLAT%R,mesh%SINLAT%R,
     &              lambd0,mesh%Y%R,npoin2)
        CALL corlat
!
        IF(ielm.EQ.11.OR.ielm.EQ.41) THEN
          DO i=1,npoin2
            t1%R(i)=mesh%X%R(i)
!            T2%R(I)=MESH%Y%R(I)*MESH%COSLAT%R(I)
          ENDDO
        ENDIF
!       COMPLETING UPPER LAYERS FOR 3D MESHES
        IF(ielm.EQ.41) THEN
          DO iplan=2,nplan
            DO i=1,npoin2
              i3d=(iplan-1)*npoin2+i
              t1%R(i3d)=mesh%X%R(i3d)
              t2%R(i3d)=t2%R(i)
            ENDDO
          ENDDO
        ENDIF
!
!       CONVERTS TO COORDINATES BY ELEMENTS (STARTING WITH X AND Y)
        CALL pttoel(mesh%XEL,t1,mesh)
        CALL pttoel(mesh%YEL,t2,mesh)
!
      ELSE
!
!       NOTE: IN 3D MESH%X AND MESH%Y FULLY BUILT IN ALMESH
!
!       CONVERTS TO COORDINATES BY ELEMENTS (STARTING WITH X AND Y)
!
        CALL pttoel(mesh%XEL,mesh%X,mesh)
        CALL pttoel(mesh%YEL,mesh%Y,mesh)
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!     CONVERTS TO A LOCAL SYSTEM IN X AND Y, WITH POINT 1 AT ORIGIN
!
      DO idp=2,ndp
        CALL ov_2('X=X-Y   ',mesh%XEL%R,idp,
     &                       mesh%XEL%R,1  ,
     &                       mesh%XEL%R,1  , 0.d0 , nelmax , nelem )
        CALL ov_2('X=X-Y   ',mesh%YEL%R,idp,
     &                       mesh%YEL%R,1  ,
     &                       mesh%YEL%R,1  , 0.d0 , nelmax , nelem )
      ENDDO
!
      CALL ov('X=C     ', x=mesh%XEL%R, c=0.d0, dim1=nelem)
      CALL ov('X=C     ', x=mesh%YEL%R, c=0.d0, dim1=nelem)
!
!
      IF(spheri) THEN
!
        IF(PRESENT(mesh2d)) THEN
          CALL longitu(mesh2d%XEL%R,mesh%COSLAT%R,mesh2d%IKLE%I,
     &                 nelmax2,nelem2)
          nelcou = (nplan-1)*nelmax2
          DO iplan=1,nplan-1
            DO ielem=1,nelem2
              ielem3d=(iplan-1)*nelmax2+ielem
!             3 POINTS OF THE BOTTOM
              mesh%XEL%R(ielem3d)         =mesh2d%XEL%R(ielem)
              mesh%XEL%R(ielem3d+nelcou  )=mesh2d%XEL%R(ielem+nelmax2)
              mesh%XEL%R(ielem3d+2*nelcou)=mesh2d%XEL%R(ielem+2*nelmax2)
!             3 POINTS OF THE TOP
              mesh%XEL%R(ielem3d+3*nelcou)=mesh%XEL%R(ielem3d)
              mesh%XEL%R(ielem3d+4*nelcou)=mesh%XEL%R(ielem3d+nelcou)
              mesh%XEL%R(ielem3d+5*nelcou)=mesh%XEL%R(ielem3d+2*nelcou)
            ENDDO
          ENDDO
        ELSE
          CALL longitu(mesh%XEL%R,mesh%COSLAT%R,mesh%IKLE%I,
     &                 nelmax,nelem)
        ENDIF
!
      ENDIF
!
!     IF DONE FOR Z (BUT IN MOVING MESHES SHOULD NOT BE USED !!!!)
!
!     IF(MESH%DIM1.EQ.3) THEN
!       CALL PTTOEL(MESH%ZEL,MESH%Z,MESH)
!       DO IDP=2,NDP
!         CALL OV_2('X=X-Y   ',MESH%ZEL%R,IDP,
!    &                         MESH%ZEL%R,1  ,
!    &                         MESH%ZEL%R,1  , 0.D0 , NELMAX , NELEM )
!       ENDDO
!       CALL OV('X=C     ', X=MESH%ZEL%R, C=0.D0, NELEM )
!     ENDIF
!
!-----------------------------------------------------------------------
!
! COMPUTES THE GEOMETRICAL COEFFICIENTS FOR EACH ELEMENT
!
      IF(ielm.EQ.11) THEN
!
        CALL geoelt(mesh%SURDET%R,mesh%SURFAC%R,
     &              mesh%XEL%R   ,mesh%YEL%R   ,nelem,nelmax,ielm)
!
! FOR THE TIME BEING, SURDET IS ONLY USED BY CARACT, WHICH DOES NOT
! WORK ON THE MESH IN SPHERICAL COORDINATES.
! ERASES SURDET COMPUTED BY GEOELE FROM XEL AND YEL
!
        IF(spheri) THEN
!
          DO ielem = 1 , nelem
            i1 = mesh%IKLE%I(ielem)
            i2 = mesh%IKLE%I(ielem+nelmax)
            i3 = mesh%IKLE%I(ielem+2*nelmax)
            x2 = - mesh%X%R(i1) + mesh%X%R(i2)
            x3 = - mesh%X%R(i1) + mesh%X%R(i3)
            y2 = - mesh%Y%R(i1) + mesh%Y%R(i2)
            y3 = - mesh%Y%R(i1) + mesh%Y%R(i3)
            mesh%SURDET%R(ielem) = 1.d0 / (x2*y3 - x3*y2)
          ENDDO
!
        ENDIF
!
      ELSEIF(ielm.EQ.41.OR.ielm.EQ.51.OR.ielm.EQ.31) THEN
!
!       FOR PRISMS, SURFAC IS THE SURFACE OF THE TRIANGLES
!       FOR ELEMENTS 51 AND 31 ???????? SHOULD NOT BE USED...
!
        DO ielem = 1 , nelem
          x2 = mesh%XEL%R(ielem+nelmax)
          x3 = mesh%XEL%R(ielem+2*nelmax)
          y2 = mesh%YEL%R(ielem+nelmax)
          y3 = mesh%YEL%R(ielem+2*nelmax)
          mesh%SURFAC%R(ielem) = 0.5d0 * (x2*y3 - x3*y2)
        ENDDO
!
      ELSE
        WRITE(lu,*) 'UNEXPECTED ELEMENT IN INBIEF:',ielm
        CALL plante(1)
        stop
      ENDIF
!
!-----------------------------------------------------------------------
!
! DEFINES THE OUTGOING NORMALS AT THE BOUNDARIES
!         AND THE DISTANCES TO THE BOUNDARY
!
      IF(ielm.EQ.11) THEN
!
      CALL normab(mesh%XNEBOR%R,mesh%YNEBOR%R,
     &            mesh%XSGBOR%R,mesh%YSGBOR%R,
     &            mesh%DISBOR%R,mesh%SURFAC%R,nelmax,mesh%NELBOR%I,
     &            mesh%NULONE%I,mesh%LGSEG%R,nptfr,mesh,
     &            mesh%XEL%R,mesh%YEL%R,mesh%IKLBOR%I,
     &            mesh%NELEBX,mesh%NELEB)
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!  DATA STRUCTURE FOR EDGE-BASED STORAGE (FROM 5.9 ON ALWAYS DONE IN 2D)
!                                        (FROM 6.2 ON ALWAYS DONE IN 3D)
!  SEE CALL TO COMP_SEG BELOW TO COMPLETE THE STRUCTURE
!
      IF(ielm.EQ.11) THEN
!
      CALL stoseg(mesh%IFABOR%I,nelem,nelmax,nelmax,ielmx,mesh%IKLE%I,
     &            mesh%NBOR%I,nptfr,mesh%GLOSEG%I,mesh%GLOSEG%MAXDIM1,
     &            mesh%ELTSEG%I,mesh%ORISEG%I,mesh%NSEG,
     &            mesh%NELBOR%I,mesh%NULONE%I,
     &            mesh%KNOLG%I,mesh%IKLBOR%I,mesh%NELEBX,mesh%NELEB)
!
      ELSEIF(ielm.EQ.41) THEN
!
      CALL stoseg41(mesh%IFABOR%I,nelmax,ielmx,mesh%IKLE%I,mesh%NBOR%I,
     &              mesh%GLOSEG%I,mesh%GLOSEG%MAXDIM1,
     &              mesh%ELTSEG%I,mesh%ORISEG%I,
     &              mesh%NELBOR%I,mesh%NULONE%I,
     &              nelmax2,nelem2,nptfr2,npoin2,nplan,mesh%KNOLG%I,
     &              bief_nbseg(11,mesh),
     &              mesh%IKLBOR%I,mesh%NELEBX,mesh%NELEB)
!
      ELSEIF(ielm.EQ.51) THEN
!
      IF(PRESENT(mesh2d)) THEN
!       NOTE THE USE OF MESH2D FOR NELBOR AND NULONE
!       THIS IS FOR CALLING STOSEG
        CALL stoseg51(mesh%IFABOR%I,nelmax,ielmx,
     &                mesh%IKLE%I,mesh%NBOR%I,
     &                mesh%GLOSEG%I,mesh%GLOSEG%MAXDIM1,
     &                mesh%ELTSEG%I,mesh%ORISEG%I,
     &                mesh2d%NELBOR%I,mesh2d%NULONE%I,
     &                nelmax2,nelem2,nptfr2,npoin2,nplan,mesh%KNOLG%I,
     &                mesh2d%NSEG,mesh2d%IKLBOR%I,mesh2d%NELEB,
     &                mesh2d%NELEBX)
      ELSE
        WRITE(lu,*) 'ARGUMENT MESH2D SHOULD BE ADDED TO INBIEF'
        WRITE(lu,*) 'FOR A CALL WITH IELM=51'
        CALL plante(1)
        stop
      ENDIF
!
      ELSE
!
        WRITE(lu,*) 'ELEMENT ',ielm,' NOT IMPLEMENTED FOR SEGMENTS'
        CALL plante(1)
        stop
!
      ENDIF
!
!     NOW THE 3D VALUE OF IFABOR IS BUILT FOR PRISMS CUT INTO
!     TETRAHEDRA (UP TO STOSEG51 A 2D VALUE WAS USED)
!
      IF(ielm.EQ.51) THEN
        CALL voisin31(mesh%IFABOR%I,nelem,nelmax,ielm,mesh%IKLE%I,
     &                mesh%IKLE%MAXDIM1,npoin,mesh%NBOR%I,nptfr,
     &                lihbor,klog,mesh%INDPU%I,iklestr,1)
      ENDIF
!
!-----------------------------------------------------------------------
!
      IF(ncsize.GT.1.AND.ielm.EQ.11) THEN
!
!       COMPLETES NH_COM_SEG WITH SEGMENT NUMBERS ONCE ELTSEG IS KNOWN
!
        CALL comp_nh_com_seg(mesh%ELTSEG%I,mesh%ELTSEG%DIM1,
     &                       mesh%NH_COM_SEG%I,
     &                       mesh%NH_COM_SEG%DIM1,mesh%NB_NEIGHB_SEG,
     &                       mesh%NB_NEIGHB_PT_SEG%I,
     &                       mesh%GLOSEG%I,mesh%GLOSEG%DIM1,
     &                       mesh%KNOLG%I,npoin)
!
!       COMPLETES FAC ONCE IFABOR AND ELTSEG ARE KNOWN
!
        IF(ielm.EQ.11.AND.ielmx.EQ.13) THEN
          CALL comp_fac(mesh%ELTSEG%I,mesh%ORISEG%I,mesh%IFABOR%I,nelem,
     &                  npoin,mesh%IFAC)
        ENDIF
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!  DATA STRUCTURE FOR EDGE-BASED STORAGE
!
      IF(ielm.EQ.11.AND.produc.EQ.2) THEN
!
        CALL fropro(mesh%NBOR%I,mesh%IKLE%I,
     &              nelem,nelmax,npoin,mesh%NPMAX,nptfr,ielm,
     &              mesh%IKLEM1%I,mesh%LIMVOI%I,optass,produc,mxptvs,
     &              it1%I,mesh%GLOSEG%I,mesh%GLOSEG%DIM1,mesh%NSEG)
!
      ENDIF
!
!-----------------------------------------------------------------------
!
!  COMPLEMENTS IKLE AND NBOR BEYOND LINEAR ELEMENTS
!
      IF(ielm.EQ.11.AND.ielm.NE.ielmx) THEN
        IF(mesh%IKLE%DIM2.NE.bief_nbpel(ielmx,mesh)) THEN
          WRITE(lu,101) ielmx
101       FORMAT(1x,'INBIEF (BIEF): WRONG DIMENSION OF IKLE',/,1x,
     &              'FOR AN ELEMENT WITH TYPE :',1i6)
          CALL plante(1)
          stop
        ENDIF
        CALL comp_ikle(mesh%IKLE%I,mesh%IKLBOR%I,
     &                 mesh%ELTSEG%I,mesh%NBOR%I,mesh%NELBOR%I,
     &                 mesh%NULONE%I,ielmx,nelem,nelmax,npoin,nptfr,
     &                 mesh%NELEB,mesh%NELEBX)
      ENDIF
!
!-----------------------------------------------------------------------
!
! COMPLEMENTS THE SEGMENT STRUCTURE BEYOND THE LINEAR ELEMENTS
!
      IF(ielm.NE.ielmx) THEN
        CALL comp_seg(nelem,nelmax,ielmx,mesh%IKLE%I,mesh%GLOSEG%I,
     &                mesh%GLOSEG%MAXDIM1,mesh%ELTSEG%I,mesh%ORISEG%I,
     &                mesh%NSEG)
      ENDIF
!
!-----------------------------------------------------------------------
!
!     V6P3 FOR NEW DATA STRUCTURE OF FINITE VOLUMES
!     COMPUTE THE COORDINATES OF CENTRE OF GRAVITY FOR ELEMENTS RIGHT AND
!     LEFT OF EDGES
!
      IF(equa(1:15).EQ.'SAINT-VENANT VF') THEN
!
        CALL centre_mass_seg(mesh%X%R,mesh%Y%R,mesh%COORDG%R,
     &                       mesh%IKLE%I,npoin,mesh%ELTSEG%I,
     &                       mesh%ORISEG%I,nelem,mesh%NSEG,
     &                       mesh%JMI%I,mesh%CMI%R,mesh%GLOSEG%I,
     &                       mesh%IFABOR%I,mesh)
!
      ENDIF
!
!-----------------------------------------------------------------------
!
! COMPLEMENTS THE DATA STRUCTURE FOR FINITE VOLUMES
!
      IF(equa(1:15).EQ.'SAINT-VENANT VF') THEN
!
        CALL infcel(mesh%X%R,mesh%Y%R,
     &              mesh%NUBO%I,mesh%VNOIN%R,npoin,
     &              nelem,mesh%NSEG,mesh%CMI%R,
     &              mesh%AIRST%R,mesh%GLOSEG%I,
     &              mesh%COORDG%R,mesh%ELTSEG%I,
     &              mesh%ORISEG%I,mesh%IFABOR%I)
!
!       COMPUTES THE SURFACE OF THE CELLS
!
        CALL vector(t1,'=','MASBAS          ',11,
     &              1.d0,t2,t2,t2,t2,t2,t2,mesh,.false.,t2)
        IF(ncsize.GT.1) CALL parcom(t1,2,mesh)
!
!       COMPUTES THE LOCAL SPACE STEP PER CELL
!
        CALL hloc(npoin,mesh%NSEG,nelem,mesh%NUBO%I,mesh%VNOIN%R,t1%R,
     &            mesh%DTHAUT%R,mesh,mesh%ELTSEG%I,mesh%IFABOR%I)
!
!       COMPUTES THE GRADIENTS OF THE BASE FUNCTIONS
!
        CALL gradp(npoin,mesh%NELMAX,mesh%IKLE%I,mesh%SURFAC%R,
     &               mesh%X%R,mesh%Y%R,mesh%DPX%R,mesh%DPY%R)
!
      ENDIF
!
!-----------------------------------------------------------------------
!
! COMPUTES THE STARTING ELEMENT FOR THE METHOD OF CHARACTERISTICS
!
      CALL make_eltcar(mesh%ELTCAR%I,mesh%IFAC%I,
     &                 mesh%IKLE%I,npoin2,nelem2,
     &                 nelmax,mesh%KNOLG%I,it1%I,mesh,nplan,ielmx)
!
!-----------------------------------------------------------------------
!
!     CHECKING THE MESH
!
      IF(check_mesh) THEN
        CALL checkmesh(mesh,npoin)
      ENDIF
!
!-----------------------------------------------------------------------
!
      RETURN
      END