ini_seatbelt.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/.
!||====================================================================
!||    ini_seatbelt   ../starter/source/tools/seatbelts/ini_seatbelt.F
!||--- called by ------------------------------------------------------
!||    initia         ../starter/source/elements/initia/initia.F
!||--- calls      -----------------------------------------------------
!||    ancmsg         ../starter/source/output/message/message.F
!||    ifrontplus     ../starter/source/spmd/node/frontplus.F
!||    new_seatbelt   ../starter/source/tools/seatbelts/new_seatbelt.F
!||    nlocal         ../starter/source/spmd/node/ddtools.F
!||--- uses       -----------------------------------------------------
!||    message_mod    ../starter/share/message_module/message_mod.F
!||====================================================================
      SUBROUTINE ini_seatbelt(IPARG,ELBUF_TAB,KNOD2EL1D,NOD2EL1D,IXR,
     .                        X,ITAB,IPM,ALEA,KNOD2ELC,
     .                        NOD2ELC,IXC)
C-----------------------------------------------
C   M o d u l e s
C-----------------------------------------------
      USE my_alloc_mod
      USE message_mod
      USE seatbelt_mod
      USE elbufdef_mod  
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      "random_c.inc"
#include      "tabsiz_c.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER IPARG(NPARG,NGROUP),KNOD2EL1D(*),NOD2EL1D(*),IXR(NIXR,*),ITAB(*),IPM(NPROPMI,*),
     .        KNOD2ELC(*),NOD2ELC(*),IXC(NIXC,*)
      my_real               :: x(3,*),alea(*)
      TYPE (ELBUF_STRUCT_), TARGET, DIMENSION(NGROUP) :: ELBUF_TAB
C-----------------------------------------------
C   F u n c t i o n
C-----------------------------------------------
      INTEGER  NLOCAL
      EXTERNAL NLOCAL 
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I,J,K,L,M,N,NG,NEL,ITY,IAD,LFT,LLT,NFT,MTN,NODE,NODE_NEXT,
     .        ELEM_NEXT,N1,N2,NTOOL,KK,NOT_USED,FLAG,ELEM_CUR,NN,ID,NNOD,
     .        mid,mtyp,p,nb_shell,nfram,n3,n4,ms,nodes,iseatbelt,nsl,idrb,
     .        nod,nfound_rby,nfound_bcs,bcs_x,bcs_y,bcs_z,ic,ic1,ic2,n1sp,
     .        n2sp,elem,orient
C
      my_real dist1,dist2,alea_max,tole_2
C
      INTEGER , DIMENSION(:), ALLOCATABLE:: TAG_RETRACTOR,TAGN_RETRACTOR,TAG_RES
C
      TYPE(g_bufel_),POINTER :: GBUF
C-----------------------------------------------
C   S o u r c e   L i n e s
C-----------------------------------------------
C
C-----Reset of TAG_RES from SEATBELT_TAB structure----------------------------------
      CALL my_alloc(tag_res,numelr)
      tag_res(1:numelr) = 0
      DO i=1,n_seatbelt
        DO j=1,seatbelt_tab(i)%NSPRING
          tag_res(seatbelt_tab(i)%SPRING(j)) = i
        ENDDO 
      ENDDO
C
C-----Preparation for SPMD computation  ---------------------------------------------
      IF (nspmd > 1) THEN
C
        DO i=1,nslipring
          DO j=1,slipring(i)%NFRAM
            n2 = slipring(i)%FRAM(j)%NODE(2)
            DO p=1,nspmd
              IF ((nlocal(n2,p)==1).AND.(nlocal(slipring(i)%FRAM(j)%ANCHOR_NODE,p)==0)) THEN
C--             Anchor node and orientation node must be stick on the proc of the slipring
                CALL ifrontplus(slipring(i)%FRAM(j)%ANCHOR_NODE,p)
                IF (slipring(i)%FRAM(j)%ORIENTATION_NODE > 0)
     .             CALL ifrontplus(slipring(i)%FRAM(j)%ORIENTATION_NODE,p)
              ENDIF
            ENDDO
          ENDDO
        ENDDO
C
        DO i=1,nretractor
          n2 = retractor(i)%NODE(2)
          DO p=1,nspmd
            IF ((nlocal(n2,p)==1).AND.(nlocal(retractor(i)%ANCHOR_NODE,p)==0)) THEN
C--           Anchor node must be stick on the proc of the retractor
              CALL ifrontplus(retractor(i)%ANCHOR_NODE,p)
            ENDIF
          ENDDO
        ENDDO
C
      ENDIF
 
