## UANISOHYPER STRAIN help

## UANISOHYPER STRAIN help

(OP)

Hi Guys,

First time poster here. Im writing a UMAT to implement the fung elastic model into abaqus. The equation is U=c/2(e^^Q -1) and Q is int terms of the material coefficients of the strain components Eij.

Im using the UANISOHYPER Strain in abaqus, and what i need to define is U - the strain energy density function, UA(1)- Derivatives of strain energy potential with respect to the components of the modified Green strain tensor, UA(2)-2nd Derivatives of strain energy potential with respect to the components of the modified Green strain tensor.

Im modelling an incompressible soft tissue so the third derivative and the derivatives with respect to J ( volume ratio) are ignored.

I'm running a very simple job, 1 element, fixed at one end and being pulled at the opposite end. The umat compiles properly but i get an unexpected error 193 running the job. I have attached the UMAT if anyone has any experience writing something similar. Any help guys would be greatly appreciated.

SUBROUTINE UANISOHYPER_STRAIN (EBAR, AJ, UA, DU1, DU2, DU3,

1 TEMP, NOEL, CMNAME, INCMPFLAG, IHYBFLAG, NDI, NSHR, NTENS,

2 NUMSTATEV, STATEV, NUMFIELDV, FIELDV, FIELDVINC,

3 NUMPROPS, PROPS)

C

INCLUDE 'ABA_PARAM.INC'

C

DIMENSION EBAR(NTENS), UA(2), DU1(NTENS+1),

2 DU2((NTENS+1)*(NTENS+2)/2),

3 DU3((NTENS+1)*(NTENS+2)/2),

4 STATEV(NUMSTATEV), FIELDV(NUMFIELDV),

5 FIELDVINC(NUMFIELDV), PROPS(NUMPROPS)

C

DIMENSION D(3,3)

DIMENSION ET(3,3)

C where ET = elasticity tensor and D = term1*ET

CHARACTER*80 FUNG_HYPERELASTIC

C User Code to define

C

C UA(1)=Strain energy density function

C UA(2)=The deviatoric part of the strain energy desnity function

C DU1(NTENS+1)=Derivatives of the strain enery potential wrt to modified green strain

C DU2=Second Derivatives of the strain enery potential wrt to modified green strain

C Parameters

C PROPS(X) = Used to determine constant value and location in Abaqus

C Material Property Interface

zero = 0.d0

half = 0.5d0

one = 1.d0

two = 2.d0

twothirds = 2.d0/3.d0

fourthirds = 4.d0/3.d0

oneninth = 1.d0/9.d0

twoninths = 2.d0/9.d0

fourninths = 4.d0/9.d0

print*, 'in'

C Constants

C=PROPS(1)

A_1=PROPS(2)

A_2=PROPS(3)

A_3=PROPS(4)

A_4=PROPS(5)

A_5=PROPS(6)

A_6=PROPS(7)

C converting Green strain to modified green strain

C xpow= j**2/3= determination of deformation gradient (volume ratio)

xpow=exp(log(aj)*twothirds)

E11 = xpow * ebar(1) + half * ( xpow - one )

E22 = xpow * ebar(2) + half * ( xpow - one )

E33 = xpow * ebar(3) + half * ( xpow - one )

E12 = xpow * ebar(4)

E13 = xpow * ebar(5)

E23 = xpow * ebar(6)

C term 1 is c/2*expQ and term2 is expQ

term1=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22)))

term2=(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22)))

C let z1 - z3 be the derivative of the u substition of the first derivative.

C i.e z1 wrt to E11 z2 wrt E22 etc

z1=two*A_1*E11+two*A_3*E22+two*A_5*E12

z2=two*A_2*E22+two*A_3*E11+two*A_6*E12

z3=two*A_4*E12+two*A_5*E11+two*A_6*E22

C

ET(1,1)=2*A_1+(z1**2)

ET(1,2)=2*A_3+(z1*z2)

ET(2,2)=2*A_2+(z2**2)

ET(1,3)=2*A_5+(z1*z2)

ET(2,3)=2*A_6+(z2*z3)

ET(3,3)=2*A_4+(z3**2)

C Symmetry

ET(2,1)=ET(1,2)

ET(3,1)=ET(1,3)

ET(3,2)=ET(2,3)

C

D1111=term1*ET(1,1)

D1122=term1*ET(1,2)

D2222=term1*ET(2,2)

D1133=term1*ET(1,3)

D2233=term1*ET(2,3)

D3333=term1*ET(3,3)

D2211=D1122

D3311=D1133

D3322=D2233

C

C write(6,*)aj,twothirds,xpow

C

C write(6,*)NOEL,E11

C Strain Energy Density Function

UA(1)=U

U=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))-1)

C

C UA(2)=not needed

C

C let z1 - z3 be the derivative of the u substition of the first derivative.

C i.e z1 wrt to E11 z2 wrt E22 etc

z1=two*A_1*E11+two*A_3*E22+two*A_5*E12

z2=two*A_2*E22+two*A_3*E11+two*A_6*E12

z3=two*A_4*E12+two*A_5*E11+two*A_6*E22

C

dUde11=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))*(z1))

