SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,
     1 RPL,DDSDDT,DRPLDE,DRPLDT,
     2 STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME,
     3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT,
     4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC)
C
      INCLUDE 'ABA_PARAM.INC'
C
      CHARACTER*80 CMNAME
      DIMENSION STRESS(NTENS),STATEV(NSTATV),
     1 DDSDDE(NTENS,NTENS),
     2 DDSDDT(NTENS),DRPLDE(NTENS),
     3 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),
     4 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3)
C
C
      INTEGER i, j
      REAL Y, n, eta
      REAL lambda, thetalambda, mu, thetamu
      REAL A, B, C, D, F, G
C
      REAL, DIMENSION(NTENS,NTENS):: DDSDE
      REAL, DIMENSION(NTENS,NTENS):: DDSDS
      REAL, DIMENSION(NTENS):: SIGold
C
      Y = PROPS(1)
      n = PROPS(2)
      eta = PROPS(3) ! shear viscosity
C
      lambda = Y*n/((1.0d0+n)*(1.0d0-2.0d0*n))
      mu = Y/(2.0d0*(1.0d0+n))
      thetalambda = 3.0d0*eta/Y
      thetamu = eta/mu
C
      A = lambda*(1.0d0 + thetalambda/DTIME) +
     1    2.0d0*mu*(1.0d0 + thetamu/DTIME)
      B = lambda*(1.0d0 + thetalambda/DTIME)
      C = lambda + 2.0d0*mu
      D = lambda
      F = 2.0d0*mu*(1.0d0 + thetamu/DTIME)
      G = 2.0d0*mu
C
      DO i = 1,NTENS
        DO j = 1,NTENS
          DDSDDE(i,j) = 0.0d0
          DDSDE(i,j) = 0.0d0
          DDSDS(i,j) = 0.0d0
        END DO
      END DO
C
C DEFINE TENSOR COEFFICIENTS RELATED TO NORMAL STRESSES
C
      DO i = 1,NDI
         DO j = 1,NDI
            DDSDDE(i,j) = B
            DDSDE(i,j) = D
            DDSDS(i,j) = 0.0d0
         END DO
         DDSDDE(i,i) = A
         DDSDE(i,i) = C
         DDSDS(i,i) = -1.0d0
      END DO
C
C DEFINE TENSOR COEFFICIENTS RELATED TO SHEAR STRESSES
C
      DO i = NDI+1,NTENS
       DDSDDE(i,i) = F/2.0d0
       DDSDE(i,i) = G/2.0d0
       DDSDS(i,i) = - 1.0d0
      END DO
C
C STORE PRESENT STRESS STATE IN ARRAY SIGold
C
      DO i= 1,NTENS
        SIGold(i) = STRESS(i)
      END DO
C
C UPDATE STRESSES
C
      DO i = 1,NTENS
        DO j = 1,NTENS
          STRESS(i) = STRESS(i) + DDSDDE(i,j)*DSTRAN(j) +
     1                DDSDE(i,j)*STRAN(j) +
     2                DDSDS(i,j)*SIGold(j)
        END DO
      END DO
C
      RETURN
      END