i20sti3.F

Go to the documentation of this file.

Copyright>        OpenRadioss
Copyright>        Copyright (C) 1986-2025 Altair Engineering Inc.
Copyright>
Copyright>        This program is free software: you can redistribute it and/or modify
Copyright>        it under the terms of the GNU Affero General Public License as published by
Copyright>        the Free Software Foundation, either version 3 of the License, or
Copyright>        (at your option) any later version.
Copyright>
Copyright>        This program is distributed in the hope that it will be useful,
Copyright>        but WITHOUT ANY WARRANTY; without even the implied warranty of
Copyright>        MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
Copyright>        GNU Affero General Public License for more details.
Copyright>
Copyright>        You should have received a copy of the GNU Affero General Public License
Copyright>        along with this program.  If not, see <https://www.gnu.org/licenses/>.
Copyright>
Copyright>
Copyright>        Commercial Alternative: Altair Radioss Software
Copyright>
Copyright>        As an alternative to this open-source version, Altair also offers Altair Radioss
Copyright>        software under a commercial license.  Contact Altair to discuss further if the
Copyright>        commercial version may interest you: https://www.altair.com/radioss/.
!||====================================================================
!||    i20sti3                ../starter/source/interfaces/inter3d1/i20sti3.F
!||--- called by ------------------------------------------------------
!||    i20ini3                ../starter/source/interfaces/inter3d1/i20ini3.F
!||--- calls      -----------------------------------------------------
!||    ancmsg                 ../starter/source/output/message/message.f
!||    bitget                 ../starter/source/interfaces/inter3d1/bitget.F
!||    bitset                 ../starter/source/interfaces/inter3d1/bitget.F
!||    bitunset               ../starter/source/interfaces/inter3d1/bitget.F
!||    i20nelts               ../starter/source/interfaces/inter3d1/inelt.F
!||    i4gmx3                 ../starter/source/interfaces/inter3d1/i4gmx3.F
!||    incoq3                 ../starter/source/interfaces/inter3d1/incoq3.F
!||    ineltc                 ../starter/source/interfaces/inter3d1/inelt.F
!||    insol3                 ../starter/source/interfaces/inter3d1/insol3.F
!||    volint                 ../starter/source/interfaces/inter3d1/volint.F
!||--- uses       -----------------------------------------------------
!||    message_mod            ../starter/share/message_module/message_mod.F
!||====================================================================
      SUBROUTINE i20sti3(
     1 PM        ,GEO       ,X         ,MS        ,
     2 IXS       ,IXC       ,IXTG      ,IXT       ,
     3 IXP       ,WA        ,NINT      ,NTY       ,
     4 NOINT     ,NRT       ,NSN       ,IRECT     ,
     5 NSV       ,INACTI    ,GAP       ,IGAP      ,
     6 GAP_S     ,GAP_M     ,GAPMIN    ,GAPINF    ,
     7 GAPMAX    ,STFAC     ,STF       ,STFN      ,
     8 KNOD2ELS  ,KNOD2ELC  ,KNOD2ELTG ,NOD2ELS   ,
     9 NOD2ELC   ,NOD2ELTG  ,IGRSURF1  ,IFS2      ,
     A IGRSURF2  ,INTTH     ,IELES     ,
     B IELEC     ,AREAS     ,IPARTC    ,IPARTTG   ,
     C THK       ,THK_PART  ,GAP_SH    ,XANEW     ,
     D GAPSHMAX  ,NBINFLG   ,MBINFLG   ,NLN       ,
     E NLG       ,GAPSOL    ,IXS10     ,IXS16     ,
     F IXS20     ,ID,TITR   ,IGEO      ,PM_STACK  ,
     G IWORKSH )
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE groupdef_mod
      USE message_mod
      USE names_and_titles_mod , ONLY : nchartitle
      use element_mod , only : nixs,nixc,nixtg,nixt,nixp
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "param_c.inc"
#include      "com01_c.inc"
#include      "com04_c.inc"
#include      "scr05_c.inc"
#include      "scr08_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER NRT, NINT, NTY, NOINT,NSN,IGAP,
     .        INACTI,IFS2,NLN      
      INTEGER IRECT(4,*), IXS(NIXS,*), IXC(NIXC,*),
     .   NSV(*), IXTG(NIXTG,*), IXT(NIXT,*), IXP(NIXP,*),
     .   KNOD2ELS(*), KNOD2ELC(*), KNOD2ELTG(*), NOD2ELS(*), NOD2ELC(*), 
     .   NOD2ELTG(*),IELES(*),INTTH,IELEC(*),
     .   IPARTC(*), IPARTTG(*),NBINFLG(*),MBINFLG(*),NLG(*) ,
     .   IXS10(6,*), IXS16(*), IXS20(*),  IGEO(NPROPGI,*),IWORKSH(3,*)
C     REAL
      my_real
     .   STFAC, GAP,GAPMIN,GAPINF, GAPMAX,GAPSHMAX,GAPSOLIDMAX,GAPSOL
C     REAL
      my_real
     .   X(3,*), STF(*), PM(NPROPM,*), GEO(NPROPG,*),  STFN(*),
     .   MS(*),WA(*),GAP_S(*),GAP_M(*),GAP_SH(*),AREAS(*),
     .   THK(*),THK_PART(*),XANEW(3,*),PM_STACK(20,*)
      INTEGER ID
      CHARACTER(LEN=NCHARTITLE) :: TITR
      TYPE (SURF_) :: IGRSURF1
      TYPE (SURF_) :: IGRSURF2
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER NDX, I, J, INRT, NELS, MT, JJ, JJJ, NELC,
     .   MG, L, NELTG,IE,IP,NM1,
     .  igtyp,ipgmat,igmat,isubstack
C     REAL
      my_real
     .   dxm, gapmx, gapmn, area, vol, dx,gaps1,gaps2, gapm, 
     .   gaptmp, gapscale,sx1,sy1,sz1,sx2,sy2,sz2,sx3,sy3,sz3,
     .   slsfac,gapinfs,gapinfm,gapsups,gapsupm,st
      INTEGER TAG(NUMNOD)
      INTEGER BITUNSET,BITGET,BITSET
      EXTERNAL bitunset,bitget,bitset
C--------------------------------------------------------------
C     CALCULATION OF SEGMENT STIFFNESSES
C     V16 : IN THE CASE WHERE ONE SEGMENT BELONGS BOTH
C           TO A BRICK AND A SHELL, THE SHELL STIFFNESS IS CHOSEN
C           OF THE SHELL UNLESS THE SHELL MATERIAL IS NULL.
C---------------------------------------------------------------
 
C---------------------------------------------------------------
      slsfac = one
