Creating Your Own Material Model in Ansys

The UserMat subroutine enables you to write your own material constitutive equations within a general material framework using current-technology elements.

Overview

The UserMat subroutine defines the material stress-strain relationship of a material and applies to time-domain and full-harmonic analysis types.
The subroutine is called at every material integration point of the elements during the solution phase. The program passes in stresses, strains, and state variable values at the beginning of the increment (and for time-domain analyses, the strain increment at the current increment). The UserMat subroutine then updates the stresses and state variables to current values.

Time Domain

The two most important UserMat-calculated quantities are stress and the consistent tangent stiffness:
The measure of the stress is Cauchy stress (true stress), .
The consistent tangent stiffness is defined as , where is the stress increment and is the strain increment. The measure of strain is logarithmic strain (true strain), .
Generally, the material model can be formulated in either incremental (rate) or total forms.
The incremental form is suitable for material models such as plasticity. The program uses a co-rotated framework for such material formulation. All variables are defined and updated in the co-rotated material coordinate system. The implementation of a constitutive model looks exactly the same for both small and large deformation (NLGEOM,ON).
The total form is commonly used with hyperelastic material, and the material response is characterized by a potential relating to the invariants of the deformation gradient. The discrete equations are formulated with respect to the reference configuration. The deformation gradient is defined in the global coordinate system (with the exception of shell elements, where the deformation gradient is expressed on the co-rotated element coordinate system).

Harmonic

For full-harmonic analyses, the incoming strain argument is a two-dimensional vector containing the real (material stiffness) and imaginary (damping) parts of the harmonic strain. The subroutine must calculate and return the real and imaginary parts of the stress and material Jacobian matrix for the current frequency.
The complex stress returned by the subroutine has no effect on the harmonic solution and is used for postprocessing only.

Stress, Strain, and Material Jacobian Matrix

For nonlinear geometry analysis, the stress measure () used by the subroutine is the Cauchy stress (true stress), and the strain measure () is the logarithmic strain (true strain). The strains and incremental strains passed into UserMat are the total mechanical strains from which the thermal strains (if they exist) are subtracted.
UserMat must also provide the material Jacobian matrix defined as . is the stress increment, and is the strain increment.
UserMat is based on the current configuration for nonlinear geometry analysis (NLGEOM,ON). The program uses a co-rotational approach to account for rigid body rotation. Because the program already accounts for the strains passed into UserMat for the rigid body rotation, there is no need to apply additional rotation within UserMat.
Stress, strain, and the material Jacobian tensors are stored in a vector or matrix format.

The UserMat API

*deck,usermat                      USERDISTRIB
    subroutine usermat(
    &                   matId, elemId,kDomIntPt, kLayer, kSectPt,
    &                   ldstep,isubst,keycut,
    &                   nDirect,nShear,ncomp,nStatev,nProp,
    &                   Time,dTime,Temp,dTemp,
    &                   stress,ustatev,dsdePl,sedEl,sedPl,epseq,
    &                   Strain,dStrain, epsPl, prop, coords, 
    &                   var0, defGrad_t, defGrad,
    &                   tsstif, epsZZ,
    &                   cutFactor, pVolDer, hrmflg, var3, var4,
    &                   var5, var6, var7)