C-----------------------------------------------------------------------------------
 
      CALL my_alloc(tag_retractor,numelr)
      CALL my_alloc(tagn_retractor,numnod)
      tag_retractor(1:numelr) = 0
      tagn_retractor(1:numnod) = 0
      flag = 1
      not_used = 0
C
C----- Loop on retractors to identify elements in retractor to deactivate (same algorithm as in creat_seatbelt) --
C
      DO i=1,nretractor
C
        n1 = retractor(i)%NODE(1)
        n2 = retractor(i)%NODE(2)
C
C--     loop of elements initially inside the retractor - tag with negative value
C
        nnod = 2
C
        DO k=knod2el1d(n1)+1,knod2el1d(n1+1)
          IF (nod2el1d(k) > numelt+numelp) THEN
            elem_cur = nod2el1d(k)-numelt-numelp
            mid = ixr(5,elem_cur)
            IF (mid > 0) THEN
              mtyp = ipm(2,mid)
              IF (((ixr(2,elem_cur)==n2).OR.(ixr(3,elem_cur)==n2)).AND.(mtyp == 114)) THEN
C--             Loop on belt inside the retractor
                CALL new_seatbelt(ixr,itab,knod2el1d,nod2el1d,n1,
     .                            elem_cur,tag_retractor,tagn_retractor,i,flag,
     .                            not_used,ipm,not_used,not_used,not_used,
     .                            not_used)
              ENDIF
            ENDIF
          ENDIF
        ENDDO
C
      ENDDO
C
C--   Loop to initialise elements with seatbelt material 
C
      DO ng=1,ngroup
C
        mtn = iparg(1,ng)
        nel = iparg(2,ng)
        ity = iparg(5,ng)
        nft = iparg(3,ng)
        iad = iparg(4,ng)
        iseatbelt = iparg(91,ng)
        lft = 1
        llt = nel
C
        IF ((ity==6).AND.(mtn==114)) THEN
C
C---      1D seatbelts springs 
C
          gbuf => elbuf_tab(ng)%GBUF
C
          DO i=1,nel
C
            j = nft + i
            n1 = ixr(2,j)
            n2 = ixr(3,j)          
            ntool = 0
C
C----------------------------------------------------------------------------
C-          Detection if element is in slipring
C----------------------------------------------------------------------------
C
            DO k=1,nslipring
              DO l=1,slipring(k)%NFRAM
                IF (((n1==slipring(k)%FRAM(l)%NODE(1)).AND.(n2==slipring(k)%FRAM(l)%NODE(2))).OR.
     .              ((n1==slipring(k)%FRAM(l)%NODE(2)).AND.(n2==slipring(k)%FRAM(l)%NODE(1)))) THEN
CC-               element is strand number 1
                  ntool = ntool + 1
                  gbuf%SLIPRING_ID(i) = k
                  gbuf%SLIPRING_FRAM_ID(i) = l
                  gbuf%SLIPRING_STRAND(i) = 1
C
CC-               determinaton of slipring direction for strand1
                  IF (n1==slipring(k)%FRAM(l)%NODE(1)) THEN
                    slipring(k)%FRAM(l)%STRAND_DIRECTION(1) = 1
                  ELSE
                    slipring(k)%FRAM(l)%STRAND_DIRECTION(1) = -1
                  ENDIF
C
                ELSEIF (((n1==slipring(k)%FRAM(l)%NODE(2)).AND.(n2==slipring(k)%FRAM(l)%NODE(3))).OR.
     .                  ((n1==slipring(k)%FRAM(l)%NODE(3)).AND.(n2==slipring(k)%FRAM(l)%NODE(2)))) THEN
CC-               element is strand number 2
                  ntool = ntool + 1
                  gbuf%SLIPRING_ID(i) = k
                  gbuf%SLIPRING_FRAM_ID(i) = l
                  gbuf%SLIPRING_STRAND(i) = 2
C
CC-               determinaton of slipring direction for strand2
                  IF (n1==slipring(k)%FRAM(l)%NODE(2)) THEN
                    slipring(k)%FRAM(l)%STRAND_DIRECTION(2) = 1
                  ELSE
                    slipring(k)%FRAM(l)%STRAND_DIRECTION(2) = -1
                  ENDIF
C
                ENDIF
              ENDDO
            ENDDO
C
C----------------------------------------------------------------------------
C-          Detection if element is in retractor
C----------------------------------------------------------------------------
C
            DO k=1,nretractor
              IF (((n1==retractor(k)%NODE(1)).AND.(n2==retractor(k)%NODE(2))).OR.
     .            ((n2==retractor(k)%NODE(1)).AND.(n1==retractor(k)%NODE(2)))) THEN