dUde22=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))*(z2))

dUde33=zero

dUde12=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))*(z3))

dUde13=zero

dUde23=zero

dUdeJ=zero

C

dUde11=du1(1)

dUde22=du1(2)

dUde33=du1(3)

dUde12=du1(4)

dUde13=du1(5)

dUde23=du1(6)

dUdeJ =du1(7)

C Second derivatives with respect to modified green strain tensors

C du2(X)

d2Ude11de11=dUde11*z1

d2Ude11de22=dUde11*z2

d2Ude22de22=dUde22*z2

d2Ude11de12=dUde11*z3

d2Ude22de12=dUde22*z3

d2Ude12de12=dUde12*z3

C

du2(indx(1,1)) = d2Ude11de11

du2(indx(1,2)) = d2Ude11de22

du2(indx(2,2)) = d2Ude22de22

du2(indx(1,3)) = zero

du2(indx(2,3)) = zero

du2(indx(3,3)) = zero

du2(indx(1,4)) = d2Ude11de12

du2(indx(2,4)) = d2Ude22de12

du2(indx(3,4)) = zero

du2(indx(4,4)) = d2Ude12de12

du2(indx(1,5)) = zero

du2(indx(2,5)) = zero

du2(indx(3,5)) = zero

du2(indx(4,5)) = zero

du2(indx(5,5)) = zero

return

end

* Maps index from Square to Triangular storage

* of symmetric matrix

*

integer function indx(i,j)

*

include 'aba_param.inc'

C write(6,*)'in'

*

ii = min(i,j)

jj = max(i,j)

*

indx = ii + jj*(jj-1)/2

C write(6,*)'out'

*

return

Thanks Guys

First time poster here. Im writing a UMAT to implement the fung elastic model into abaqus. The equation is U=c/2(e^^Q -1) and Q is int terms of the material coefficients of the strain components Eij.

Im using the UANISOHYPER Strain in abaqus, and what i need to define is U - the strain energy density function, UA(1)- Derivatives of strain energy potential with respect to the components of the modified Green strain tensor, UA(2)-2nd Derivatives of strain energy potential with respect to the components of the modified Green strain tensor.

Im modelling an incompressible soft tissue so the third derivative and the derivatives with respect to J ( volume ratio) are ignored.

I'm running a very simple job, 1 element, fixed at one end and being pulled at the opposite end. The umat compiles properly but i get an unexpected error 193 running the job. I have attached the UMAT if anyone has any experience writing something similar. Any help guys would be greatly appreciated.

SUBROUTINE UANISOHYPER_STRAIN (EBAR, AJ, UA, DU1, DU2, DU3,

1 TEMP, NOEL, CMNAME, INCMPFLAG, IHYBFLAG, NDI, NSHR, NTENS,

2 NUMSTATEV, STATEV, NUMFIELDV, FIELDV, FIELDVINC,

3 NUMPROPS, PROPS)

C

INCLUDE 'ABA_PARAM.INC'

C

DIMENSION EBAR(NTENS), UA(2), DU1(NTENS+1),

2 DU2((NTENS+1)*(NTENS+2)/2),

3 DU3((NTENS+1)*(NTENS+2)/2),

4 STATEV(NUMSTATEV), FIELDV(NUMFIELDV),

5 FIELDVINC(NUMFIELDV), PROPS(NUMPROPS)

C

DIMENSION D(3,3)

DIMENSION ET(3,3)

C where ET = elasticity tensor and D = term1*ET

CHARACTER*80 FUNG_HYPERELASTIC

C User Code to define

C

C UA(1)=Strain energy density function

C UA(2)=The deviatoric part of the strain energy desnity function

C DU1(NTENS+1)=Derivatives of the strain enery potential wrt to modified green strain

C DU2=Second Derivatives of the strain enery potential wrt to modified green strain

C Parameters

C PROPS(X) = Used to determine constant value and location in Abaqus

C Material Property Interface

zero = 0.d0

half = 0.5d0

one = 1.d0

two = 2.d0

twothirds = 2.d0/3.d0

fourthirds = 4.d0/3.d0

oneninth = 1.d0/9.d0

twoninths = 2.d0/9.d0

fourninths = 4.d0/9.d0

print*, 'in'

C Constants

C=PROPS(1)

A_1=PROPS(2)

A_2=PROPS(3)

A_3=PROPS(4)

A_4=PROPS(5)

A_5=PROPS(6)

A_6=PROPS(7)

C converting Green strain to modified green strain

C xpow= j**2/3= determination of deformation gradient (volume ratio)

xpow=exp(log(aj)*twothirds)

E11 = xpow * ebar(1) + half * ( xpow - one )

E22 = xpow * ebar(2) + half * ( xpow - one )

E33 = xpow * ebar(3) + half * ( xpow - one )

E12 = xpow * ebar(4)

E13 = xpow * ebar(5)

E23 = xpow * ebar(6)

C term 1 is c/2*expQ and term2 is expQ

term1=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22)))

term2=(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22)))

C let z1 - z3 be the derivative of the u substition of the first derivative.

C i.e z1 wrt to E11 z2 wrt E22 etc