c*************************************************************************
c     *** primary function ***
c
c           user defined material constitutive model
c
c      Attention:
c           User must define material constitutive law properly
c           according to the stress state such as 3D, plane strain
c           and axisymmetry, plane stress and 3D/1D beam.
c
c           A 3D material constitutive model can be used for
c           plane strain and axisymmetry cases.
c
c           When using shell elements, a plane stress algorithm
c           must be used.
c
c                                             gal July, 1999
c
c       The following demonstrates a USERMAT subroutine for 
c       a plasticity model, which is the same as TB, BISO,
c       for different stress states. 
c       See "ANSYS user material subroutine USERMAT" for detailed
c       description of how to write a USERMAT routine.
c
c       This routine calls four routines,
c       usermat3d.F, usermatps.F usermatbm.F and usermat1d.F, w.r.t.
c       the corresponding stress states.
c       Each routine can be also a usermat routine for the specific 
c       element.
c
c*************************************************************************
c Copyright ANSYS.  All Rights Reserved.
c
c     input arguments
c     ===============
c      matId     (int,sc,i)               material #
c      elemId    (int,sc,i)               element #
c      kDomIntPt (int,sc,i)               "k"th domain integration point
c      kLayer    (int,sc,i)               "k"th layer
c      kSectPt   (int,sc,i)               "k"th Section point
c      ldstep    (int,sc,i)               load step number
c      isubst    (int,sc,i)               substep number
c      nDirect   (int,sc,in)              # of direct components
c      nShear    (int,sc,in)              # of shear components
c      ncomp     (int,sc,in)              nDirect + nShear
c      nstatev   (int,sc,i)               Number of state variables
c      nProp     (int,sc,i)               Number of material constants
c
c      Temp      (dp,sc,in)               temperature at beginning of
c                                         time increment
c      dTemp     (dp,sc,in)               temperature increment 
c      Time      (dp,sc,in)               time at beginning of increment (t)
c      dTime     (dp,sc,in)               current time increment (dt)
c
c      Strain   (dp,ar(ncomp),i)          Strain at beginning of time increment
c      dStrain  (dp,ar(ncomp),i)          Strain increment
c      prop     (dp,ar(nprop),i)          Material constants defined by TB,USER
c      coords   (dp,ar(3),i)              current coordinates
c      defGrad_t(dp,ar(3,3),i)            Deformation gradient at time t
c      defGrad  (dp,ar(3,3),i)            Deformation gradient at time t+dt
c      hrmflg   (dp,sc,io)                flag to indicate harmonic analysis 
c
c     input output arguments              
c     ======================             
c      stress   (dp,ar(ncomp),io)         stress
c      ustatev   (dp,ar(nstatev),io)      user state variables
c      sedEl    (dp,sc,io)                elastic work
c      sedPl    (dp,sc,io)                plastic work
c      epseq    (dp,sc,io)                equivalent plastic strain
c      epsPl   (dp,ar(ncomp),io)          plastic strain
c      var?     (dp,sc,io)                not used, they are reserved arguments 
c                                         for further development
c
c     output arguments
c     ================
c      keycut   (int,sc,o)                loading bisect/cut control
c                                         0 - no bisect/cut
c                                         1 - bisect/cut 
c                                         (factor will be determined by solution control)
c      dsdePl   (dp,ar(ncomp,ncomp),o)    material jacobian matrix
c      pVolDer  (dp,ar(3),o)              derivatives of volumetric potential wrt to J
c                                         pVolDer(1) = dU/dJ
c                                         pVolDer(2) = d^2U/dJ^2
c                                         pVolDer(3) = d^3U/dJ^3
c      tsstif   (dp,ar(2),o)              transverse shear stiffness
c                                         tsstif(1) - Gxz
c                                         tsstif(2) - Gyz
c                                         tsstif(1) is also used to calculate hourglass
c                                         stiffness, this value must be defined when low
c                                         order element, such as 181, 182, 185 with uniform 
c                                         integration is used.
c      epsZZ    (dp,sc,o)                 strain epsZZ for plane stress,
c                                         define it when accounting for thickness change
c                                         in shell and plane stress states
c      cutFactor(dp,sc,o)                 time step size cut-back factor 
c                                         define it if a smaller step size is wished
c                                         recommended value is 0~1
c
c*************************************************************************
c
c      ncomp   6   for 3D  (nshear=3)
c      ncomp   4   for plane strain or axisymmetric (nShear = 1)
c      ncomp   3   for plane stress (nShear = 1)
c      ncomp   3   for 3d beam      (nShear = 2)
c      ncomp   1   for 1D (nShear = 0)
c
c      stresses and strains, plastic strain vectors
c          11, 22, 33, 12, 23, 13    for 3D
c          11, 22, 33, 12            for plane strain or axisymmetry
c          11, 22, 12                for plane stress
c          11, 13, 12                for 3d beam
c          11                        for 1D
c
c      material jacobian matrix
c        3D
c           dsdePl    |  1111   1122   1133   1112   1123   1113 |
c           dsdePl    |  2211   2222   2233   2212   2223   2213 |
c           dsdePl    |  3311   3322   3333   3312   3323   3313 |
c           dsdePl    |  1211   1222   1233   1212   1223   1213 |
c           dsdePl    |  2311   2322   2333   2312   2323   2313 |
c           dsdePl    |  1311   1322   1333   1312   1323   1313 |
c        plane strain or axisymmetric (11, 22, 33, 12)
c           dsdePl    |  1111   1122   1133   1112 |
c           dsdePl    |  2211   2222   2233   2212 |
c           dsdePl    |  3311   3322   3333   3312 |
c           dsdePl    |  1211   1222   1233   1212 |
c        plane stress (11, 22, 12)
c           dsdePl    |  1111   1122   1112 |
c           dsdePl    |  2211   2222   2212 |
c           dsdePl    |  1211   1222   1212 |
c        3d beam (11, 13, 12)
c           dsdePl    |  1111   1113   1112 |
c           dsdePl    |  1311   1313   1312 |
c           dsdePl    |  1211   1213   1212 |
c        1d
c           dsdePl    |  1111 |
c
c*************************************************************************