CC-             element is mouth element of retractor - 1st direction
                ntool = ntool + 1
                gbuf%RETRACTOR_ID(i) = k
                gbuf%SLIPRING_STRAND(i) = -1
              ENDIF
            ENDDO
C
            IF (tag_retractor(j) > 0) THEN
C-            element initially in the retractor
              gbuf%OFF(i) = zero
              gbuf%RETRACTOR_ID(i) = -k
              k = tag_retractor(j)
              id = retractor(k)%ID
              nn = retractor(k)%ANCHOR_NODE 
              dist1 = (x(1,nn)-x(1,n1))**2+(x(2,nn)-x(2,n1))**2+(x(3,nn)-x(3,n1))**2
              dist2 = (x(1,nn)-x(1,n2))**2+(x(2,nn)-x(2,n2))**2+(x(3,nn)-x(3,n2))**2
C
C--           default tolerance
              tole_2 = em10*retractor(k)%ELEMENT_SIZE*retractor(k)%ELEMENT_SIZE
C--           compatibility with random noise
              IF (nrand > 0) THEN
                alea_max = zero
                DO j=1,nrand
                  alea_max = max(alea_max,alea(j))
                ENDDO 
                tole_2 = max(tole_2,ten*alea_max*alea_max)
              ENDIF
 
C--           tolerance if nodes are very close to anchorage node
              IF (dist1 <= tole_2) THEN
                  x(1,n1) = x(1,nn)
                  x(2,n1) = x(2,nn)
                  x(3,n1) = x(3,nn)
                  dist1 = zero
              ENDIF
              IF (dist2 <= tole_2) THEN
                  x(1,n2) = x(1,nn)
                  x(2,n2) = x(2,nn)
                  x(3,n2) = x(3,nn)
                  dist2 = zero
              ENDIF
C
              IF(dist2 + dist1 > em30) CALL ancmsg(msgid=2011,
     .                                 msgtype=msgerror,
     .                                 anmode=aninfo,
     .                                 i1=id,i2=ixr(nixr,j),i3=id)
C
              IF (tagn_retractor(n1) > 0) THEN
                retractor(k)%INACTI_NNOD = retractor(k)%INACTI_NNOD + 1
                retractor(k)%INACTI_NODE(retractor(k)%INACTI_NNOD) = n1
                tagn_retractor(n1) = 0
              ENDIF
C
              IF (tagn_retractor(n2) > 0) THEN
                retractor(k)%INACTI_NNOD = retractor(k)%INACTI_NNOD + 1
                retractor(k)%INACTI_NODE(retractor(k)%INACTI_NNOD) = n2
                tagn_retractor(n2) = 0
              ENDIF
C               
            ENDIF
C
C----------------------------------------------------------------------------
C
            IF(ntool > 1) CALL ancmsg(msgid=2006,
     .                                msgtype=msgerror,
     .                                anmode=aninfo,
     .                                i1=ixr(nixr,j))
C
CCC-        Find node before node 1
            node = n1
            node_next = 0
            DO k=knod2el1d(node)+1,knod2el1d(node+1)
              IF ((nod2el1d(k) /= j + numelt+numelp).AND.(nod2el1d(k) > numelt+numelp)) THEN
                elem_next = nod2el1d(k)-numelt-numelp
                mid = ixr(5,elem_next)
                IF (mid > 0) THEN
                  mtyp = ipm(2,mid)
                  IF (mtyp == 114) THEN
                    IF (ixr(2,elem_next) == node) THEN
                      node_next = ixr(3,elem_next)
                    ELSE
                      node_next = ixr(2,elem_next)
                    ENDIF
                  ENDIF
                ENDIF                
              ENDIF
            ENDDO
            gbuf%ADD_NODE(i) = node_next
 
CCC-        Find node after node 2
            node = n2
            node_next = 0
            DO k=knod2el1d(node)+1,knod2el1d(node+1)
              IF ((nod2el1d(k) /= j + numelt+numelp).AND.(nod2el1d(k) > numelt+numelp)) THEN
                elem_next = nod2el1d(k)-numelt-numelp
                mid = ixr(5,elem_next)
                IF (mid > 0) THEN
                  mtyp = ipm(2,mid)
                  IF (mtyp == 114) THEN
                    IF (ixr(2,elem_next) == node) THEN
                      node_next = ixr(3,elem_next)
                    ELSE
                      node_next = ixr(2,elem_next)
                    ENDIF
                  ENDIF
                ENDIF                
              ENDIF
            ENDDO
            gbuf%ADD_NODE(nel+i) = node_next
