VUCREEPNETWORK

Model Description

The user subroutine allows a creep law of the following general form to be defined:

{\dot{\bar{ε}}}^{c r} = g^{c r} ({\bar{ε}}^{c r}, I_{1}^{c r}, {\bar{I}}_{1}, {\bar{I}}_{2}, J, p, q, t, θ, F V),

where

I_{1}^{c r} = I : C^{c r},

and

$I$: is the identity tensor,
$C^{c r}$: is the right Cauchy-Green creep strain tensor,
${\dot{\bar{ε}}}^{c r}$: is the equivalent creep strain rate,
${\bar{ε}}^{c r}$: is the equivalent creep strain,
${\bar{I}}_{1}$: is the first invariant of $\bar{B}$ ,
${\bar{I}}_{2}$: is the second invariant of $\bar{B}$ ,
$J$: is the determinant of the deformation gradient, $F$ ,
$p$: is the Kirchhoff pressure,
$q$: is the equivalent deviatoric Kirchhoff stress,
t: is the time,
$θ$: is the temperature, and
FV: are field variables.

The left Cauchy-Green strain tensor, $\bar{B}$ , is defined as

\bar{B} = \bar{F} {\bar{F}}^{T},

where $\bar{F}$ is the deformation gradient with volume change eliminated, which is computed using

\bar{F} = J^{- \frac{1}{3}} F .

The user subroutine must define the increment of creep equivalent strain, $Δ {\bar{ε}}^{c r}$ , as a function of the time increment, $Δ t$ , and the variables used in the definition of $g^{c r}$ , as well as the derivatives of the equivalent creep strain increment with respect to those variables. If any solution-dependent state variables are included in the definition of $g^{c r}$ , they must also be integrated forward in time in this user subroutine.

User Subroutine Interface

     subroutine vucreepnetwork (
C Read only -
     *   nblock, networkid, nstatev, nfieldv,
     *   nprops, nDg, stepTime, totalTime, dt,
     *   jElem, kIntPt, kLayer, kSecPt, cmname,
     *   props, coordMp, tempOld, fieldOld,
     *   stateOld, tempNew, fieldNew,
     *   nIarray, i_array, nRarray, r_array,
     *   q, p, eqcs, TrCc,
C Write only -
     *   dg, stateNew )
C
      include 'vaba_param.inc'
C
C     indices for equivalent creep strain and its derivatives
      parameter( i_deqcs       = 1,
     *           i_DdeqcsDq    = 2,
     *           i_DdeqcsDeqcs = 3,
     *           i_DdeqcsDi1c  = 4 )
C
C     indices for strain invariants
      parameter( i_I1 = 1,
     *           i_I2 = 2,
     *           i_J  = 3 )
C
      dimension props(nprops),
     *  tempOld(nblock),
     *  fieldOld(nblock,nfieldv),
     *  stateOld(nblock,nstatev),
     *  tempNew(nblock),
     *  fieldNew(nblock,nfieldv),
     *  coordMp(nblock,*),
     *  jElem(nblock),
     *  i_array(nblock,nIarray),
     *  r_array(nblock,nRarray),
     *  q(nblock),
     *  p(nblock),
     *  eqcs(nblock),
     *  TrCc(nblock),
     *  stateNew(nblock,nstatev),
     *  dg(nblock,nDg)

      character*80 cmname
C
      do 100 km = 1,nblock
        user coding
  100 continue

      return
      end

Variables to Be Defined

dg(nblock,i_deqcs): Equivalent creep strain increment, $Δ {\bar{ε}}^{c r}$ .
dg(nblock,i_DdeqcsDq): The derivative: $\partial Δ {\bar{ε}}^{c r} / \partial q$ .
dg(nblock,i_DdeqcsDeqcs): The derivative: $\partial Δ {\bar{ε}}^{c r} / \partial {\bar{ε}}^{c r}$ .
dg(nblock,i_DdeqcsDi1c): The first invariant, $I_{1}^{c r}$ , of the right Cauchy-Green creep strain tensor, $C^{cr}$ .

Variables That Can Be Updated

stateNew(nblock,nstatev): Array containing the user-defined solution-dependent state variables at this point.

Variables Passed in for Information

