s6cinit3.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/.
!||====================================================================
!||    s6cinit3               ../starter/source/elements/thickshell/solide6c/s6cinit3.F
!||--- called by ------------------------------------------------------
!||    initia                 ../starter/source/elements/initia/initia.F
!||--- calls      -----------------------------------------------------
!||    ancmsg                 ../starter/source/output/message/message.F
!||    atheri                 ../starter/source/ale/atheri.F
!||    dtmain                 ../starter/source/materials/time_step/dtmain.F
!||    failini                ../starter/source/elements/solid/solide/failini.F
!||    fretitl2               ../starter/source/starter/freform.F
!||    matini                 ../starter/source/materials/mat_share/matini.F
!||    s6ccoor3               ../starter/source/elements/thickshell/solide6c/s6ccoor3.F
!||    s6cderi3               ../starter/source/elements/thickshell/solide6c/s6cderi3.F
!||    s6mass3                ../starter/source/elements/thickshell/solide6c/s6mass3.F
!||    sbulk3                 ../starter/source/elements/solid/solide/sbulk3.F
!||    scmorth3               ../starter/source/elements/thickshell/solidec/scmorth3.F
!||    sczero3                ../starter/source/elements/thickshell/solidec/scinit3.F
!||    sdlensh3n              ../starter/source/elements/thickshell/solide6c/s6cinit3.F
!||    sigin20b               ../starter/source/elements/solid/solide20/s20mass3.F
!||    svalue0                ../starter/source/elements/thickshell/solidec/scinit3.F
!||--- uses       -----------------------------------------------------
!||    defaults_mod           ../starter/source/modules/defaults_mod.F90
!||    detonators_mod         ../starter/share/modules1/detonators_mod.F
!||    message_mod            ../starter/share/message_module/message_mod.F
!||====================================================================
      SUBROUTINE s6cinit3(
     .                    ELBUF_STR,MAS     ,IXS     ,PM       ,X     ,
     .                    DETONATORS,GEO     ,VEUL    ,ALE_CONNECTIVITY,IPARG ,
     .                    DTELEM   ,SIGI    ,NEL     ,SKEW     ,IGEO  ,
     .                    STIFN    ,PARTSAV ,V       ,IPARTS   ,MSS   ,
     .                    IPART    ,GLOB_THERM,
     .                    SIGSP    ,NSIGI   ,IPM     ,IUSER    ,NSIGS ,
     .                    VOLNOD   ,BVOLNOD ,VNS     ,BNS      ,PTSOL ,
     .                    BUFMAT   ,MCP     ,MCPS    ,MCPSX    ,TEMP  ,
     .                    NPF      ,TF      ,STRSGLOB,STRAGLOB ,MSSA  ,
     .                    ORTHOGLOB,FAIL_INI,ILOADP  ,FACLOAD  ,RNOISE,
     .                    PERTURB  ,MAT_PARAM,DEFAULTS_SOLID)