Table(TB) Commands for UserMat

When creating your own material model, first define the material by specifying input values for the UserMat subroutine (TB,USER).

TB,USER Command

Issue the TB command using the following syntax:

TB,USER,MAT,NTEMPS,NPTS,TBOPT

where

MAT = User material ID number  
NTEMPS = Number of temperature points.  
NPTS = Number of material constants at a given temperature point.  
TBOPT = NONLINEAR (default), LINEAR, or MXUP  

The material properties at an intermediate temperature point are interpolated and passed to the UserMat subroutine.
If applicable, specify the number of state variables used (TB,STATE).
Define temperatures and material constants via TBTEMP and TBDATA commands, respectively.

Example

tb,user,1,2,4                   ! Define material 1 as a user
                                ! material with two temperatures
                                ! and four data points at each
                                ! temperature point.
tbtemp,1.0                      ! first temp. 
tbdata,1,19e5, 0.3, 1e3,100,    ! Four mat. constants for one temp.
tbtemp,2.0                      ! Second temp.
tbdata,1,21e5, 0.3, 2e3,100,    ! Four mat. constants for two temps.

TB,STATE Command

If you intend to use state variables with the UserMat subroutine, first specify the number of state variables. Issue the TB command using the following syntax:

TB,STATE,MAT, ,NPTS

where

MAT = User material ID number  
NPTS = Number of state variables that you intend to use (maximum 1000).  

The command defines only the number of state variables and must always be associated with a user material ID. No temperatures or data are associated with the command.
By default, the program initializes state variables to zero at the beginning of an analysis. Use the TBDATA command to initialize your own values for state variables.

Example

tb,state,1,,8                         ! Define material 1 with eight state variables
tbdata,1,c1,c2,c3,c4,c5,c6,c7,c8      ! Initialize the eight state variables.

Utility Functions for UserMat

The following functions are available for use with UserMat.

Initializes array a to zero.
The value n is the array dimension.

Multiplies vector a by constant c and outputs ( b = a * c) as vector b.
The value n is the dimension for both arrays.

Multiplies vector a by constant c and outputs the result as itself (that is, a = a * c).
The value n represents the array dimension.

Multiplies two double-precision matrices and outputs the result as c (that is, c = a * b ).
The value na is number of rows in matrix a , nb the number of rows in matrix b , and nc the number of rows in matrix c.
The n1 value is the number of rows in matrix c to fill, and n2 the number of columns in matrix c to fill.
The value n3 is the number of columns in matrix a and the number of rows in matrix b to work with. (The number of columns in a and rows in b is the same in order to generate the inner product correctly.)