Sources Esope | Cast3M

gurson
C GURSON    SOURCE    CHAT      05/01/13    00:22:25     5004
c*****************************************************************************
c       sous rouine utilisee dans l'ecoulement plastique de la loi
c       de gurson
c**********************************************************************
        subroutine gurson (s,p,depsp,gamma,poro,densit,epsi,dtime,
     &  wrkgur)
c
c       this subroutine is the interface between the finite element program
c       and the rest of the material update program
c       given the history variables and the starin rate it updates the
c       state
c               by arnaud de gayffier 06/13/93 in san diego california
c
c---    arguments
c               s........deviatoric stress  in gpa
c                               s(1)=sxx|
c                               s(2)=syy| s(1)+s(2)+s(3)=0
c                               s(3)=szz|
c                               s(4)=sxy
c                               s(5)=sxz
c                               s(6)=syz
c               p........pressure (if >0 compression)  in gpa
c               gamma....cumulative plastic strain
c               poro.....porosity or void fraction
c               densit...material density  in gr/cm3
c               epsi.....total strain rate same conventions as "s"
c               depsp....total plasti strain rate same convention
c               dtime....time increment
c
c         warning!! the argument epsi might be modified
c
c---    variables
c
      IMPLICIT INTEGER(I-N)
        IMPLICIT REAL*8 (A-H,O-Z)
        segment wrkgur
         real*8 sigbar, sy0,phi0,rho0,g,b,h
         real*8 epn,phin,sqrtj2,rho,sig(6)
         real*8 e(7),dt
         real*8 conv,tol1,tol2
        endsegment
 
*        COMMON/PROP/SIGBAR,SY0,PHI0,RHO0,G,B,H
*        COMMON/STATE/EPN,PHIN,SQRTJ2,RHO,SIG(6)
c               e is the strain rate:
c                       e(1-6) is the deviatoric part
c                       e(7) is the trace/3
c               dt ist the time increment
*        COMMON /DELTA/ E(7),DT
        dimension s(6),epsi(6),depsp(6)
c
c---    identification
        DPHI=0.D0
        SQRTJ2=0.D0
        epn=gamma
        phin=poro
        rho=densit
        do 100 i=1,6
          e(i)=epsi(i)
          SIG(I)=S(I)
  100   CONTINUE
        dt=dtime
c
c---    partition in deviatoric and hydrostatic
        E(7)=(E(1)+E(2)+E(3))/3.D0
        E(1)=E(1)-E(7)
        E(2)=E(2)-E(7)
        E(3)=E(3)-E(7)
        E(4)=E(4)/2.D0
        E(5)=E(5)/2.D0
        E(6)=E(6)/2.D0
c
  101   FORMAT (A,7(1PE12.5,1X))
  102   FORMAT (5(A,1PE12.5))
c
c---    material update
        call newsta(wrkgur)
c
c---    arguments update
        DPHI = PHIN-PORO
        DEPN = EPN-GAMMA
c
        P=(1.D0-PHIN)*B*(1.D0-(1.D0-PHIN)*RHO0/((1.D0-PHI0)*RHO))
        gamma=epn
        poro=phin
        densit=rho
        do 200 i=1,6
  200   S(I)=SIG(I)
c                plastic strain
        DO 300 I=1,6
        IF (SQRTJ2 .GE. 1.D-6 ) THEN
          DEPSP(I)=1.5D0*SIG(I)/SQRTJ2*DEPN
        ELSE
          DEPSP(I) =0.D0
        ENDIF
  300   CONTINUE
c
        DEPSP(4)=2.D0*DEPSP(4)
        DEPSP(5)=2.D0*DEPSP(5)
        DEPSP(6)=2.D0*DEPSP(6)
c
        DO 400 I=1,3
         DEPSP(I)=DEPSP(I)+DPHI/(1.D0-PHIN)/3.D0
  400   CONTINUE
c
        RETURN
        end
c------------------------------------------------------------------------------