C MATER3TG  SOURCE    AM        12/02/28    00:28:36     7289
      SUBROUTINE mater3tg(cmat,ntot,idimpara3,
     .              parahot3,npttot3,lcp,NOEL,NPT)
c END THG add variables


C     This routine is run at the end of every time step
c       when convergence has been obtained.

c       Ntot is equal to ntot1 or ntot2 depending on the finite element type

C
C     *****************************************************************
C     CALCULATION OF THE PLASTIC STRAIN FOR ALL THE MATERIAL LAWS
C     *****************************************************************
C

      IMPLICIT REAL*8 (A-B,D-H,O-Z)
      implicit integer (I-K,M,N)
      implicit logical (L)
      implicit character*10 (C)

      dimension parahot3(idimpara3,npttot3)
      dimension Hinv3(6,6),trav3(6),esigma6v(6)
      dimension epssig6(6),epsipl6(6)

      i1 = 1
      i2 = 2
      i3 = 3
      i4 = 4
      i5 = 5
      i6 = 6
      i7 = 7
      i8 = 8
      i9 = 9
      i10=10
      i11=11
      i12=12
      i13=13
      i14=14
      i15=15
      i16=16
      i17=17
      i18=18
      i19=19
      i20=20
      i21=21
      i22=22
      i23=23
      i24=24
      i25=25
      i26=26
      i27=27
      i28=28
      i29=29
      i30=30
      i31=31
      i32=32
      i33=33
      i34=34
      i35=35
      i36=36
      i37=37
      i38=38
      i39=39

      r0 = 0.
      r1 = 1.


c     3D concrete model (lcp = .false.)
c     and plane stress version (lcp = .true.)

c     Store the mechanical strain at the end of this step
      do iloc=i1,i6
        PARAHOT3(idimpara3-31+iloc,ntot) =
     *          PARAHOT3(idimpara3-19+iloc,ntot)
      end do
c      if (lcp) then
c        PARAHOT2(idimpara2-14,ntot)=parahot3(idimpara3-30,ntot)
c        PARAHOT2(idimpara2-13,ntot)=parahot3(idimpara3-29,ntot)
c        PARAHOT2(idimpara2-12,ntot)=parahot3(idimpara3-27,ntot)
c      endif

c     Store the equivalent plastic strain
      parahot3(idimpara3-36,ntot) = parahot3(idimpara3-34,ntot)
      parahot3(idimpara3-37,ntot) = parahot3(idimpara3-35,ntot)

c     Store the damage in compression
      parahot3(idimpara3-40,ntot) = parahot3(idimpara3-32,ntot)
c     Store the damage in tension
      parahot3(idimpara3-41,ntot) = parahot3(idimpara3-33,ntot)

c     Store the thermal strain at the end of this step
      parahot3(idimpara3-31,ntot) = parahot3(idimpara3,ntot)
c      if (lcp) then
c        parahot2(idimpara2-15,ntot) =  parahot3(idimpara3-31,ntot)
c      endif

c     calcul des deformations plastiques.
      Eo=parahot3(i1,ntot)
      poison=parahot3(i2,ntot)
      rloc=r1/Eo
c     calcul de la matrice D-1 pour relation sigma effectif - eps et calcul
c     de eps elastique en fonction de sigma effectif
      do iloc=i1,i6
        do jloc=i1,i6
          HINV3(iloc,jloc) = r0
        end do
      end do

      Hinv3(i1,i1)=rloc
      Hinv3(i2,i2)=rloc
      Hinv3(i3,i3)=rloc
      Hinv3(i1,i2)=-rloc*poison
      Hinv3(i2,i1)=-rloc*poison
      Hinv3(i1,i3)=-rloc*poison
      Hinv3(i3,i1)=-rloc*poison
      Hinv3(i2,i3)=-rloc*poison
      Hinv3(i3,i2)=-rloc*poison
      Hinv3(i4,i4)=2.0d0*(r1+poison)*rloc
      Hinv3(i5,i5)=2.0d0*(r1+poison)*rloc
      Hinv3(i6,i6)=2.0d0*(r1+poison)*rloc

c     Effective stress in compression at the end of the step
      do iloc=i1,i6
        parahot3(29+iloc,ntot) = parahot3(23+iloc,ntot)
      end do
c     Effective stress at the end of the step
      do iloc=i1,i6
        PARAHOT3(idimpara3-48+iloc,ntot) =
     .    PARAHOT3(idimpara3-i7+iloc,ntot)
      end do
      do iloc=i1,i6
        esigma6v(iloc) = PARAHOT3(idimpara3-i7+iloc,ntot)
      end do

c     eps,el = sigma,eff / C0     and not      sigma / C0      because damage
      call mulAB(Hinv3,esigma6v,trav3,i6,i6,i6,i1)
c                [6;6] x [6;1] = [6;1]
      do iloc=i1,i6
        EPSSIG6(iloc) = PARAHOT3(idimpara3-19+iloc,ntot) -
     .                     PARAHOT3((17+iloc),ntot)
        EPSIPL6(iloc) = EPSSIG6(iloc) - TRAV3(iloc)
      end do

c     Here is pseudo incremental plasticity for variable temperatures.
c     ----------------------------------------------------------------
c     The starting point is the plastic strain vector,
C       corresponding to a zero stress vector

c     initial strain = plastic strain
      do iloc=i1,i6
        PARAHOT3(idimpara3-25+iloc,ntot) = EPSIPL6(iloc)
      end do
c      if (lcp) then
c        parahot2(idimpara2-10,ntot) = EPSIPL6(1)
c        parahot2(idimpara2-9,ntot) = EPSIPL6(2)
c        parahot2(idimpara2-8,ntot) = EPSIPL6(4)
c      endif


      return
      end