z1=two*A_1*E11+two*A_3*E22+two*A_5*E12

z2=two*A_2*E22+two*A_3*E11+two*A_6*E12

z3=two*A_4*E12+two*A_5*E11+two*A_6*E22

C

ET(1,1)=2*A_1+(z1**2)

ET(1,2)=2*A_3+(z1*z2)

ET(2,2)=2*A_2+(z2**2)

ET(1,3)=2*A_5+(z1*z2)

ET(2,3)=2*A_6+(z2*z3)

ET(3,3)=2*A_4+(z3**2)

C Symmetry

ET(2,1)=ET(1,2)

ET(3,1)=ET(1,3)

ET(3,2)=ET(2,3)

C

D1111=term1*ET(1,1)

D1122=term1*ET(1,2)

D2222=term1*ET(2,2)

D1133=term1*ET(1,3)

D2233=term1*ET(2,3)

D3333=term1*ET(3,3)

D2211=D1122

D3311=D1133

D3322=D2233

C

C write(6,*)aj,twothirds,xpow

C

C write(6,*)NOEL,E11

C Strain Energy Density Function

UA(1)=U

U=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))-1)

C

C UA(2)=not needed

C

C let z1 - z3 be the derivative of the u substition of the first derivative.

C i.e z1 wrt to E11 z2 wrt E22 etc

z1=two*A_1*E11+two*A_3*E22+two*A_5*E12

z2=two*A_2*E22+two*A_3*E11+two*A_6*E12

z3=two*A_4*E12+two*A_5*E11+two*A_6*E22

C

dUde11=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))*(z1))

dUde22=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))*(z2))

dUde33=zero

dUde12=half*C*(exp((A_1*(E11)**2)+(A_2*(E22)**2)

* +(two*A_3*E11*E22)+(A_4*(E12)**2)

* +(two*A_5*E12*E11)+(two*A_6*E12*E22))*(z3))

dUde13=zero

dUde23=zero

dUdeJ=zero

C

dUde11=du1(1)

dUde22=du1(2)

dUde33=du1(3)

dUde12=du1(4)

dUde13=du1(5)

dUde23=du1(6)

dUdeJ =du1(7)

C Second derivatives with respect to modified green strain tensors

C du2(X)

d2Ude11de11=dUde11*z1

d2Ude11de22=dUde11*z2

d2Ude22de22=dUde22*z2

d2Ude11de12=dUde11*z3

d2Ude22de12=dUde22*z3

d2Ude12de12=dUde12*z3

C

du2(indx(1,1)) = d2Ude11de11

du2(indx(1,2)) = d2Ude11de22

du2(indx(2,2)) = d2Ude22de22

du2(indx(1,3)) = zero

du2(indx(2,3)) = zero

du2(indx(3,3)) = zero

du2(indx(1,4)) = d2Ude11de12

du2(indx(2,4)) = d2Ude22de12

du2(indx(3,4)) = zero

du2(indx(4,4)) = d2Ude12de12

du2(indx(1,5)) = zero

du2(indx(2,5)) = zero

du2(indx(3,5)) = zero

du2(indx(4,5)) = zero

du2(indx(5,5)) = zero

return

end

* Maps index from Square to Triangular storage

* of symmetric matrix

*

integer function indx(i,j)

*

include 'aba_param.inc'

C write(6,*)'in'

*

ii = min(i,j)

jj = max(i,j)

*

indx = ii + jj*(jj-1)/2

C write(6,*)'out'

*

return

Thanks Guys

## RE: UANISOHYPER STRAIN help

## RE: UANISOHYPER STRAIN help

I suppose that your problems have been resolved,

however, I am doing the same subroutine in Abaqus for the first time (I am a Ansys user)

-) I am not sure if this is correct:

dUde11=du1(1)

dUde22=du1(2)

dUde33=du1(3)

dUde12=du1(4)

dUde13=du1(5)

dUde23=du1(6)

dUdeJ =du1(7)

why not the inverse?...

du1(1)=dUde11

...

...

-) additionally, did you try to define UA2 and don't use it?

These my two hypotheses arise from the example "UANISOHYPER_STRAIN" reported in the Abaqus manual

best

## RE: UANISOHYPER STRAIN help

Cheers for the reply, I'm still stuck on this! I noticed your first correction straight away and changed it as you said, but still to no avail.

I didnt use UA2. My reason being the documentation says this is the deviatoric part of the strain energy function,I however am just differentiating the whole strain energy function. Perhaps this could be a source error?

My umat is running and I'm using print statements "in" at the start of the code and "out" at the end of the code and I can see that the umat is being called 10 times after which an "unexpected error 193" appears. I am also printing my arrays and none of them are filling.

How are you getting on implementing it into Ansys? There is an existing fung model in abaqus if you dont already know, may be of use to you. I need to changed it slightly hence why im doing a umat.

Cheers for the tips, any other help would be great

## RE: UANISOHYPER STRAIN help

Is there any particular reason why ebar is not updating from its original value of 0?

## RE: UANISOHYPER STRAIN help

Ebar probably isn't updating because you're probably not returning realistic strain energy density values to Abaqus, they might be very small etc.