C---------------------------------------------------------------
      ipgmat = 700
      igmat = 0
      DO i=1,numnod
        xanew(1,i)=x(1,i)
        xanew(2,i)=x(2,i)
        xanew(3,i)=x(3,i)
        tag(i)=0
      ENDDO
      dxm=0.
      ndx=0
      gapsolidmax=ep30
      gapmx=ep30
      gapmn=ep30
      gaps1=zero
      gaps2=zero
      IF(igap==2)THEN
        igap = 1
        gapscale = gapmin
        gapmin   = zero
      ELSE
        gapscale = one
      ENDIF  
C------------------------------------
C     GAP NODES SECONDS
C------------------------------------
      IF(igap>=1)THEN
       DO i=1,numnod
        wa(i)=zero
       ENDDO
       DO i=1,numelc
              mg=ixc(6,i)
              igtyp = igeo(11,mg)
                 ip = ipartc(i)
                 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
                   dx=half*thk_part(ip)
                 ELSEIF ( thk(i) /= zero .AND. iintthick == 0) THEN
                   dx=half*thk(i)
                 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR.igtyp == 52)THEN
                   dx=half*thk(i)
               ELSE
                 dx=half*geo(1,mg)
                 ENDIF
        wa(ixc(2,i))=max(wa(ixc(2,i)),dx)
        wa(ixc(3,i))=max(wa(ixc(3,i)),dx)
        wa(ixc(4,i))=max(wa(ixc(4,i)),dx)
        wa(ixc(5,i))=max(wa(ixc(5,i)),dx)
       ENDDO
       DO i=1,numeltg
        mg=ixtg(5,i)
        igtyp = igeo(11,mg)
        ip = iparttg(i)
        IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
          dx=half*thk_part(ip)
        ELSEIF (thk(numelc+i)/=zero .AND. iintthick==0) THEN
          dx=half*thk(numelc+i)
        ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
          dx=half*thk(numelc+i)
        ELSE
          dx=half*geo(1,mg)
        ENDIF
        wa(ixtg(2,i))=max(wa(ixtg(2,i)),dx)
        wa(ixtg(3,i))=max(wa(ixtg(3,i)),dx)
        wa(ixtg(4,i))=max(wa(ixtg(4,i)),dx)
       ENDDO
       DO i=1,numelt
        mg=ixt(4,i)
        dx=half*sqrt(geo(1,mg))
        wa(ixt(2,i))=max(wa(ixt(2,i)),dx)
        wa(ixt(3,i))=max(wa(ixt(3,i)),dx)
       ENDDO
       DO i=1,numelp
        mg=ixp(5,i)
        dx=0.5*sqrt(geo(1,mg))
        wa(ixp(2,i))=max(wa(ixp(2,i)),dx)
        wa(ixp(3,i))=max(wa(ixp(3,i)),dx)
       ENDDO
       DO i=1,nsn
        gap_s(i)=gapscale * wa(nsv(i))
        gaps1=max(gaps1,gap_s(i))
       ENDDO
      ENDIF
C
C calculation of the second surface. ---
      IF(intth > 0 ) THEN
       DO i = 1,nsn    
          areas(i) = zero
          DO j= knod2elc(nsv(i))+1,knod2elc(nsv(i)+1)
              ie = nod2elc(j)
              sx1 = x(1,ixc(4,ie)) - x(1,ixc(2,ie))
              sy1 = x(2,ixc(4,ie)) - x(2,ixc(2,ie))
              sz1 = x(3,ixc(4,ie)) - x(3,ixc(2,ie))
              sx2 = x(1,ixc(5,ie)) - x(1,ixc(3,ie))
              sy2 = x(2,ixc(5,ie)) - x(2,ixc(3,ie))
              sz2 = x(3,ixc(5,ie)) - x(3,ixc(3,ie))
              sx3  = sy1*sz2 - sz1*sy2
              sy3  = sz1*sx2 - sx1*sz2
              sz3  = sx1*sy2 - sy1*sx2
              areas(i) = areas(i) + one_over_8*sqrt(sx3*sx3+sy3*sy3+sz3*sz3)
          ENDDO
             ielec(i) = ixc(1,ie)
        ENDDO
      ENDIF
C
C------------------------------------
C     STIF NODES SECONDS
C------------------------------------
      IF(slsfac >= zero)THEN
       DO i=1,numelc
        mg=ixc(6,i)
        igtyp = igeo(11,mg)
                 ip = ipartc(i)
                 IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
                   dx=half*thk_part(ip)
                 ELSEIF ( thk(i) /= zero .AND. iintthick == 0) THEN
                   dx=half*thk(i)
                 ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52)THEN
                   dx=half*thk(i)
               ELSE
                 dx=half*geo(1,mg)
                 ENDIF
        wa(ixc(2,i))=max(wa(ixc(2,i)),dx)
        wa(ixc(3,i))=max(wa(ixc(3,i)),dx)
        wa(ixc(4,i))=max(wa(ixc(4,i)),dx)
        wa(ixc(5,i))=max(wa(ixc(5,i)),dx)
       ENDDO
       DO i=1,numeltg
        mg=ixtg(5,i)
        igtyp = igeo(11,mg)
        ip = iparttg(i)
        IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
          dx=half*thk_part(ip)
        ELSEIF (thk(numelc+i)/=zero .AND. iintthick==0) THEN
          dx=half*thk(numelc+i)
        ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52)THEN
          dx=half*thk(numelc+i)
        ELSE
          dx=half*geo(1,mg)
        ENDIF
        wa(ixtg(2,i))=max(wa(ixtg(2,i)),dx)
        wa(ixtg(3,i))=max(wa(ixtg(3,i)),dx)
        wa(ixtg(4,i))=max(wa(ixtg(4,i)),dx)
       ENDDO
       DO i=1,numelt
        mg=ixt(4,i)
        dx=half*sqrt(geo(1,mg))
        wa(ixt(2,i))=max(wa(ixt(2,i)),dx)
        wa(ixt(3,i))=max(wa(ixt(3,i)),dx)
       ENDDO
       DO i=1,numelp
        mg=ixp(5,i)
        dx=0.5*sqrt(geo(1,mg))
        wa(ixp(2,i))=max(wa(ixp(2,i)),dx)
        wa(ixp(3,i))=max(wa(ixp(3,i)),dx)
       ENDDO
       ! new interface buffer development : GAP_S not sized when IGAP=0 (cf bufintr)
       
c       DO I=1,NSN
c        GAP_S(I)=GAPSCALE * WA(NSV(I))
c        GAPS1=MAX(GAPS1,GAP_S(I))
c       ENDDO
      ENDIF
 
C------------------------------------
C     SHELL OR SOLID SURFACE
C------------------------------------
C------------------------------------
C     STIF FACES MAIN
C------------------------------------
 
      DO 500 i=1,nrt
      stf(i)=zero
      IF(intth > 0 ) ieles(i) = 0
      IF(slsfac<zero)THEN
        stf(i)=slsfac
      ENDIF
      gap_sh(i)=zero
      gapm  =zero
      inrt=i
      CALL i4gmx3(x,irect,inrt,gapmx)
