!DEC$ FREEFORM
!------------------------------------------------------------------------
! User hardening subroutine for multilinear hardening curves
!   according to UHARD interface of Abaqus
!   
!   Input PROPS(*): Table of yield stress and values of eq. plastic strain in ascending order 
!
! distributed under CC BY-NC-SA 4.0 license
! (https://creativecommons.org/licenses/by-nc-sa/4.0/)
! Reference: 
!   N. Lange, G. Huetter, B. Kiefer: "An efficient monolithic solution scheme for FE2 problems",
!   https://arxiv.org/abs/2101.01802
!
! Nils Lange, 2020-11-11
! modified: Geralf Hütter, 2021-02-05
!------------------------------------------------------------------------ 
      SUBROUTINE UHARD(SYIELD,HARD,EQPLAS,EQPLASRT,TIME,DTIME,TEMP,&
                       DTEMP,NOEL,NPT,LAYER,KSPT,KSTEP,KINC,CMNAME,NSTATV,&
                       STATEV,NUMFIELDV,PREDEF,DPRED,NUMPROPS,PROPS)
!
      USE ABQINTERFACE
!
      IMPLICIT NONE
!
      CHARACTER(len=80):: CMNAME
      REAL(KIND=AbqRK)::EQPLAS,EQPLASRT,DTIME,TEMP,DTEMP,TIME(2)
      REAL(KIND=AbqRK)::SYIELD,HARD(3)
      INTEGER(KIND=AbqIK)::NSTATV,NOEL,NPT,LAYER,KSPT,KSTEP,KINC,NUMPROPS,NUMFIELDV
      REAL(KIND=AbqRK)::PREDEF(NUMFIELDV),DPRED(NUMFIELDV),PROPS(NUMFIELDV)
      REAL(KIND=AbqRK)::STATEV(NSTATV)
!     Internal variables of UHARD	  
	  REAL(KIND=AbqRK)::SYIEL0,EQPL1,EQPL0
      REAL(KIND=AbqRK),DIMENSION(2,NUMPROPS/2)::TABLE
      INTEGER::I,NENTRIES
	  LOGICAL::FOUND
!
      IF (NUMPROPS==1) THEN
          SYIELD=PROPS(1)
          HARD(1)=0.0
          RETURN
      END IF
!      
	  NENTRIES=NUMPROPS/2  
	  IF (MOD(NUMPROPS,2)==1) WRITE(*,*) 'Odd number of UHARD properties specified. Last entry is ignored.'    
      TABLE=RESHAPE(PROPS,[2,NENTRIES])
!     YIELD STRESS = LAST VALUE OF TABLE
      SYIELD=TABLE(1,NENTRIES)
      HARD(1)=0!
!     SEARCH AND INTERPOLATE TABLE, IF EXISTING
      IF(NENTRIES>1) THEN
	    I=1
!       FIRST EQPL ALWAYS SET TO ZERO, EVEN IF SPECIFIED INCORRECTLY		
		EQPL1=0 
		FOUND=.FALSE.
	    DO WHILE (.NOT.(FOUND).AND.(I<NENTRIES))
		  I=I+1
		  EQPL0=EQPL1
		  EQPL1=TABLE(2,I)
		  IF(EQPL1<=EQPL0) THEN
                WRITE(6,*) '***ERROR in UHARD: plastic strains not entered in ascending order!'
                CALL XIT
          ENDIF
		  FOUND=(EQPLAS<EQPL1)
		END DO
		IF (FOUND) THEN 
	      SYIEL0=TABLE(1,I-1)
          HARD(1)=(TABLE(1,I)-SYIEL0)/(EQPL1-EQPL0)
          SYIELD=SYIEL0+(EQPLAS-EQPL0)*HARD(1)
        END IF
      END IF
!
      RETURN
      END SUBROUTINE UHARD
