#include "implicit_f.inc"
#include "mvsiz_p.inc"
#include "param_c.inc"

Functions/Subroutines

subroutine sigeps91gpinch (jft, jlt, nuvar, uparam, rho0, thk, thk0, nel, ssp, area, depsxx, depsyy, depszz, depsxy, depsyz, depszx, depbxx, depbyy, depbxy, deppxz, deppyz, sigoxx, sigoyy, sigozz, sigoxy, sigoyz, sigozx, momoxx, momoyy, momoxy, momopxz, momopyz, signxx, signyy, signzz, signxy, signyz, signzx, momnxx, momnyy, momnxy, momnpxz, momnpyz, time, uvar, dt_inv, viscmx, aldt, vol0, ipm, mat, pla, degmb, degfx, ngl, ezzavg, areapinch)

Function/Subroutine Documentation

◆ sigeps91gpinch()

subroutine sigeps91gpinch	(	integer	jft,
		integer	jlt,
		integer	nuvar,
			uparam,
			rho0,
			thk,
			thk0,
		integer	nel,
			ssp,
			area,
			depsxx,
			depsyy,
			depszz,
			depsxy,
			depsyz,
			depszx,
			depbxx,
			depbyy,
			depbxy,
			deppxz,
			deppyz,
			sigoxx,
			sigoyy,
			sigozz,
			sigoxy,
			sigoyz,
			sigozx,
			momoxx,
			momoyy,
			momoxy,
			momopxz,
			momopyz,
			signxx,
			signyy,
			signzz,
			signxy,
			signyz,
			signzx,
			momnxx,
			momnyy,
			momnxy,
			momnpxz,
			momnpyz,
			time,
			uvar,
			dt_inv,
			viscmx,
			aldt,
			vol0,
		integer, dimension(npropmi,*)	ipm,
		integer, dimension(nel)	mat,
			pla,
			degmb,
			degfx,
		integer, dimension(*)	ngl,
			ezzavg,
			areapinch )

Definition at line 28 of file sigeps91gpinch.F.

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      "param_c.inc"
C-----------------------------------------------
C   I N P U T   A R G U M E N T S
C-----------------------------------------------
      INTEGER JFT,JLT,NEL,ISMSTR,NUVAR, IPM(NPROPMI,*),MAT(NEL),NGL(*)
C     REAL
      my_real
     .  uparam(*),uvar(nel,nuvar)
      my_real
     .  time
      my_real
     .   depsxx(nel),depsyy(nel),depszz(nel),
     .   depsxy(nel),depsyz(nel),depszx(nel),
     .   depbxx(nel),depbyy(nel),depbxy(nel),
     .   deppxz(nel),deppyz(nel),
     .   sigoxx(nel),sigoyy(nel),sigozz(nel),
     .   sigoxy(nel),sigoyz(nel),sigozx(nel),
     .   momoxx(nel),momoyy(nel),momoxy(nel),
     .   momopxz(nel),momopyz(nel),
     .   area(nel),thk(nel),thk0(nel),rho0(nel),dt_inv(nel),aldt(mvsiz),
     .   vol0(mvsiz),pla(nel),degmb(mvsiz),degfx(mvsiz),ezzavg(mvsiz),
     .   areapinch(mvsiz)