nblock: Number of material points to be processed in this call to user subroutine VUCREEPNETWORK.
networkid: Network identification number, which identifies the network for which creep is defined.
nstatev: Number of user-defined state variables that are associated with this material type (see Allocating Space for Solution-Dependent State Variables).
nfieldv: Number of user-defined external field variables.
nprops: User-specified number of user-defined material properties.
nDg: Size of array dg.
stepTime: Value of time since the step began.
totalTime: Value of total time. The time at the beginning of the step is given by totalTime-stepTime.
dt: Time increment size.
jElem(nblock): Array of element numbers.
kIntPt: Integration point number.
kLayer: Layer number (for composite shells).
kSecPt: Section point number within the current layer.
cmname: Material name, left justified. It is passed in as an uppercase character string. Some internal material models are given names starting with the “ABQ_” character string. To avoid conflict, “ABQ_” should not be used as the leading string for cmname.
props(nprops): User-supplied material properties.
coordMp(nblock,*): Material point coordinates. It is the midplane material point for shell elements and the centroid for beam elements.
tempOld(nblock): Temperatures at the material points at the beginning of the increment.
fieldOld(nblock,nfieldv): Values of the user-defined field variables at the material points at the beginning of the increment.
stateOld(nblock,nstatev): State variables at the material points at the beginning of the increment.
tempNew(nblock): Temperatures at the material points at the end of the increment.
fieldNew(nblock): Values of the user-defined field variables at the material points at the end of the increment.
nIarray: Size of array i_array.
i_array(nblock,nIarray): Array containing integer arguments. Currently it is not used.
nRarray: Size of array r_array.
r_array(nblock,i_I1): The first invariant, $\bar{I_{1}}$ , of the left Cauchy-Green strain tensor, $\bar{B}$ .
r_array(nblock,i_I2): The second invariant, $\bar{I_{2}}$ , of the left Cauchy-Green strain tensor, $\bar{B}$ .
r_array(nblock,i_J): The determinant of the deformation gradient, $F$ .
q(nblock): Array containing equivalent deviatoric Kirchhoff stresses.
p(nblock): Array containing Kirchhoff pressures.
eqcs(nblock): Array containing equivalent creep strains.
TrCc(nblock): Array containing the first invariants, $I_{1}^{c r}$ , of the right Cauchy-Green creep strain tensor, $C^{cr}$ .

Example: Power-Law Strain Hardening Model

As an example of the coding of user subroutine VUCREEPNETWORK, consider the power-law strain hardening model. In this case the equivalent creep strain rate is expressed as

{\dot{\bar{ε}}}^{c r} = {(A q^{n} {[(m + 1) {\bar{ε}}^{c r}]}^{m})}^{\frac{1}{m + 1}},

where

${\bar{ε}}^{c r}$: is the equivalent creep strain,
$q$: is the equivalent deviatoric Kirchhoff stress, and
A, m, and n: are material parameters.

The user subroutine would be coded as follows:

     subroutine vucreepnetwork (
C Read only -
     *   nblock, networkid, nstatev, nfieldv,
     *   nprops, nDg, stepTime, totalTime, dt,
     *   jElem, kIntPt, kLayer, kSecPt, cmname,
     *   props, coordMp, tempOld, fieldOld,
     *   stateOld, tempNew, fieldNew,
     *   nIarray, i_array, nRarray, r_array,
     *   q, p, eqcs, TrCc,
C Write only -
     *   dg, stateNew )
C
      include 'vaba_param.inc'
C
      parameter ( one = 1.d0, half = 0.5d0 )
      parameter ( eqcsSmall = 1.d-8 )
      parameter ( rMinVal = 1.d-12 )
C
      parameter( i_deqcs       = 1,
     *           i_DdeqcsDq    = 2,
     *           i_DdeqcsDeqcs = 3,
     *           i_DdeqcsDi1c  = 4 )
C
      dimension props(nprops),
     *  tempOld(nblock),
     *  fieldOld(nblock,nfieldv),
     *  stateOld(nblock,nstatev),
     *  tempNew(nblock),
     *  fieldNew(nblock,nfieldv),
     *  coordMp(nblock,*),
     *  jElem(nblock),
     *  i_array(nblock,nIarray),
     *  r_array(nblock,nRarray),
     *  q(nblock),
     *  p(nblock),
     *  eqcs(nblock),
     *  TrCc(nblock),
     *  stateNew(nblock,nstatev),
     *  dg(nblock,nDg)
C
      character*80 cmname
C
C Read properties
C
      rA = props(1)
      rN = props(2)
      rM = props(3)
C
C Update equivalent creep strain and its derivatives
C
      do k = 1, nblock
        om1 = one / ( one + rM )
        test = half - sign( half, q(k) - rMinVal )
        qInv = ( one - test ) / ( q(k) + test )
        eqcs_t = eqcs(k)
        if ( eqcs_t .le. eqcsSmall .and. q(k).gt.rMinVal ) then 
C Initial guess based on constant creep strain rate during increment
          eqcs_t = dt*(exp(log(rA)+rN*log(q(k)))*
     *             ((one+rM)*dt)**rM)
        end if
C
        test2 = half - sign( half, eqcs_t - rMinVal )
        eqcsInv = ( one - test2 ) / ( eqcs_t + test2 )
        g = dt*(exp(log(rA)+rN*log(q(k)))*
     *       ((one+rM)*(test2+eqcs_t))**rM)**om1
        dg(k,i_deqcs) = g
        dg(k,i_DdeqcsDq) = qInv * rN * om1 * g
        dg(k,i_DdeqcsDeqcs) = eqcsInv * rM * om1 * g
      end do
C
      return
      end