sigeps34t.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/.
!||====================================================================
!||    sigeps34t   ../engine/source/materials/mat/mat034/sigeps34t.F
!||--- called by ------------------------------------------------------
!||    tforc3      ../engine/source/elements/truss/tforc3.F
!||====================================================================
      SUBROUTINE sigeps34t(NEL     ,NGL     ,MAT     ,PID     ,UPARAM  ,
     .                     IPM     ,GEO     ,OFF     ,FOR     ,STI     ,
     .                     EINT    ,AREA    ,AL0     ,AL      ,EPSP    ,
     .                     NUVAR   ,UVAR    )
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      "com04_c.inc"
#include      "com08_c.inc"
#include      "comlock.inc"
C-----------------------------------------------
C   D u m m y   A r g u m e n t s
C-----------------------------------------------
      INTEGER ,INTENT(IN)  :: NEL,NUVAR
      INTEGER ,DIMENSION(NEL) ,INTENT(IN)  :: MAT,PID,NGL
      INTEGER ,DIMENSION(NPROPMI,NUMMAT) ,INTENT(IN) :: IPM
      my_real ,DIMENSION(NPROPG ,NUMGEO) ,INTENT(IN) :: geo
      my_real ,DIMENSION(*) ,INTENT(IN) :: uparam
      my_real ,DIMENSION(NEL) :: off,for,eint,area,al0,al,sti,epsp
      my_real ,DIMENSION(NEL,NUVAR) :: uvar
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER :: I,J,IADBUF,NINDX
      INTEGER ,DIMENSION(NEL) :: INDX,ICC,ISRATE,VFLAG
      my_real :: K3,NU2,DSIG,H01,H02,H1,H2,SV,DDEXX,DEPSDXX,
     .  dexx,depsvxx,dp
      my_real, DIMENSION(NEL) :: bulk,g_ins,g_inf,ge,ge2,gv,gv2,beta,
     .   deps,gap,eps,c1,c2
c=======================================================================     
      DO i=1,nel
        iadbuf   = ipm(7,mat(i))-1
        bulk(i)  = uparam(iadbuf+1)
        g_ins(i) = uparam(iadbuf+2)
        g_inf(i) = uparam(iadbuf+3)
        beta(i)  = uparam(iadbuf+4)
        gap(i)   = geo(2,pid(i))
      ENDDO
c        
      ge(1:nel)  = g_inf(1:nel)                ! elastic part of shear modulus
      gv(1:nel)  = g_ins(1:nel)-g_inf(1:nel)   ! viscous part of shear modulus
      ge2(1:nel) = ge(1:nel) * two
      gv2(1:nel) = gv(1:nel) * two 
      c1(1:nel)  = one - exp(-beta(1:nel)*dt1)
      c2(1:nel)  =-c1(1:nel) / beta(1:nel)
c      
      deps(1:nel) = epsp(1:nel) * dt1   ! strain increment
      eps(1:nel) = uvar(1:nel,2) + deps(1:nel) ! total normal strain
      uvar(1:nel,2) = eps(1:nel) 
c
      DO i=1,nel
        IF (gap(i) > zero .AND. al(i) <= (al0(i)-gap(i))) off(i)=one
      ENDDO
c
      DO i=1,nel
        k3  = bulk(i)*three
        nu2 = (k3 - ge2(i)) / (k3 + ge(i))
        nu2 = max(nu2, one)
        eint(i) = eint(i) + for(i)*deps(i)*al(i)*half
        area(i) = area(i)*(one - nu2*deps(i)*off(i))
        sti(i)  = k3
      ENDDO
c
      DO i=1,nel
        ddexx = deps(i)*two_third   ! strain increment deviator    
        depsdxx = two_third*epsp(i)
        dexx = two_third*eps(i)
        !
        depsvxx = c1(i)*(dexx - uvar(i,1))     + c2(i)*depsdxx
        dp  = bulk(i)*deps(i)
        !!        
        dsig = ge2(i)*ddexx - gv2(i)*depsvxx + dp
        for(i) = for(i) + dsig * area(i)
        sti(i) = max(sti(i), abs(dsig / max(ddexx,em20)))
        sti(i) = sti(i)*off(i)
        for(i) = for(i)*off(i)  
        uvar(i,1) =  uvar(i,1) + depsvxx + ddexx
      END DO
c
      DO i=1,nel
        eint(i) = eint(i) + for(i)*deps(i)*al(i)*half
      ENDDO
c-----------
      RETURN
      END