C-----------------------------------------------
C   O U T P U T   A R G U M E N T S
C-----------------------------------------------
      my_real
     .    signxx(nel),signyy(nel),signzz(nel),
     .    signxy(nel),signyz(nel),signzx(nel),
     .    momnxx(nel),momnyy(nel),momnxy(nel),
     .    momnpxz(nel),momnpyz(nel),
     .    ssp(nel),rho(nel)
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER I,MX,J,IADBUF
C     REAL
      my_real
     .   pa1,pa2,pa3,pa4,pa5,
     .   e,nu,g,gs,sigy,hm,sigy0,
     .   ms,fs,d1,d2,d3,c1,c2,c3, 
     .   aaa,bbb,ccc,p,bulk
      my_real
     .   vol(mvsiz),spheps(mvsiz), 
     .   devdepsxx(mvsiz),devdepsyy(mvsiz),devdepszz(mvsiz),
     .   devdepsxy(mvsiz),devdepsyz(mvsiz),devdepsxz(mvsiz),
     .   b1(mvsiz),b2(mvsiz),b3(mvsiz),b4(mvsiz),b5(mvsiz),
     .   thk08(mvsiz),viscmx(mvsiz),thkx(mvsiz),mu(mvsiz),
     .   svonm(mvsiz),crit(mvsiz),sphsig(mvsiz),magdev(mvsiz),invmagdev(mvsiz),
     .   devsxx(mvsiz),devsyy(mvsiz),devszz(mvsiz),
     .   devsxy(mvsiz),devsyz(mvsiz),devszx(mvsiz),
     .   n1(mvsiz),n2(mvsiz),n3(mvsiz),n4(mvsiz),n5(mvsiz),n6(mvsiz),
     .   rr(mvsiz),unsyeq(mvsiz),hh(mvsiz),etse(mvsiz),
     .   degsh_loc(mvsiz),degmb_loc(mvsiz),degfx_loc(mvsiz),dwelm(mvsiz),
     .   dwelf(mvsiz),dwpla(mvsiz),dpla(mvsiz),
     .   sdevxx(mvsiz),sdevyy(mvsiz),sdevzz(mvsiz),
     .   signdevxx(mvsiz),signdevyy(mvsiz),signdevzz(mvsiz),pnew(mvsiz),ptrial(mvsiz),
     .   pold(mvsiz),sigodevxx(mvsiz),sigodevyy(mvsiz),sigodevzz(mvsiz),dd(mvsiz)
C-----------------------------------------------
C
C initialize
      IF (time == zero) THEN
        DO i=1,nel
          uvar(i,1) = areapinch(i)*thk(i) 
          uvar(i,2) = thk(i)
        ENDDO
      ENDIF
C retrieve material parameters from the buffer       
        mx = mat(1)
        iadbuf = ipm(7,mx)
        e = uparam(iadbuf)
        nu  = uparam(iadbuf+1)
        sigy0    = uparam(iadbuf+2)
        hm      = uparam(iadbuf+3)
        g    = half*e/(one+nu)
        gs   = 5.d0/6.d0*g
        bulk = e/(three*(one-two*nu))
C other derived material parameters
        pa1  = e*(one-nu)/(one+nu)/(one-two*nu)
        pa2  = e*nu/(one+nu)/(one-two*nu)
        pa3  = g
        pa4  = (one+nu)*(one-two*nu)/(one-nu**2)/(one-nu)*pa1
        pa5  = (one+nu)*(one-two*nu)/(one-nu**2)*pa2 
C C1/C2/C3 are PM(28:30) [ref:sigeps02g]
        c1 = one/e
        c2 = -nu*c1
        c3 = one/g
C
      DO i=1,nel
        thk08(i)= thk0(i)*one_over_12
        b1(i)   = pa1*thk08(i)
        b2(i)   = pa2*thk08(i)
        b3(i)   = pa3*thk08(i)
        b4(i)   = pa4*thk08(i)
        b5(i)   = pa5*thk08(i)
      ENDDO
C
      DO i=1,nel
C        THKX(I) = UVAR(I,2)*(1+HALF*DEPSZZ(I))/(1-HALF*DEPSZZ(I))
         thkx(i) = uvar(i,2)*(1+half*ezzavg(i))/(1-half*ezzavg(i))
 
        uvar(i,2) = thkx(i)
        vol(i) = areapinch(i)*thkx(i)
        rho(i) = uvar(i,1)*rho0(i)/vol(i)
        mu(i) = rho(i)/rho0(i)-one
        pnew(i) = bulk*mu(i)
C        print*, NGL(I),PNEW(I)
      ENDDO
 
        
      DO i=jft,jlt
C       energy of transverse shear + pinching
        degsh_loc(i) = sigoyz(i)*depsyz(i)+sigozx(i)*depszx(i)
C       energy of only in plane components, xx,yy,xy
        degmb_loc(i) = degmb(i) - degsh_loc(i) 
C       bending energy
        degfx_loc(i) = degfx(i) 
C       calculate old pressure and old derivatoric stress
        pold(i) = -third*(sigoxx(i)+sigoyy(i)+sigozz(i))
        sigodevxx(i) = sigoxx(i)+pold(i)
        sigodevyy(i) = sigoyy(i)+pold(i)
        sigodevzz(i) = sigozz(i)+pold(i)
C       calculate change in spherical strain = DD
        dd(i) = third*(depsxx(i)+depsyy(i)+depszz(i))