C-----------------------------------------------
C   D e s c r i p t i o n
C   Initialization of thick shell PA6 element
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE elbufdef_mod            
      USE message_mod
      USE detonators_mod      
      USE ale_connectivity_mod
      USE matparam_def_mod
      USE defaults_mod
      USE names_and_titles_mod , ONLY : nchartitle
      use glob_therm_mod
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   C o m m o n   B l o c k s
C-----------------------------------------------
#include      "com04_c.inc"
#include      "param_c.inc"
#include      "scr12_c.inc"
#include      "scr17_c.inc"
#include      "scry_c.inc"
#include      "vect01_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER IXS(NIXS,*), IPARG(*),IPARTS(*), 
     .    NEL, IPART(LIPART1,*),PERTURB(NPERTURB), 
     .    IPM(NPROPMI,*), PTSOL(*), NSIGI, IUSER, NSIGS, NPF(*)
      INTEGER IGEO(NPROPGI,*),STRSGLOB(*),STRAGLOB(*),ORTHOGLOB(*),
     .    FAIL_INI(*)
      my_real
     .   MAS(*), PM(NPROPM,*), X(*), GEO(NPROPG,*),
     .   VEUL(LVEUL,*), DTELEM(*),SIGI(NSIGS,*),SKEW(LSKEW,*),STIFN(*),
     .   PARTSAV(20,*), V(*), MSS(8,*),SIGSP(NSIGI,*),
     .    VOLNOD(*), BVOLNOD(*), VNS(8,*), BNS(8,*),BUFMAT(*),MCP(*),
     .    MCPS(8,*), MCPSX(12,*),TEMP(*), TF(*), MSSA(*),RNOISE(NPERTURB,*)
      TYPE(elbuf_struct_), TARGET :: ELBUF_STR
      INTEGER,INTENT(IN) :: ILOADP(SIZLOADP,*)
      my_real,INTENT(IN) :: FACLOAD(LFACLOAD,*)
      TYPE(detonators_struct_) :: DETONATORS
      TYPE(t_ale_connectivity), INTENT(INOUT) :: ALE_CONNECTIVITY
      TYPE (MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(INOUT) :: MAT_PARAM
      TYPE(SOLID_DEFAULTS_), INTENT(IN) :: DEFAULTS_SOLID
      type (glob_therm_) ,intent(in)      :: glob_therm
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I,NF1,IBID,IGTYP,IREP,IP,ILAY,NLAY,NUVAR,NCC,JHBE,
     .   nuvarr,idef,ipang,ipthk,ippos,ipmat,ig,im,mtn0,nlymax,
     .   ipid1,nptr,npts,nptt,l_pla,l_sigb,imas_ds
      INTEGER MAT(MVSIZ), PID(MVSIZ), NGL(MVSIZ), MAT0(MVSIZ)
      INTEGER IX1(MVSIZ), IX2(MVSIZ), IX3(MVSIZ), IX4(MVSIZ),
     .        ix5(mvsiz), ix6(mvsiz)
      my_real
     .     x1(mvsiz), x2(mvsiz), x3(mvsiz), x4(mvsiz), x5(mvsiz), x6(mvsiz), 
     .     y1(mvsiz), y2(mvsiz), y3(mvsiz), y4(mvsiz), y5(mvsiz), y6(mvsiz), 
     .     z1(mvsiz), z2(mvsiz), z3(mvsiz), z4(mvsiz), z5(mvsiz), z6(mvsiz)
      CHARACTER(LEN=NCHARTITLE)::TITR1
      my_real
     .   bid, fv, sti, zi,wi
      my_real
     .   v8loc(51,mvsiz),volu(mvsiz),dtx(mvsiz),vzl(mvsiz),vzq(mvsiz),
     .   rx(mvsiz) ,ry(mvsiz) ,rz(mvsiz) ,sx(mvsiz) ,
     .   sy(mvsiz) ,sz(mvsiz) ,tx(mvsiz) ,ty(mvsiz) ,tz(mvsiz) ,
     .   e1x(mvsiz),e1y(mvsiz),e1z(mvsiz),
     .   e2x(mvsiz),e2y(mvsiz),e2z(mvsiz),
     .   e3x(mvsiz),e3y(mvsiz),e3z(mvsiz),
     .   f1x(mvsiz) ,f1y(mvsiz) ,f1z(mvsiz) ,llsh(mvsiz) ,
     .   f2x(mvsiz) ,f2y(mvsiz) ,f2z(mvsiz) ,rhocp(mvsiz),temp0(mvsiz), deltax(mvsiz), aire(mvsiz)
      my_real :: tempel(nel)
C-----------------------------------------------
      TYPE(g_bufel_) ,POINTER :: GBUF     
      TYPE(BUF_LAY_) ,POINTER :: BUFLY
      TYPE(L_BUFEL_) ,POINTER :: LBUF     
      TYPE(BUF_MAT_) ,POINTER :: MBUF
C-----------------------------------------------
      my_real
     .  W_GAUSS(9,9),A_GAUSS(9,9),ANGLE(MVSIZ),DTX0(MVSIZ)
      DATA W_GAUSS /
     1 2.               ,0.               ,0.               ,
     1 0.               ,0.               ,0.               ,
     1 0.               ,0.               ,0.               ,
     2 1.               ,1.               ,0.               ,
     2 0.               ,0.               ,0.               ,
     2 0.               ,0.               ,0.               ,
     3 0.555555555555556,0.888888888888889,0.555555555555556,
     3 0.               ,0.               ,0.               ,
     3 0.               ,0.               ,0.               ,
     4 0.347854845137454,0.652145154862546,0.652145154862546,
     4 0.347854845137454,0.               ,0.               ,
     4 0.               ,0.               ,0.               ,
     5 0.236926885056189,0.478628670499366,0.568888888888889,
     5 0.478628670499366,0.236926885056189,0.               ,
     5 0.               ,0.               ,0.               ,
     6 0.171324492379170,0.360761573048139,0.467913934572691,
     6 0.467913934572691,0.360761573048139,0.171324492379170,
     6 0.               ,0.               ,0.               ,
     7 0.129484966168870,0.279705391489277,0.381830050505119,
     7 0.417959183673469,0.381830050505119,0.279705391489277,
     7 0.129484966168870,0.               ,0.               ,
     8 0.101228536290376,0.222381034453374,0.313706645877887,
     8 0.362683783378362,0.362683783378362,0.313706645877887,
     8 0.222381034453374,0.101228536290376,0.               ,
     9 0.081274388361574,0.180648160694857,0.260610696402935,
     9 0.312347077040003,0.330239355001260,0.312347077040003,
     9 0.260610696402935,0.180648160694857,0.081274388361574/
      DATA a_gauss /
     1 0.               ,0.               ,0.               ,
     1 0.               ,0.               ,0.               ,
     1 0.               ,0.               ,0.               ,
     2 -.577350269189626,0.577350269189626,0.               ,
     2 0.               ,0.               ,0.               ,
     2 0.               ,0.               ,0.               ,
     3 -.774596669241483,0.               ,0.774596669241483,
     3 0.               ,0.               ,0.               ,
     3 0.               ,0.               ,0.               ,
     4 -.861136311594053,-.339981043584856,0.339981043584856,
     4 0.861136311594053,0.               ,0.               ,
     4 0.               ,0.               ,0.               ,
     5 -.906179845938664,-.538469310105683,0.               ,
     5 0.538469310105683,0.906179845938664,0.               ,
     5 0.               ,0.               ,0.               ,
     6 -.932469514203152,-.661209386466265,-.238619186083197,
     6 0.238619186083197,0.661209386466265,0.932469514203152,
     6 0.               ,0.               ,0.               ,
     7 -.949107912342759,-.741531185599394,-.405845151377397,
     7 0.               ,0.405845151377397,0.741531185599394,
     7 0.949107912342759,0.               ,0.               ,
     8 -.960289856497536,-.796666477413627,-.525532409916329,
     8 -.183434642495650,0.183434642495650,0.525532409916329,
     8 0.796666477413627,0.960289856497536,0.               ,
     9 -.968160239507626,-.836031107326636,-.613371432700590,
     9 -.324253423403809,0.               ,0.324253423403809,
     9 0.613371432700590,0.836031107326636,0.968160239507626/
C-----------------------------------------------
C   S o u r c e  L i n e s
C=======================================================================
      gbuf => elbuf_str%GBUF
      lbuf  => elbuf_str%BUFLY(1)%LBUF(1,1,1)
      mbuf  => elbuf_str%BUFLY(1)%MAT(1,1,1)
      bufly => elbuf_str%BUFLY(1)
      nptr  =  elbuf_str%NPTR
      npts  =  elbuf_str%NPTS
      nptt  =  elbuf_str%NPTT
      nlay  =  elbuf_str%NLAY
C
      jhbe  = iparg(23)
      irep  = iparg(35)
      igtyp = iparg(38)
      nf1=nft+1
      idef =0
      ibid = 0
      bid   = zero
      IF (igtyp /= 22) THEN
       isorth = 0
      END IF
      imas_ds = defaults_solid%IMAS
C
      DO i=1,nel
        rhocp(i) =  pm(69,ixs(1,nft+i))
        temp0(i) =  pm(79,ixs(1,nft+i))
      ENDDO
C
      CALL s6ccoor3(x     ,ixs(1,nf1)   ,geo  ,ngl  ,mat  ,pid  ,
     .     rx   ,ry   ,rz   ,sx   ,sy   ,sz   ,tx   ,ty   ,tz   ,
     .     e1x  ,e1y  ,e1z  ,e2x  ,e2y  ,e2z  ,e3x  ,e3y  ,e3z  ,
     .     f1x  ,f1y  ,f1z  ,f2x  ,f2y  ,f2z  ,temp0, temp,glob_therm%NINTEMP, 
     .     ix1, ix2, ix3, ix4, ix5, ix6, 
     .     x1, x2, x3, x4, x5, x6, 
     .     y1, y2, y3, y4, y5, y6,
     .     z1, z2, z3, z4, z5, z6)
      IF (igtyp == 21 .OR. igtyp == 22) THEN
        DO i=1,nel
         angle(i) =  geo(1,pid(i))
        END DO
        CALL scmorth3(pid  ,geo  ,igeo ,skew ,irep ,gbuf%GAMA   ,
     .         rx   ,ry   ,rz   ,sx   ,sy   ,sz   ,tx   ,ty   ,tz   ,
     .         e1x  ,e1y  ,e1z  ,e2x  ,e2y  ,e2z  ,e3x  ,e3y  ,e3z  ,
     .         ngl  ,angle,nsigi,sigsp,nsigs,sigi ,ixs  ,1    ,
     .         orthoglob,ptsol,nel)
       IF (igtyp == 22) THEN
        nlymax= 200
        ipang = 200
        ipthk = ipang+nlymax
        ippos = ipthk+nlymax
        ipmat = 100
        ig=pid(1)
        mtn0=mtn
        DO i=1,nel
         mat0(i)=mat(i)
         dtx0(i) = ep20
        ENDDO
       END IF
      END IF
      CALL s6cderi3(nel,gbuf%VOL,geo,vzl,ngl,deltax,volu , 
     .     x1, x2, x3, x4, x5, x6, 
     .     y1, y2, y3, y4, y5, y6,
     .     z1, z2, z3, z4, z5, z6)
      IF (idttsh > 0) THEN
         CALL sdlensh3n(nel,llsh,
     .                  x1, x2, x3, x4, x5, x6,
     .                  y1, y2, y3, y4, y5, y6,
     .                  z1, z2, z3, z4, z5, z6)
        DO i=1,nel
          IF (gbuf%IDT_TSH(i)>0) 
     .     deltax(i)=max(llsh(i),deltax(i))
        ENDDO
      END IF        
!
!     Initialize element temperature from /initemp
!      
      IF (jthe == 0 .and. glob_therm%NINTEMP > 0) THEN
        DO i=1,nel
          tempel(i) = one_over_8 *(temp(ixs(2,i)) + temp(ixs(3,i))
     .                           + temp(ixs(4,i)) + temp(ixs(5,i))        
     .                           + temp(ixs(6,i)) + temp(ixs(7,i))        
     .                           + temp(ixs(8,i)) + temp(ixs(9,i))) 
        ENDDO
      ELSE
        tempel(1:nel) = temp0(1:nel)
      END IF
!
      ip=0
      CALL matini(pm      ,ixs    ,nixs       ,x      ,
     .            geo     ,ale_connectivity  ,detonators ,iparg  ,
     .            sigi    ,nel    ,skew       ,igeo   ,
     .            ipart   ,iparts ,
     .            mat     ,ipm    ,nsigs  ,numsol     ,ptsol  ,
     .            ip      ,ngl    ,npf    ,tf         ,bufmat , 
     .            gbuf    ,lbuf   ,mbuf   ,elbuf_str  ,iloadp ,
     .            facload, deltax ,tempel )                         
C
      IF (igtyp == 22) CALL sczero3(gbuf%RHO,gbuf%SIG,gbuf%EINT,nel)
C----------------------------------------
C Thermal initialization
      IF(jthe /=0) CALL atheri(mat,pm,gbuf%TEMP)
C-----------------------------
C Loop on integration points
      DO ilay=1,nlay
        lbuf => elbuf_str%BUFLY(ilay)%LBUF(1,1,1)
        mbuf => elbuf_str%BUFLY(ilay)%MAT(1,1,1)
        l_pla = elbuf_str%BUFLY(ilay)%L_PLA
        l_sigb= elbuf_str%BUFLY(ilay)%L_SIGB
c
       IF (igtyp == 22) THEN
          zi = geo(ippos+ilay,ig)
          wi = geo(ipthk+ilay,ig)
          im=igeo(ipmat+ilay,ig)
         mtn=nint(pm(19,im))
         DO i=1,nel
          mat(i)=im
            angle(i) = geo(ipang+ilay,pid(i))
         ENDDO
       ELSE
          zi = a_gauss(ilay,nlay)
          wi = w_gauss(ilay,nlay)
       ENDIF
c
        DO i=1,nel
          lbuf%VOL0DP(i)= half*wi*(gbuf%VOL(i)+vzl(i)*zi)
          lbuf%VOL(i)= lbuf%VOL0DP(i)
        ENDDO
        IF (igtyp == 22)
     .  CALL scmorth3(pid  ,geo  ,igeo ,skew ,irep ,lbuf%GAMA   ,
     .         rx   ,ry   ,rz   ,sx   ,sy   ,sz   ,tx   ,ty   ,tz   ,
     .         e1x  ,e1y  ,e1z  ,e2x  ,e2y  ,e2z  ,e3x  ,e3y  ,e3z  ,
     .         ngl  ,angle,nsigi,sigsp,nsigs,sigi ,ixs,ilay,
     .         orthoglob,ptsol,nel)
!
!     Initialize element temperature from /initemp
!      
        IF (jthe == 0 .and. glob_therm%NINTEMP > 0) THEN
          DO i=1,nel
            tempel(i) = one_over_8 *(temp(ixs(2,i)) + temp(ixs(3,i))
     .                             + temp(ixs(4,i)) + temp(ixs(5,i))      
     .                             + temp(ixs(6,i)) + temp(ixs(7,i))      
     .                             + temp(ixs(8,i)) + temp(ixs(9,i))) 
          ENDDO
        ELSE
          tempel(1:nel) = temp0(1:nel)
        END IF
!
        CALL matini(pm      ,ixs    ,nixs   ,x         ,
     .              geo     ,ale_connectivity  ,detonators,iparg  ,
     .              sigi    ,nel    ,skew      ,igeo   ,
     .              ipart   ,iparts ,
     .              mat     ,ipm    ,nsigs  ,numsol    ,ptsol  ,
     .              ilay    ,ngl    ,npf    ,tf        ,bufmat ,
     .              gbuf    ,lbuf   ,mbuf   ,elbuf_str ,iloadp ,
     .              facload, deltax ,tempel )
        IF (mtn >= 28) THEN
           nuvar = ipm(8,ixs(1,nft+1))
           idef =1
        ELSE
           nuvar = 0
            IF(mtn == 14 .OR. mtn == 12)THEN
               idef =1
            ELSEIF(mtn == 24)THEN
              idef =1
            ELSEIF(istrain == 1)THEN
             IF(mtn == 1)THEN
               idef =1
             ELSEIF(mtn == 2)THEN
               idef =1
             ELSEIF(mtn == 4)THEN
               idef =1
            ELSEIF(mtn == 3.OR.mtn == 6.OR.mtn == 10
     .       .OR.mtn == 21.OR.mtn == 22.OR.mtn == 23.OR.mtn == 49)THEN
               idef =1
             ENDIF
            ENDIF
        ENDIF
        CALL sigin20b(
     .     lbuf%SIG,pm      ,lbuf%VOL ,sigsp    ,
     .     sigi    ,lbuf%EINT,lbuf%RHO,mbuf%VAR ,lbuf%STRA,
     .     ixs     ,nixs     ,nsigi   ,ilay     ,nuvar    , 
     .     nel     ,iuser    ,idef    ,nsigs    ,strsglob ,
     .     straglob,jhbe     ,igtyp   ,x        ,lbuf%GAMA,
     .     mat     ,lbuf%PLA ,l_pla   ,ptsol    ,lbuf%SIGB,
     .     l_sigb  ,ipm      ,bufmat  ,lbuf%VOL0DP)
c
        IF(igtyp == 22) THEN
C         moyene density,sig,...---
          aire(:) = zero
          CALL dtmain(geo       ,pm        ,ipm         ,pid     ,mat     ,fv    ,
     .         lbuf%EINT ,lbuf%TEMP ,lbuf%DELTAX ,lbuf%RK ,lbuf%RE ,bufmat, deltax, aire, 
     .         volu, dtx , igeo,igtyp)
C
          CALL svalue0(
     .         lbuf%RHO,lbuf%VOL,lbuf%OFF,lbuf%SIG,lbuf%EINT,dtx,
     .         gbuf%RHO,gbuf%VOL,gbuf%OFF,gbuf%SIG,gbuf%EINT,dtx0,
     .         nel     )
        ENDIF
      ENDDO  ! ILAY = 1,NLAY
C----------------------------------------
      IF(igtyp == 22) THEN
       mtn=mtn0
       DO i=1,nel
         mat(i)=mat0(i)
       ENDDO
      ENDIF
C----------------------------------------
C Mass initialization
      CALL s6mass3(gbuf%RHO,mas,partsav,x,v,iparts(nf1),mss(1,nf1),
     .     rhocp,mcp ,mcps(1,nf1),mssa(nf1),gbuf%FILL, volu, 
     .     ix1, ix2, ix3, ix4, ix5, ix6,imas_ds)
C----------------------------------------
C Failure model initialization
      CALL failini(elbuf_str,nptr,npts,nptt,nlay,
     .             ipm,sigsp,nsigi,fail_ini ,
     .             sigi,nsigs,ixs,nixs,ptsol,
     .             rnoise,perturb,mat_param)
C------------------------------------------ 
C  Assemble nodal volumes and moduli for interface stiffness
C  Warning : IX1, IX2 ... IX6 <=> NC(MVSIZ,6)
      IF(i7stifs/=0)THEN
        ncc=6
        CALL sbulk3(volu  ,ix1    ,ncc,mat,pm ,
     2              volnod,bvolnod,vns(1,nf1),bns(1,nf1),bid,
     3              bid ,gbuf%FILL)
      ENDIF
C------------------------------------------
C Element time step
        aire(:) = zero
        CALL dtmain(geo       ,pm        ,ipm         ,pid     ,mat     ,fv    ,
     .       lbuf%EINT ,lbuf%TEMP ,lbuf%DELTAX ,lbuf%RK ,lbuf%RE ,bufmat, deltax, aire, 
     .       volu, dtx, igeo,igtyp)
C------------------------------------------
       IF(igtyp == 22) THEN
        DO i=1,nel
         dtx(i)=dtx0(i)
        ENDDO
       ENDIF
c
      DO i=1,nel
        IF(ixs(10,i+nft) /= 0) THEN
          IF (igtyp < 20 .OR. igtyp > 22) THEN
             ipid1=ixs(nixs-1,i+nft)
             CALL fretitl2(titr1,igeo(npropgi-ltitr+1,ipid1),ltitr)
             CALL ancmsg(msgid=226,
     .                   msgtype=msgerror,
     .                   anmode=aninfo_blind_1,
     .                   i1=igeo(1,ipid1),
     .                   c1=titr1,
     .                   i2=igtyp)
          ENDIF
        ENDIF
        dtelem(nft+i)=dtx(i)
        sti = fourth * gbuf%FILL(i) * gbuf%RHO(i) * volu(i) /
     .        max(em20,dtx(i)*dtx(i))
        stifn(ixs(2,i+nft))=stifn(ixs(2,i+nft))+sti
        stifn(ixs(3,i+nft))=stifn(ixs(3,i+nft))+sti
        stifn(ixs(4,i+nft))=stifn(ixs(4,i+nft))+sti
        stifn(ixs(5,i+nft))=stifn(ixs(5,i+nft))+sti
        stifn(ixs(6,i+nft))=stifn(ixs(6,i+nft))+sti
        stifn(ixs(7,i+nft))=stifn(ixs(7,i+nft))+sti
        stifn(ixs(8,i+nft))=stifn(ixs(8,i+nft))+sti
        stifn(ixs(9,i+nft))=stifn(ixs(9,i+nft))+sti
      ENDDO
C-----------
      RETURN
      END SUBROUTINE s6cinit3
!||====================================================================
!||    sdlensh3n   ../starter/source/elements/thickshell/solide6c/s6cinit3.F
!||--- called by ------------------------------------------------------
!||    s6cinit3    ../starter/source/elements/thickshell/solide6c/s6cinit3.F
!||====================================================================
      SUBROUTINE sdlensh3n(NEL,LLSH3N,
     .                     X1, X2, X3, X4, X5, X6,
     .                     Y1, Y2, Y3, Y4, Y5, Y6,
     .                     Z1, Z2, Z3, Z4, Z5, Z6)
C-----------------------------------------------
C   I m p l i c i t   T y p e s
C-----------------------------------------------
#include      "implicit_f.inc"
C-----------------------------------------------
C   G l o b a l   P a r a m e t e r s
C-----------------------------------------------
#include      "mvsiz_p.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER :: NEL
      my_real,DIMENSION(MVSIZ),INTENT(OUT) :: LLSH3N
      my_real,DIMENSION(MVSIZ),INTENT(IN) :: 
     .   X1, X2, X3, X4, X5, X6,
     .   Y1, Y2, Y3, Y4, Y5, Y6,  
     .   Z1, Z2, Z3, Z4, Z5, Z6
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I
      my_real
     .        E1X(MVSIZ), E1Y(MVSIZ), E1Z(MVSIZ),
     .        E2X(MVSIZ), E2Y(MVSIZ), E2Z(MVSIZ),
     .        E3X(MVSIZ), E3Y(MVSIZ), E3Z(MVSIZ),
     .        x31(mvsiz), y31(mvsiz), z31(mvsiz),
     .        x32(mvsiz), y32(mvsiz), z32(mvsiz),
     .        x21(mvsiz), y21(mvsiz), z21(mvsiz), area(mvsiz),
     .        x2l(mvsiz), x3l(mvsiz), y3l(mvsiz),
     .        xn(mvsiz,3) , yn(mvsiz,3) , zn(mvsiz,3) 
      my_real
     .   al1,al2,al3,almax,sum
C=======================================================================
       DO i=1,nel
          xn(i,1) = half*(x1(i)+x4(i))
          yn(i,1) = half*(y1(i)+y4(i))
          zn(i,1) = half*(z1(i)+z4(i))
          xn(i,2) = half*(x2(i)+x5(i))
          yn(i,2) = half*(y2(i)+y5(i))
          zn(i,2) = half*(z2(i)+z5(i))
          xn(i,3) = half*(x3(i)+x6(i))
          yn(i,3) = half*(y3(i)+y6(i))
          zn(i,3) = half*(z3(i)+z6(i))
       ENDDO 
      DO i=1,nel
        x21(i)=xn(i,2)-xn(i,1)
        y21(i)=yn(i,2)-yn(i,1)
        z21(i)=zn(i,2)-zn(i,1)
        x31(i)=xn(i,3)-xn(i,1)
        y31(i)=yn(i,3)-yn(i,1)
        z31(i)=zn(i,3)-zn(i,1)
        x32(i)=xn(i,3)-xn(i,2)
        y32(i)=yn(i,3)-yn(i,2)
        z32(i)=zn(i,3)-zn(i,2)
      ENDDO
C
      DO i=1,nel
        e1x(i)= x21(i)
        e1y(i)= y21(i)
        e1z(i)= z21(i)
        x2l(i) = sqrt(e1x(i)*e1x(i)+e1y(i)*e1y(i)+e1z(i)*e1z(i))
        e1x(i)=e1x(i)/x2l(i)
        e1y(i)=e1y(i)/x2l(i)
        e1z(i)=e1z(i)/x2l(i)
      ENDDO
C
      DO i=1,nel
        e3x(i)=y31(i)*z32(i)-z31(i)*y32(i)
        e3y(i)=z31(i)*x32(i)-x31(i)*z32(i)
        e3z(i)=x31(i)*y32(i)-y31(i)*x32(i)
        sum = sqrt(e3x(i)*e3x(i)+e3y(i)*e3y(i)+e3z(i)*e3z(i))
        e3x(i)=e3x(i)/sum
        e3y(i)=e3y(i)/sum
        e3z(i)=e3z(i)/sum
        area(i) = half * sum
      ENDDO
C
      DO i=1,nel
        e2x(i)=e3y(i)*e1z(i)-e3z(i)*e1y(i)
        e2y(i)=e3z(i)*e1x(i)-e3x(i)*e1z(i)
        e2z(i)=e3x(i)*e1y(i)-e3y(i)*e1x(i)
        sum = sqrt(e2x(i)*e2x(i)+e2y(i)*e2y(i)+e2z(i)*e2z(i))
        e2x(i)=e2x(i)/sum
        e2y(i)=e2y(i)/sum
        e2z(i)=e2z(i)/sum
        y3l(i)=e2x(i)*x31(i)+e2y(i)*y31(i)+e2z(i)*z31(i)
        x3l(i)=e1x(i)*x31(i)+e1y(i)*y31(i)+e1z(i)*z31(i)
      ENDDO
C
      DO i=1,nel
        al1 = x2l(i) * x2l(i)      
        al2 = (x3l(i)-x2l(i)) * (x3l(i)-x2l(i)) + y3l(i) * y3l(i)
        al3 = x3l(i) * x3l(i) + y3l(i) * y3l(i)
        almax = max(al1,al2,al3)
        llsh3n(i)= two*area(i) / sqrt(almax)
      ENDDO
C
      RETURN
      END SUBROUTINE sdlensh3n