anchor.diff:  !-------------------------------------------------------------------------------
   !> This subroutine provides to update stress for viscoelastic material
     subroutine UpdateViscoelastic(matl, sectType, eps, sig, vsig, dt, temp, tempn)
  +      use m_ElasticLinear
         TYPE( tMaterial ), INTENT(IN)   :: matl      !< material properties
         INTEGER, INTENT(IN)             :: sectType  !< not used currently
         REAL(KIND=kreal), INTENT(IN)    :: eps(6)    !< strain after this step
  
positive.code:       real(kind=kreal) :: G,Gg,K,Kg,Kth, exp_n,mu_0,mu_n,dq_n,dtau, theta
        real(kind=kreal) :: ina(1), outa(2), EE, PP, ddt
  
        real(kind=kreal) :: devstrain(6), en(6)
        type(tTable), pointer  :: dicval
        logical :: ierr
  
  
negative.code: 
        if(mtype==VISCOELASTIC+2) then   ! Arrhenius
  
          if( tn1==tn ) then
            hsn = mvar(2)*( 1.d0/(tn-mvar(3))-1.d0/(mvar(1)-mvar(3)) )/mvar(4)
            trsinc = dexp(hsn)
          else
            hsn1 = mvar(2)*( 1.d0/(tn1-mvar(3))-1.d0/(mvar(1)-mvar(3)) )/mvar(4)
            hsn = mvar(2)*( 1.d0/(tn-mvar(3))-1.d0/(mvar(1)-mvar(3)) )/mvar(4)
            asn1 = dexp(-hsn1)
            asn = dexp(-hsn)
            trsinc = (asn1-asn)/(hsn1-hsn)
          endif
  
        else   ! WLF
  
          if( mvar(3)+tn1-mvar(1)<=0.d0 ) then
             trsinc= -1.d0; return
          endif
          if( mvar(3)+tn-mvar(1)<=0.d0 ) then
             trsinc= -1.d0; return
          endif
          if( tn1==tn ) then
            hsn = mvar(2)*( tn-mvar(1) )/ (mvar(3)+tn-mvar(1) )*dlog(10.d0)
            trsinc = dexp(hsn)
          else
            hsn1 = mvar(2)*( tn1-mvar(1) )/ (mvar(3)+tn1-mvar(1) )*dlog(10.d0)
            hsn = mvar(2)*( tn-mvar(1) )/ (mvar(3)+tn-mvar(1) )*dlog(10.d0)
            asn1 = dexp(hsn1)
            asn = dexp(hsn)
            trsinc = (asn1-asn)/(hsn1-hsn)
          endif
  
        endif
  
    end function
  