C       update trial deviatoric stress to calculate trial deviatoric stress
        signdevxx(i)=sigodevxx(i)+two*g*(depsxx(i)-dd(i))
        signdevyy(i)=sigodevyy(i)+two*g*(depsyy(i)-dd(i))
        signdevzz(i)=sigodevzz(i)+two*g*(depszz(i)-dd(i))
C       shear stresses
        signxy(i)=sigoxy(i)+g*depsxy(i)
        signyz(i)=sigoyz(i)+gs*(depsyz(i)+zero*deppyz(i))
        signzx(i)=sigozx(i)+gs*(depszx(i)+zero*deppxz(i))
C       moments
        momnxx(i)=momoxx(i)+b4(i)*depbxx(i)+b5(i)*depbyy(i)
        momnyy(i)=momoyy(i)+b5(i)*depbxx(i)+b4(i)*depbyy(i)
        momnxy(i)=momoxy(i)+b3(i)*depbxy(i)
        momnpxz(i)=momopxz(i)+b3(i)*deppxz(i)
        momnpyz(i)=momopyz(i)+b3(i)*deppyz(i)
C       Sound speed like solids
        ssp(i) = sqrt(pa1/rho0(i))
 
 
 
 
 
 
      ENDDO
C     calculate equivalent stress 2D
      DO i=jft,jlt
        ms = momnxx(i)+momnyy(i)
C     theory manual (global plasticity algorithm)
        unsyeq(i) = one/
     .    sqrt(max(sixteen*(ms*ms + three*(momnxy(i)*momnxy(i) - momnxx(i)*momnyy(i)))
     .         +three*half*(signdevxx(i)**2+signdevyy(i)**2+signdevzz(i)**2+two*signxy(i)**2),em20))
        sigy = sigy0 + hm*pla(i)
        rr(i)  = min(one,sigy*unsyeq(i))
      ENDDO
C
      DO i=jft,jlt
        IF (rr(i) < one) THEN
          hh(i) = zero
          etse(i) = zero
        ENDIF
      ENDDO
C     energies calculation
      DO i=jft,jlt
C
          signxx(i) = signdevxx(i)*rr(i)-pnew(i)
          signyy(i) = signdevyy(i)*rr(i)-pnew(i)
          signxy(i) = signxy(i)*rr(i)
          signzz(i) = signdevzz(i)*rr(i)-pnew(i)
C
          d1 = signxx(i)-sigoxx(i)
          d2 = signyy(i)-sigoyy(i)
          d3 = signzz(i)-sigozz(i)
C
          dwelm(i) = (signxx(i)+sigoxx(i))*(c1*d1+c2*d2+c2*d3)+
     .               (signyy(i)+sigoyy(i))*(c2*d1+c1*d2+c2*d3)+
     .               (signzz(i)+sigozz(i))*(c2*d1+c2*d2+c1*d3)+
     .               (signxy(i)+sigoxy(i))*(c3*(signxy(i)-sigoxy(i)))
          degmb_loc(i) = degmb_loc(i)+signxx(i)*depsxx(i)+signyy(i)*depsyy(i)
     .                               +signzz(i)*depszz(i)
     .                               +signxy(i)*depsxy(i)
!
          momnxx(i) = momnxx(i)*rr(i)
          momnyy(i) = momnyy(i)*rr(i)
          momnxy(i) = momnxy(i)*rr(i)
          d1 = momnxx(i)-momoxx(i)
          d2 = momnyy(i)-momoyy(i)
          dwelf(i) = twelve*(
     .              (momnxx(i)+momoxx(i))*(c1*d1+c2*d2)
     .             +(momnyy(i)+momoyy(i))*(c2*d1+c1*d2)
     .             +(momnxy(i)+momoxy(i))*(c3*(momnxy(i)-momoxy(i))) )
          degfx_loc(i) = degfx_loc(i)+ momnxx(i)*depbxx(i)+momnyy(i)*depbyy(i)
     .                                +momnxy(i)*depbxy(i)
 
      ENDDO
C     delta of plastic work  
      DO i=jft,jlt
          dwpla(i) = degmb_loc(i)+degfx_loc(i)*thk0(i)-dwelm(i)-dwelf(i)
      ENDDO
C     equivalent plastic strain 
      DO i=jft,jlt
          dpla(i) = max(zero,half*dwpla(i)/sigy)
          pla(i) = pla(i) + dpla(i)
      ENDDO  
      RETURN

OpenRadioss 2025.1.11 OpenRadioss project