C----------------------
      nm1=igrsurf1%NSEG
      IF(inrt <= nm1)THEN
        CALL i20nelts(x            ,irect(1,inrt),ixs  ,nint,nels          ,
     .                inrt         ,area         ,noint,0   ,igrsurf1%ELTYP,
     .                igrsurf1%ELEM)
      ELSE
        CALL i20nelts(x            ,irect(1,inrt),ixs  ,nint,nels          ,
     .                inrt-nm1     ,area         ,noint,0   ,igrsurf2%ELTYP,
     .                igrsurf2%ELEM)
      ENDIF
      IF(nels /= 0)THEN
        mt=ixs(1,nels)
        IF(mt>0)THEN
          DO jj=1,8
            jjj=ixs(jj+1,nels)
            xc(jj)=x(1,jjj)
            yc(jj)=x(2,jjj)
            zc(jj)=x(3,jjj)
          END DO
          CALL volint(vol)
          stf(i)=slsfac*area*area*pm(100,mt)/vol
        ELSE
          IF(nint>=0) THEN
             CALL ancmsg(msgid=95,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixs(nixs,nels),
     .                   c2='SOLID',
     .                   i3=i)
          ENDIF
          IF(nint<0) THEN 
             CALL ancmsg(msgid=96,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixs(nixs,nels),
     .                   c2='SOLID',
     .                   i3=i)
          ENDIF 
        ENDIF
        IF(igap/=0)THEN
          gap_sh(i)=min(vol/area,sqrt(area))/six
          gapsolidmax = min(gapsolidmax,vol/(area*four))
          gapmn=min(gapmn,half*gap_sh(i))
          gap_m(i)=zero
          tag(irect(1,inrt)) = 1
          tag(irect(2,inrt)) = 1
          tag(irect(3,inrt)) = 1
          tag(irect(4,inrt)) = 1