C
CCC-        Compute Fram factor for 2D seatbelt
C
            nb_shell = 0
            IF (gbuf%ADD_NODE(i) > 0) THEN
              node = n1
            ELSE
              node = n2
            ENDIF
            DO k=knod2elc(node)+1,knod2elc(node+1)
              elem_cur = nod2elc(k)
              mid = ixc(1,elem_cur)
              mtyp = ipm(2,mid)
              IF (mtyp == 119) nb_shell = nb_shell + 1                
            ENDDO
C
            nfram = seatbelt_tab(tag_res(j))%NFRAM
            IF ((nfram > 1).AND.(nb_shell==4)) THEN
C--           spring is on the edge of the 2D belt
              gbuf%FRAM_FACTOR(i) = one/(nfram-1)
            ELSEIF ((nfram > 1).AND.(nb_shell==2)) THEN
C--           spring is inside the 2D belt
              gbuf%FRAM_FACTOR(i) = half/(nfram-1)
            ELSE
C--           1D seatbelt
              gbuf%FRAM_FACTOR(i) = one
            ENDIF
C--         element mass and inertia scaled by frame factor for elementary time step
            gbuf%MASS(i) = gbuf%MASS(i)*gbuf%FRAM_FACTOR(i)
            gbuf%INTVAR(i) = gbuf%INTVAR(i)*gbuf%FRAM_FACTOR(i)
C
          ENDDO  
C
        ELSEIF ((ity==3).AND.(iseatbelt==1)) THEN
C
C---      2D seatbelts shells
C
          gbuf => elbuf_tab(ng)%GBUF
C
          DO i=1,nel
C
            j = nft + i        
C
            n1 = ixc(2,j)
            n2 = ixc(3,j)
            n3 = ixc(4,j)
            n4 = ixc(5,j)
C
C--         Scale factor on stress initialised 
            gbuf%INTVAR(i) = one
C
C--         Initial length of fram 1/2 and fram 4/3
            dist1 = (x(1,n2)-x(1,n1))**2+(x(2,n2)-x(2,n1))**2+(x(3,n2)-x(3,n1))**2
            dist2 = (x(1,n3)-x(1,n4))**2+(x(2,n3)-x(2,n4))**2+(x(3,n3)-x(3,n4))**2
            gbuf%INTVAR(i+2*nel) = sqrt(dist1)
            gbuf%INTVAR(i+3*nel) = sqrt(dist2)
C
CCC-        Find node defining 1rst orthotropy direction with N1 using springs : N2 (default) or N4 
            gbuf%ADD_NODE(i) = n2 
            DO k=knod2el1d(n1)+1,knod2el1d(n1+1)
              IF (nod2el1d(k) > numelt+numelp) THEN          
                elem = nod2el1d(k)-numelt-numelp
                mid = ixr(5,elem)
                IF (mid > 0) THEN
                  mtyp = ipm(2,mid)
                  IF ((mtyp==114).AND.((ixr(2,elem)==n1).AND.(ixr(3,elem) == n4))
     .                            .OR.((ixr(3,elem)==n1).AND.(ixr(2,elem) == n4))) THEN
                    gbuf%ADD_NODE(i) = n4
                  ENDIF
                ENDIF              
              ENDIF
            ENDDO
C
CCC-        Find node after node i
            DO m=2,5
              node = ixc(m,j)
              node_next = 0
              DO k=knod2el1d(node)+1,knod2el1d(node+1)
                IF (nod2el1d(k) > numelt+numelp) THEN          
                  elem_next = nod2el1d(k)-numelt-numelp
                  mid = ixr(5,elem_next)
                  IF (mid > 0) THEN
                    mtyp = ipm(2,mid)
                    IF (mtyp==114) THEN
                      n1sp=ixr(2,elem_next)
                      n2sp=ixr(3,elem_next)
                      IF ((n1sp==node).AND.(n2sp/=n1).AND.(n2sp/=n2).AND.(n2sp/=n3).AND.(n2sp/=n4)) THEN    
                        node_next = n2sp
                      ELSEIF ((n2sp==node).AND.(n1sp/=n1).AND.(n1sp/=n2).AND.(n1sp/=n3).AND.(n1sp/=n4)) THEN 
                        node_next = n1sp
                      ENDIF
                    ENDIF  
                  ENDIF              
                ENDIF
              ENDDO
              gbuf%ADD_NODE((m-1)*nel+i) = node_next
C
            ENDDO
C
          ENDDO
C
        ENDIF
C
      ENDDO         
C
      DEALLOCATE(tag_res,tag_retractor,tagn_retractor)
C
      END SUBROUTINE ini_seatbelt