c          NBINFLG(IRECT(1,INRT))=BITUNSET(NBINFLG(IRECT(1,INRT)),7)
c          NBINFLG(IRECT(2,INRT))=BITUNSET(NBINFLG(IRECT(2,INRT)),7)
c          NBINFLG(IRECT(3,INRT))=BITUNSET(NBINFLG(IRECT(3,INRT)),7)
c          NBINFLG(IRECT(4,INRT))=BITUNSET(NBINFLG(IRECT(4,INRT)),7)
        ENDIF
        mbinflg(i)=bitset(mbinflg(i),8)
        GO TO 500
      ELSE
        IF(inrt <= nm1)THEN
          CALL ineltc(nelc ,neltg ,inrt    ,igrsurf1%ELTYP,igrsurf1%ELEM)
        ELSE
          CALL ineltc(nelc ,neltg ,inrt-nm1,igrsurf2%ELTYP,igrsurf2%ELEM)
        ENDIF
        IF(neltg/=0) THEN
          mt=ixtg(1,neltg)
          mg=ixtg(5,neltg)
          igtyp = igeo(11,mg)
          ip = iparttg(neltg)
          IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
            dx=thk_part(ip)*gapscale
          ELSEIF(thk(numelc+neltg)/=zero.AND.iintthick==0)THEN
            dx=thk(numelc+neltg)*gapscale
          ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52)THEN
            dx=thk(numelc+neltg)*gapscale
          ELSE
            dx=geo(1,mg)*gapscale
          ENDIF
          gapm=half*dx
          gaps2=max(gaps2,gapm)
          gapmn = min(gapmn,dx)
          dxm=dxm+dx
          ndx=ndx+1
          igmat = igeo(98,mg)
          IF(mt>0)THEN
            IF(igtyp == 11 .AND. igmat > 0) THEN
               IF ( thk(numelc+neltg) /=zero.AND.iintthick==0)THEN
                 stf(i)=slsfac*thk(numelc+neltg)*geo(ipgmat + 2 ,mg)
                ELSE
                 stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
               ENDIF
             ELSEIF(igtyp == 52 .OR. 
     .            ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN  
              isubstack = iworksh(3,numelc+neltg)
              st=pm_stack(2,isubstack) 
              stf(i)=slsfac*thk(numelc+neltg)*st
            ELSE
               IF ( thk(numelc+neltg) /=zero.AND.iintthick==0)THEN
                 stf(i)=slsfac*thk(numelc+neltg)*pm(20,mt)
               ELSEIF(igtyp == 17 .OR. igtyp == 51) THEN
                 stf(i)=slsfac*thk(numelc+neltg)*pm(20,mt)
                ELSE
                 stf(i)=slsfac*geo(1,mg)*pm(20,mt)
               ENDIF
           ENDIF
          ELSE
            IF(nint>=0) THEN
               CALL ancmsg(msgid=95,
     .                     msgtype=msgwarning,
     .                     anmode=aninfo_blind_2,
     .                     i1=id,
     .                     c1=titr,
     .                     i2=ixtg(nixtg,neltg),
     .                     c2='SHELL',
     .                     i3=i)
            END IF
            IF(nint<0) THEN
               CALL ancmsg(msgid=96,
     .                     msgtype=msgwarning,
     .                     anmode=aninfo_blind_2,
     .                     i1=id,
     .                     c1=titr,
     .                     i2=ixtg(nixtg,neltg),
     .                     c2='SHELL',
     .                     i3=i)
            END IF
          END IF
          IF(igap/=0) gap_m(i)=gapm
          mbinflg(i)=bitset(mbinflg(i),3)
          GO TO 500
        ELSEIF(nelc/=0) THEN
          mt=ixc(1,nelc)
          mg=ixc(6,nelc)
          igtyp = igeo(11,mg)
          ip = ipartc(nelc)
          igmat = igeo(99,mg)
          IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
            dx=thk_part(ip)*gapscale
          ELSEIF ( thk(nelc) /= zero .AND. iintthick == 0) THEN
            dx=thk(nelc)*gapscale
          ELSEIF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52) THEN
            dx=thk(nelc)*gapscale
          ELSE
            dx=geo(1,mg)*gapscale
          ENDIF
          gapm=half*dx
          gaps2=max(gaps2,gapm)
          gapmn = min(gapmn,dx)
          dxm=dxm+dx
          ndx=ndx+1
          IF(mt>0)THEN
            IF(igtyp == 11 .AND. igmat > 0) THEN
              IF ( thk(nelc) /= zero .AND. iintthick == 0) THEN 
                stf(i)=slsfac*thk(nelc)*geo(ipgmat + 2 ,mg)
              ELSE
                stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
              ENDIF            
            ELSEIF(igtyp==52 .OR.
     .            ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
              isubstack = iworksh(3,nelc)
               st=pm_stack(2,isubstack) 
               stf(i)=slsfac*thk(nelc)*st
            ELSE
              IF ( thk(nelc) /= zero .AND. iintthick == 0) THEN 
                stf(i)=slsfac*thk(nelc)*pm(20,mt)
              ELSEIF(igtyp == 17) THEN
                stf(i)=slsfac*thk(nelc)*pm(20,mt)
              ELSE
                stf(i)=slsfac*geo(1,mg)*pm(20,mt)
              ENDIF
           ENDIF 
          ELSE
            IF(nint>=0) THEN
               CALL ancmsg(msgid=95,
     .                     msgtype=msgwarning,
     .                     anmode=aninfo_blind_2,
     .                     i1=id,
     .                     c1=titr,
     .                     i2=ixc(nixc,nelc),
     .                     c2='SHELL',
     .                     i3=i)
            END IF
            IF(nint<0) THEN
               CALL ancmsg(msgid=96,
     .                     msgtype=msgwarning,
     .                     anmode=aninfo_blind_2,
     .                     i1=id,
     .                     c1=titr,
     .                     i2=ixc(nixc,nelc),
     .                     c2='SHELL',
     .                     i3=i)
            END IF
          END IF
          IF(igap/=0) gap_m(i)=gapm
          mbinflg(i)=bitset(mbinflg(i),4)
          GO TO 500
        END IF
      END IF
C----------------------
C     SEGMENT SURFACE
C----------------------
C----------------------
C     ELEMENTS SOLIDES
C----------------------
      CALL insol3(x,irect,ixs,nint,nels,inrt,
     .            area,noint,knod2els ,nod2els ,0,ixs10,
     .            ixs16,ixs20)
      IF(nels/=0) THEN
       gapm=zero
       mt=ixs(1,nels)
       IF(intth > 0 ) ieles(i) = nels
       IF(mt>0)THEN
        DO 100 jj=1,8
        jjj=ixs(jj+1,nels)
        xc(jj)=x(1,jjj)
        yc(jj)=x(2,jjj)
        zc(jj)=x(3,jjj)
  100   CONTINUE
        CALL volint(vol)
        stf(i)=slsfac*area*area*pm(100,mt)/vol
       ELSE
          IF(nint>=0) THEN
             CALL ancmsg(msgid=95,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixs(nixs,nels),
     .                   c2='SOLID',
     .                   i3=i)
          ENDIF
          IF(nint<0) THEN 
             CALL ancmsg(msgid=96,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixs(nixs,nels),
     .                   c2='SOLID',
     .                   i3=i)
          ENDIF
       ENDIF
       IF(igap/=0)THEN
          gap_sh(i)=min(vol/area,sqrt(area))/six
          gapsolidmax = min(gapsolidmax,vol/(area*four))
          gapmn=min(gapmn,half*gap_sh(i))
          gap_m(i)=zero
          tag(irect(1,inrt)) = 1
          tag(irect(2,inrt)) = 1
          tag(irect(3,inrt)) = 1
          tag(irect(4,inrt)) = 1
c          NBINFLG(IRECT(1,INRT))=BITUNSET(NBINFLG(IRECT(1,INRT)),7)
c          NBINFLG(IRECT(2,INRT))=BITUNSET(NBINFLG(IRECT(2,INRT)),7)
c          NBINFLG(IRECT(3,INRT))=BITUNSET(NBINFLG(IRECT(3,INRT)),7)
c          NBINFLG(IRECT(4,INRT))=BITUNSET(NBINFLG(IRECT(4,INRT)),7)
       ENDIF
       mbinflg(i)=bitset(mbinflg(i),8)
      ENDIF
C---------------------
C     ELEMENTS COQUES
C---------------------
      CALL incoq3(irect,ixc ,ixtg ,nint ,nelc     ,
     .            neltg,inrt,geo  ,pm   ,knod2elc ,
     .            knod2eltg ,nod2elc ,nod2eltg,thk,nty,igeo,
     .            pm_stack , iworksh)
      IF(neltg/=0) THEN
C
        mt=ixtg(1,neltg)
        mg=ixtg(5,neltg)
        igtyp = igeo(11,mg)
        ip = iparttg(neltg)
        IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
          dx=thk_part(ip)*gapscale
        ELSEIF ( thk(numelc+neltg) /= zero .AND. iintthick == 0)THEN
          dx=thk(numelc+neltg)*gapscale
        ELSEIF(igtyp ==17) THEN
          dx=thk(numelc+neltg)*gapscale
        ELSE
          dx=geo(1,mg)*gapscale
        ENDIF
        gapm=half*dx
        gaps2=max(gaps2,gapm)
        gapmn = min(gapmn,dx)
        dxm=dxm+dx
        ndx=ndx+1
        igmat = igeo(98,mg)
        IF(mt>0)THEN
          IF(igtyp == 11 .AND. igmat > 0) THEN
            IF ( thk(numelc+neltg) /= zero .AND. iintthick == 0) THEN
              stf(i)=slsfac*thk(numelc+neltg)*geo(ipgmat + 2 ,mg)
            ELSE
              stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
            ENDIF 
           ELSEIF(igtyp==52 .OR. 
     .          ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
              isubstack = iworksh(3,numelc+neltg)
              st=pm_stack(2,isubstack) 
              stf(i)=slsfac*thk(numelc+neltg)*st
          ELSE
            IF ( thk(numelc+neltg) /= zero .AND. iintthick == 0) THEN
              stf(i)=slsfac*thk(numelc+neltg)*pm(20,mt)
            ELSEIF(igtyp == 17) THEN
              stf(i)=slsfac*thk(numelc+neltg)*pm(20,mt)
            ELSE
              stf(i)=slsfac*geo(1,mg)*pm(20,mt)
            ENDIF
          ENDIF
        ELSE
           IF(nint>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixtg(nixtg,neltg),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
           IF(nint<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixtg(nixtg,neltg),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
        ENDIF
        IF(igap/=0) gap_m(i)=gapm
        mbinflg(i)=bitset(mbinflg(i),3)
      ELSEIF(nelc/=0) THEN
        mt=ixc(1,nelc)
        mg=ixc(6,nelc)
        igtyp = igeo(11,mg)
        ip = ipartc(nelc)
        IF ( thk_part(ip) /= zero .AND. iintthick == 0) THEN
          dx=thk_part(ip)*gapscale
        ELSEIF ( thk(nelc) /= zero .AND. iintthick == 0) THEN
          dx=thk(nelc)*gapscale
        ELSEIF(igtyp ==17) THEN
          dx=thk(nelc)*gapscale
        ELSE
          dx=geo(1,mg)*gapscale
        ENDIF
        gapm=half*dx
        gaps2=max(gaps2,gapm)
        gapmn = min(gapmn,dx)
        dxm=dxm+dx
        ndx=ndx+1
        igmat = igeo(98,mg)
        IF(mt>0)THEN
         IF(igtyp == 11 .AND. igmat > 0) THEN
          IF ( thk(nelc) /= zero .AND. iintthick == 0) THEN 
            stf(i)=slsfac*thk(nelc)*geo(ipgmat + 2 ,mg)
          ELSE
            stf(i)=slsfac*geo(1,mg)*geo(ipgmat + 2 ,mg)
          ENDIF   
         ELSEIF(igtyp==52 .OR. 
     .        ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
          isubstack = iworksh(3,nelc)
          st=pm_stack(2,isubstack) 
          stf(i)=slsfac*thk(nelc)*st
         ELSE 
          IF ( thk(nelc) /= zero .AND. iintthick == 0) THEN 
            stf(i)=slsfac*thk(nelc)*pm(20,mt)
          ELSEIF(igtyp ==17) THEN
            stf(i)=slsfac*thk(nelc)*pm(20,mt)
          ELSE
            stf(i)=slsfac*geo(1,mg)*pm(20,mt)
          ENDIF
         ENDIF 
        ELSE
           IF(nint>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixc(nixc,nelc),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
           IF(nint<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixc(nixc,nelc),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
        ENDIF
        IF(igap/=0) gap_m(i)=gapm
        mbinflg(i)=bitset(mbinflg(i),4)
      ENDIF
C
      IF(nels+nelc+neltg==0)THEN
 
C      In SPMD you need an element associated with the edge otherwise error
         IF(nint>0) THEN
            CALL ancmsg(msgid=481,
     .                  msgtype=msgerror,
     .                  anmode=aninfo_blind_2,
     .                  i1=id,
     .                  c1=titr,
     .                  i2=i)
         ENDIF
         IF(nint<0) THEN
            CALL ancmsg(msgid=482,
     .                  msgtype=msgerror,
     .                  anmode=aninfo_blind_2,
     .                  i1=id,
     .                  c1=titr,
     .                  i2=i)
         ENDIF
 
      ENDIF
  500 CONTINUE
C---------------------------
C     GAP 
C---------------------------
       gapmx=sqrt(gapmx)
       IF(igap==0)THEN
C GAP FIXE
         IF(gap<=zero)THEN
           IF(ndx/=0)THEN
             gap = dxm/ndx
             gap = min(half*gapmx,gap)
           ELSE
             gap = em01 * gapmx
           ENDIF
c           WRITE(IOUT,1300)GAP
         ENDIF
         gapmin = gap
         IF(inacti/=7.AND.gap>0.5*gapmx)THEN
          gaptmp = half*gapmx
          CALL ancmsg(msgid=94,
     .                msgtype=msgwarning,
     .                anmode=aninfo_blind_2,
     .                i1=id,
     .                c1=titr,
     .                r1=gap,
     .                r2=gaptmp)
         ENDIF
       ELSE
C GAP VARIABLE :
C    - GAPMIN CONTIENT ONE GAP MINIMUM UTILISE SI GAP_S(I)+GAP_M(J) < GAPMIN
C    - GAP CONTIENT LE SUP DE (GAP_S(I)+GAP_M(J),GAPMIN) 
         IF(gap<=zero)THEN
           IF(ndx/=0)THEN
             gapmin = gapmn
             gapmin = min(half*gapmx,gapmin)
           ELSE
c             GAPMIN = EM01 * GAPMX
             gapmin = min(gapmn,em01 * gapmx)
           ENDIF
c           WRITE(IOUT,1300)GAPMIN
         ELSE
           gapmin = gap
         ENDIF
C ABOVE VARIABLE GAPS
         gap = max(gaps1+gaps2,gapmin)
         gap=min(gap,gapmax)
         IF(inacti/=7.AND.gap>half*gapmx)THEN
          gaptmp = 0.5*gapmx
          CALL ancmsg(msgid=477,
     .                msgtype=msgwarning,
     .                anmode=aninfo_blind_2,
     .                i1=id,
     .                c1=titr,
     .                r1=gap)
         ENDIF
       ENDIF
C---------------------------------------------
C     SETTING TO ONE OF THE NODAL STIFFNESS MULTIPLIER
C---------------------------------------------
c     STFN is temporarily from 1 to NSN instead of 1 to NLN
      DO l=1,nsn
         stfn(l) = 1.
      ENDDO
C---------------------------------------------
C     LIMITATION OF THE GAP OF SOLIDS
C---------------------------------------------
      IF (igap/=0) THEN
        DO i = 1, nrt
          IF(gap_m(i) == zero)THEN
            gap_sh(i) = min(gapsolidmax,gap_sh(i))
            gap_sh(i) = max(gapsol,gap_sh(i))
c GAP_M for sorting
            gap_m(i)=gap_m(i)+two*gap_sh(i)
          ENDIF
        ENDDO
      ENDIF
C
C Calculation of the real gap to be used during the retri criterion
C
      gapshmax = zero
      IF (igap==0) THEN
        gapinf=gap
      ELSE
        gapinfs=ep30
        gapinfm=ep30
        gapsups = zero
        gapsupm = zero
        DO i = 1, nsn
          gapinfs = min(gapinfs,gap_s(i))
          gapsups = max(gapsups,gap_s(i))
        ENDDO
        DO i = 1, nrt
c          GAP_M(I)=GAP_M(I)+TWO*GAP_SH(I)
          gapinfm = min(gapinfm,gap_m(i))
          gapsupm = max(gapsupm,gap_m(i))
          gapshmax = max(gapshmax,gap_sh(i))
        ENDDO
        gapinf= max(gapinfs+gapinfm,gapmin)
        gap   = min(gapsups+gapsupm,gapmax)
      ENDIF  
 
      DO i=1,nln
        IF(tag(nlg(i)) == 1)nbinflg(i)=bitunset(nbinflg(i),7)
      ENDDO    
C------------
      RETURN
 1300 FORMAT(2x,'GAP MIN = ',1pg20.13)
      END
!||====================================================================
!||    i20sti3e               ../starter/source/interfaces/inter3d1/i20sti3.f
!||--- called by ------------------------------------------------------
!||    i20ini3                ../starter/source/interfaces/inter3d1/i20ini3.F
!||--- calls      -----------------------------------------------------
!||    ancmsg                 ../starter/source/output/message/message.f
!||    get_u_geo              ../starter/source/user_interface/uaccess.F
!||    i11coq                 ../starter/source/interfaces/inter3d1/i11coq.F
!||    i11fil                 ../starter/source/interfaces/inter3d1/i11coq.F
!||    i11gmx3                ../starter/source/interfaces/inter3d1/i11gmx3.F
!||    i11sol                 ../starter/source/interfaces/inter3d1/i11sol.f
!||    my_exit                ../starter/source/output/analyse/analyse.c
!||    volint                 ../starter/source/interfaces/inter3d1/volint.F
!||--- uses       -----------------------------------------------------
!||    message_mod            ../starter/share/message_module/message_mod.F
!||====================================================================
      SUBROUTINE i20sti3e(
     1 X     ,IXLIN ,STF   ,IXS   ,PM    ,
     2 GEO   ,NRT   ,IXC   ,NINTR ,SLSFAC,
     3 NTY   ,GAPMAX,NOINT ,GAP_SM,
     4 MS    ,IXTG  ,IXT   ,IXP   ,IXR   ,
     5 IGAP  ,GAPMIN,GAP0  ,GAPINF,NSNE  ,
     6 IPARTC,IPARTTG,THK  ,THK_PART,IXS10,
     7 ID    ,TITR  ,KXX   ,IXX   ,IGEO,
     8 NOD2EL1D,KNOD2EL1D,KNOD2ELS,KNOD2ELC,KNOD2ELTG,
     9 NOD2ELS,NOD2ELC,NOD2ELTG   ,LELX ,PM_STACK,IWORKSH)
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE message_mod
      USE names_and_titles_mod , ONLY : nchartitle
      use element_mod , only : nixs,nixc,nixtg,nixt,nixp,nixr
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "units_c.inc"
#include      "param_c.inc"
#include      "com04_c.inc"
#include      "scr08_c.inc"
#include      "scr23_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER NRT, NINTR, NTY, NOINT,IGAP,NSNE
C     REAL
      my_real
     .   SLSFAC, GAPMAX,GAPMIN,GAP0
      INTEGER IXLIN(2,*), IXS(NIXS,*), IXC(NIXC,*),
     .   IXTG(NIXTG,*),IXT(NIXT,*),IXP(NIXP,*),IXR(NIXR,*),
     .   IPARTC(*), IPARTTG(*),IXS10(*),KXX(NIXX,*),IXX(*),
     .   igeo(npropgi,*),
     .   knod2els(*), knod2elc(*), knod2eltg(*), 
     .   nod2els(*), nod2elc(*), nod2eltg(*),
     .   nod2el1d(*),knod2el1d(*),iworksh(3,*)
C     REAL
      my_real
     .   x(3,*), stf(*), pm(npropm,*), geo(npropg,*),
     .   ms(*),gap_sm(*),xl2, gapinf,thk(*),thk_part(*),lelx(*),
     .   pm_stack(20,*)
      INTEGER ID
      CHARACTER(LEN=NCHARTITLE) :: TITR
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER NDX, I, INRT, NELS, MT, JJ, JJJ, NELC,
     .   MG, NELTG,NELT,NELP,NELR,
     .   IGTYP,NELX,IPGMAT,IGMAT,ISUBSTACK
C     REAL
      my_real
     .   dxm, gapmx, gapmn, area, vol, dx,gap1,gaps1,gaptmp,
     .   get_u_geo
C-----------------------------------------------
      EXTERNAL get_u_geo
C--------------------------------------------------------------
C     CALCULATION OF SEGMENT STIFFNESSES
C     V16 : IN THE CASE WHERE ONE SEGMENT BELONGS BOTH
C           TO A BRICK AND A SHELL, THE SHELL STIFFNESS IS CHOSEN
C           OF THE SHELL UNLESS THE SHELL MATERIAL IS NULL.
C---------------------------------------------------------------
      dxm=zero
      ndx=0      
      gaps1=zero
      gapmx=ep30
      gapmn=ep30
 
C
      DO 500 i=1,nrt
      stf(i)=zero
      gap_sm(i)=zero
      inrt=i
      CALL i11gmx3(x,ixlin,inrt,gapmx,xl2)
C----------------------
C     ELEMENTS SOLIDES
C----------------------
      CALL i11sol(x,ixlin,ixs,nintr,nels,inrt,
     .            area,noint,knod2els,nod2els,ixs10)
      IF(nels/=0) THEN
       mt=ixs(1,nels)
       IF(mt>0)THEN
        DO 100 jj=1,8
        jjj=ixs(jj+1,nels)
        xc(jj)=x(1,jjj)
        yc(jj)=x(2,jjj)
        zc(jj)=x(3,jjj)
  100   CONTINUE
        CALL volint(vol)
        IF(xl2>0.0)THEN
          stf(i)=slsfac*vol*pm(100,mt)/xl2
        ELSE
          stf(i)=zero
        ENDIF
       ELSE
C        IF(NINTR>=0)WRITE (IOUT,1500) IXS(11,NELS),I, NOINT
          IF(nintr>=0) THEN
             CALL ancmsg(msgid=95,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixs(nixs,nels),
     .                   c2='SOLID',
     .                   i3=i)
          ENDIF
C        IF(NINTR<0)WRITE (IOUT,1600) IXS(11,NELS),I, NOINT
C        IWARN=IWARN+1
          IF(nintr<0) THEN
             CALL ancmsg(msgid=96,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixs(nixs,nels),
     .                   c2='SOLID',
     .                   i3=i)
          ENDIF
       ENDIF
      ENDIF
C---------------------
C     ELEMENTS COQUES
C---------------------
      CALL i11coq(ixlin,ixc ,ixtg,nintr,nelc ,
     .            neltg,inrt,geo,pm,thk,igeo,
     .            knod2elc,knod2eltg,nod2elc,nod2eltg,
     .            pm_stack, iworksh )
      IF(neltg/=0) THEN
C
        mt=ixtg(1,neltg)
        mg=ixtg(5,neltg)
        igtyp = igeo(11,mg)
        dx=geo(1,mg)
        IF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp ==52)
     .                                   dx = thk(numelc + neltg)
        gap_sm(i)=half*dx
        gaps1=max(gaps1,gap_sm(i))
        gapmn = min(gapmn,dx)
        dxm=dxm+dx
        ndx=ndx+1
        igmat = igeo(98,mg)
        ipgmat=700
        IF(mt>0)THEN
         IF(igtyp == 11 .AND. igmat > 0)THEN
          stf(i)=slsfac*dx*geo(ipgmat + 2 ,mg)
         ELSEIF(igtyp==52 .OR. 
     .        ((igtyp == 17 .OR. igtyp == 51 ) .AND. igmat > 0))THEN
            isubstack = iworksh(3,neltg + numelc)
            stf(i)=slsfac*dx*pm_stack(2,isubstack)
         ELSE
           stf(i)=slsfac*dx*pm(20,mt)
         ENDIF  
        ELSE
C         IF(NINTR>=0)WRITE (IOUT,1700) IXTG(NIXTG,NELTG),I, NOINT
           IF(nintr>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixtg(nixtg,neltg),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
C         IF(NINTR<0)WRITE (IOUT,1800) IXTG(NIXTG,NELTG),I, NOINT
C         IWARN=IWARN+1
           IF(nintr<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixtg(nixtg,neltg),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
        ENDIF
      ELSEIF(nelc/=0) THEN
C
        mt=ixc(1,nelc)
        mg=ixc(6,nelc)
        igtyp = igeo(11,mg)
        dx=geo(1,mg)
        IF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp ==52)
     .                                         dx = thk(nelc)
        gap_sm(i)=half*dx
        gaps1=max(gaps1,gap_sm(i))
        gapmn = min(gapmn,dx)
        dxm=dxm+dx
        ndx=ndx+1
        igmat = igeo(98,mg)
        IF(mt>0)THEN
           IF(igtyp == 11 .AND. igmat > 0) THEN
            stf(i)=slsfac*dx*geo(ipgmat + 2 ,mg)
           ELSEIF(igtyp ==52 .OR. 
     .          ((igtyp == 17 .OR. igtyp == 51) .AND. igmat > 0))THEN
            isubstack = iworksh(3,nelc)
            stf(i)=slsfac*dx*pm_stack(2,isubstack) 
          ELSE
           stf(i)=slsfac*dx*pm(20,mt)
          ENDIF   
        ELSE
C           IF(NINTR>=0)WRITE (IOUT,1700) IXC(7,NELC),I, NOINT
           IF(nintr>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixc(nixc,nelc),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
C           IF(NINTR<0)WRITE (IOUT,1800) IXC(7,NELC),I, NOINT
C           IWARN=IWARN+1
           IF(nintr<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixc(nixc,nelc),
     .                    c2='SHELL',
     .                    i3=i)
           ENDIF
        ENDIF
      ENDIF
C---------------------
C     ELEMENTS TIGE POUTRE RESSORT
C---------------------
      CALL i11fil(ixlin,ixt,ixp,ixr,nintr,nelt ,
     .            nelp,nelr,nelx,inrt,nod2el1d,
     .            knod2el1d,kxx,ixx)
 
      IF(nelt/=0) THEN
C
        mt=ixt(1,nelt)
        mg=ixt(4,nelt)
        dx=sqrt(geo(1,mg))
        IF(igtyp == 17 .OR. igtyp == 51 .OR. igtyp == 52 )
     .                                  dx = sqrt(thk(numelc + nelt))
        gap_sm(i)=max(gap_sm(i),half*dx)
        gaps1=max(gaps1,gap_sm(i))
        gapmn = min(gapmn,dx)
        dxm=dxm+dx
        ndx=ndx+1
        igmat = igeo(98,mg)
        IF(mt>0)THEN
         IF(igtyp == 11 .AND. igmat > 0) THEN
          stf(i)=slsfac*dx*geo(ipgmat + 2 ,mg)
         ELSEIF((igtyp == 17 .OR. igtyp == 17) .AND. igmat > 0) THEN
            isubstack = iworksh(3,numelc + nelt)
             stf(i)=slsfac*dx*pm_stack(2,isubstack) 
         ELSE
            stf(i)=slsfac*dx*pm(20,mt)
         ENDIF
        ELSE
C         IF(NINTR>=0)WRITE (IOUT,1700) IXT(NIXT,NELT),I, NOINT
           IF(nintr>=0) THEN 
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixt(nixt,nelt),
     .                    c2='TRUSS',
     .                    i3=i)
           ENDIF
C         IF(NINTR<0)WRITE (IOUT,1800) IXT(NIXT,NELT),I, NOINT
C         IWARN=IWARN+1
           IF(nintr<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixt(nixt,nelt),
     .                    c2='TRUSS',
     .                    i3=i)
           ENDIF
        ENDIF
      ELSEIF(nelp/=0) THEN
C
        mt=ixp(1,nelp)
        mg=ixp(5,nelp)
        dx=sqrt(geo(1,mg))
        gap_sm(i)=max(gap_sm(i),half*dx)
        gaps1=max(gaps1,gap_sm(i))
        gapmn = min(gapmn,dx)
        dxm=dxm+dx
        ndx=ndx+1
        IF(mt>0)THEN
         stf(i)=slsfac*dx*pm(20,mt)
        ELSE
C         IF(NINTR>=0)WRITE (IOUT,1700) IXP(NIXP,NELP),I, NOINT
           IF(nintr>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixp(nixp,nelp),
     .                    c2='BEAM',
     .                    i3=i)
           ENDIF
C         IF(NINTR<0)WRITE (IOUT,1800) IXP(NIXP,NELP),I, NOINT
C         IWARN=IWARN+1
          IF(nintr<0) THEN
             CALL ancmsg(msgid=96,
     .                   msgtype=msgwarning,
     .                   anmode=aninfo_blind_2,
     .                   i1=id,
     .                   c1=titr,
     .                   i2=ixp(nixp,nelp),
     .                   c2='BEAM',
     .                   i3=i)
          ENDIF
        ENDIF
      ELSEIF(nelr/=0) THEN
C
        mg=ixr(1,nelr)
        mt = ixr(5,nelr)
        IF(mg>0)THEN
         igtyp=nint(geo(12,mg))
         IF(igtyp==4.OR.igtyp==12)THEN
           stf(i)=slsfac*geo(2,mg)
         ELSEIF(igtyp==8.OR.igtyp==13)THEN
           stf(i)=slsfac*max(geo(3,mg),geo(10,mg),geo(15,mg))
         ELSEIF(igtyp == 23)THEN
           stf(i)=slsfac*max(pm(191,mt),pm(192,mt),pm(193,mt))
         ELSEIF(igtyp==25)THEN
           stf(i)=slsfac*geo(10,mg)
         ELSEIF(igtyp>=29)THEN
           stf(i)=slsfac*geo(3,mg)
         ELSE
            WRITE(6,'(A)') 'INTERNAL ERROR 987'
            CALL my_exit(2)
C           STOP 987
         ENDIF
        ELSE
C         IF(NINTR>=0)WRITE (IOUT,1700) IXR(NIXR,NELR),I, NOINT
           IF(nintr>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixr(nixr,nelr),
     .                    c2='SPRING',
     .                    i3=i)
           ENDIF
C         IF(NINTR<0)WRITE (IOUT,1800) IXR(NIXR,NELR),I, NOINT
C         IWARN=IWARN+1
           IF(nintr<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=ixr(nixr,nelr),
     .                    c2='SPRING',
     .                    i3=i)
           ENDIF
        ENDIF
      ELSEIF(nelx/=0) THEN
C
        mg=kxx(2,nelx)
        IF(mg>0)THEN
         stf(i)=slsfac*get_u_geo(4,mg)*(kxx(3,nelx)-1)/lelx(nelx)
        ELSE
           IF(nintr>=0) THEN
              CALL ancmsg(msgid=95,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=kxx(nixx,nelx),
     .                    c2='XELEM',
     .                    i3=i)
           ENDIF
           IF(nintr<0) THEN
              CALL ancmsg(msgid=96,
     .                    msgtype=msgwarning,
     .                    anmode=aninfo_blind_2,
     .                    i1=id,
     .                    c1=titr,
     .                    i2=kxx(nixx,nelx),
     .                    c2='XELEM',
     .                    i3=i)
           ENDIF
        ENDIF
      ENDIF
C---------------------------
      IF(nels+nelc+neltg+nelt+nelp+nelr+numelx==0.)THEN
C In SPMD you need an element associated with the edge otherwise error
        IF(nintr>0) THEN
          CALL ancmsg(msgid=481,
     .                msgtype=msgerror,
     .                anmode=aninfo_blind_2,
     .                i1=id,
     .                c1=titr,
     .                i2=i)
        ENDIF
        IF(nintr<0) THEN
          CALL ancmsg(msgid=482,
     .                msgtype=msgerror,
     .                anmode=aninfo_blind_2,
     .                i1=id,
     .                c1=titr,
     .                i2=i)
        ENDIF
      ENDIF
 
  500 CONTINUE
C---------------------------
C     GAP 
C---------------------------
      gapmx=sqrt(gapmx)
      IF(igap==0)THEN
C---------------------------
C     GAP  FIXE
C---------------------------
        IF(gap0>zero)THEN
          gap1 = gap0
        ELSE
          IF(ndx/=0)THEN
           gap1 = min(half*gapmx,dxm/ndx)
          ELSE
           gap1 = em01* gapmx
          ENDIF
          IF(nintr<0)WRITE(iout,1300)half*(gapmin+gap1)
        ENDIF
C
        IF(nintr<0) gap1 = half*(gapmin+gap1)
        gapmin = gap1
        gapmax = gap1
C
        IF(gap1>half*gapmx)THEN
C          WRITE(IOUT,1400)GAP1,0.5*GAPMX
C          IWARN=IWARN+1
          gaptmp = half*gapmx
          CALL ancmsg(msgid=94,
     .                msgtype=msgwarning,
     .                anmode=aninfo_blind_2,
     .                i1=id,
     .                c1=titr,
     .                r1=gap1,
     .                r2=gaptmp)
        ENDIF
      ELSE
C---------------------------
C GAP VARIABLE
C---------------------------
        IF(gap0>zero)THEN
          gap1 = gap0
        ELSE
          IF(ndx/=0)THEN
           gap1 = min(half*gapmx,gapmn)
          ELSE
           gap1 = em01 * gapmx
          ENDIF
          IF(nintr<0)WRITE(iout,1300)half*(gapmin+gap1)
        ENDIF
C MINI GAP AND ABOVE VARIABLE GAPS
        IF(nintr>0)THEN
          gapmin = gap1
          gapmax = gaps1
        ELSE
          gapmin = half*(gapmin+gap1)
          gapmax = max(gapmax+gaps1,gapmin)
        ENDIF
C
        IF(gapmax>half*gapmx)THEN
           gaptmp = half*gapmx
           CALL ancmsg(msgid=94,
     .                 msgtype=msgwarning,
     .                 anmode=aninfo_blind_2,
     .                 i1=id,
     .                 c1=titr,
     .                 r1=gapmax,
     .                 r2=gaptmp)
        ENDIF
      ENDIF
C---------------------------
C     STIF GLOBAL
C---------------------------
      IF(slsfac<zero)THEN
        DO i=1,nrt
          stf(i)=-slsfac
        ENDDO
      ENDIF
C---------------------------------------------
C
C Calculation of the real gap to be used during the retri criterion
C
      IF (igap==0) THEN
        gapinf=gapmax/two ! multiplied by 2 in i20ini3
      ELSE
        DO i = 1, nrt
          gapinf = min(gapinf,gap_sm(i))
        ENDDO
      ENDIF      
      RETURN
 
 1300 FORMAT(2x,'COMPUTED GAP = ',1pg20.13)
 
      END
!||====================================================================
!||    i20nlg    ../starter/source/interfaces/inter3d1/i20sti3.F
!||--- called by ------------------------------------------------------
!||    i20ini3   ../starter/source/interfaces/inter3d1/i20ini3.F
!||====================================================================
      SUBROUTINE i20nlg(NLN  ,NRTM,NSN  ,NLINS ,NLINM ,
     .                  NLG  ,IRECT,NSV ,IXLINS,IXLINM,
     .                  NMN  ,NSNE ,NMNE,MSR   ,NSVE  ,
     .                  MSRE ,STFA ,DXANC,XANEW,X     ,
     .                  PENIA,ALPHAK)
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "com04_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER NLN,NRTM, NSN,NLINS ,NLINM ,NMN  ,NSNE ,NMNE     
      INTEGER IRECT(4,NRTM), NSV(NSN),IXLINS(2,NLINS),IXLINM(2,NLINM),
     .        msr(nmn),nsve(nsne),msre(nmne),nlg(nln)
      my_real
     .   stfa(*),dxanc(3,*),xanew(3,*),x(3,*),penia(5,*),alphak(3,*)
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I,J,K
      INTEGER TAG(NUMNOD)
      my_real
     .   AAA,STIF(NLN)
 
      DO i=1,nln
        j = nlg(i)
        tag(j)=i
      ENDDO
 
      DO k=1,nsn
        nsv(k)=tag(nsv(k))
      ENDDO
      DO k=1,nmn
        msr(k)=tag(msr(k))
      ENDDO
      DO k=1,nsne
        nsve(k)=tag(nsve(k))
      ENDDO
      DO k=1,nmne
        msre(k)=tag(msre(k))
      ENDDO
 
      DO k=1,nrtm
        irect(1,k)=tag(irect(1,k))
        irect(2,k)=tag(irect(2,k))
        irect(3,k)=tag(irect(3,k))
        irect(4,k)=tag(irect(4,k))
      ENDDO
      DO k=1,nlins
        ixlins(1,k)=tag(ixlins(1,k))
        ixlins(2,k)=tag(ixlins(2,k))
      ENDDO
      DO k=1,nlinm
        ixlinm(1,k)=tag(ixlinm(1,k))
        ixlinm(2,k)=tag(ixlinm(2,k))
      ENDDO
 
      DO i=1,nln
        stif(i) = one
        alphak(1,i)  = one
        alphak(2,i)  = one
        alphak(3,i)  = one
      ENDDO
 
      DO i=1,nsn
c In input stfa (1: nln) is STFN (1: nsn) Possibly put to zero if pene initial
        j = nsv(i)
        stif(j) = stfa(i)
      ENDDO
 
      DO i=1,nln
        stfa(i) = stif(i)
      ENDDO
 
c STFA sera recalcule in I20STIFN /inter3d1/i20stifn.F
      DO i=1,nln
        dxanc(1,i) = xanew(1,nlg(i))-x(1,nlg(i))
        dxanc(2,i) = xanew(2,nlg(i))-x(2,nlg(i))
        dxanc(3,i) = xanew(3,nlg(i))-x(3,nlg(i))
        penia(4,i) = sqrt(dxanc(1,i)*dxanc(1,i)
     +                   +dxanc(2,i)*dxanc(2,i)
     +                   +dxanc(3,i)*dxanc(3,i))
        penia(5,i) = penia(4,i)
        aaa = one/max(penia(4,i),em20)
        penia(1,i) = dxanc(1,i)*aaa
        penia(2,i) = dxanc(2,i)*aaa
        penia(3,i) = dxanc(3,i)*aaa
      ENDDO
 
      RETURN
      END