VUMATHT

User subroutine to define a material's thermal behavior.

Warning: The use of this user subroutine generally requires considerable expertise. You are cautioned that the implementation of any realistic constitutive model requires extensive development and testing. Initial testing on a single-element model with prescribed traction loading is strongly recommended.

User subroutine VUMATHT:

can be used to define the thermal constitutive behavior of a material as well as internal heat generation during heat transfer processes;
will be called for blocks of material calculation points for which the thermal material behavior is defined in a user subroutine;
can use and update solution-dependent state variables;
can use any field variables that are passed in;
can be used with thermally coupled continuum elements and Eulerian elements; and
is described further in User-defined thermal material behavior.

The following topics are discussed:

User subroutine interface
Variables to be defined
Variables that can be updated
Variables passed in for information
Example: Using more than one user-defined material model
Example: A simple thermal material

User subroutine interface

      subroutine vumatht (
C Read only (unmodifiable) variables -
     *     nblock, nElem, nIntPt, nLayer, nSectPt, 
     *     ntgrad, nstatev, nfieldv, nprops,  
     *     cmname, stepTime, totalTime, dt,  
     *     coordMp, density, props,  
     *     tempOld, fieldOld, stateOld, enerThermOld, 
     *     tempNew, tempgradNew, fieldNew, 
C Write only (modifiable) variables -
     *     stateNew, fluxNew, enerThermNew, dEnerThDTemp, condEff )
C
      include 'vaba_param.inc'
C
      dimension nElem(nblock)
C
      dimension coordMp(nblock,*), density(nblock), props(nprops),
     *    tempOld(nblock), fieldOld(nblock, nfieldv), 
     *    stateOld(nblock, nstatev), enerThermOld(nblock),
     *    tempNew(nblock), tempgradNew(nblock, ntgrad), 
     *    fieldNew(nblock, nfieldv), stateNew(nblock, nstatev),
     *    fluxNew(nblock, ntgrad), enerThermNew(nblock), 
     *    dEnerThDTemp(nblock,2), condEff(nblock)
C
      character*80 cmname
C
      do km = 1, nblock
C         user coding
      end do
C
      return
      end

Variables to be defined

stateNew(nblock, nstatev): State variables at each material point at the end of the increment. You define the size of this array by allocating space for it (see About user subroutines and utilities for more information).
fluxNew(nblock, ntgrad): Heat flux vector at each material point at the end of the increment.
enerThermNew(nblock): Internal thermal energy per unit mass at each material point at the end of the increment.
dEnerThDTemp(nblock,1): Variation of internal energy per unit mass with respect to temperature at each material point at the end of the increment.
condEff(nblock): Effective conductivity at each material point at the end of the increment. This variable is required to estimate the element stable time increment size in Abaqus/Explicit. You can provide the largest value in conductivity tensor at each material point for a conservative estimate of the element stable time increment size.

Variables that can be updated

dEnerThDTemp(nblock,2): Variation of internal heat energy per unit mass resulting from internal heat generation (for example, latent heat release and absorption) with respect to the temperature at each material point at the end of the increment. This variable can be defined when a more conservative estimate of the element stable time increment size is needed during internal heat generation.

Variables passed in for information

nblock: Number of material points to be processed in this call to VUMATHT.
nElem(nblock): User-defined element number.
nIntPt: Integration point number.
nLayer: Layer number.
nSectPt: Section point number in the current layer.
ntgrad: Number of spatial gradients of temperature.
nstatev: Number of user-defined state variables that are associated with this material type (see Allocating space).
nfieldv: Number of user-defined external field variables.
nprops: User-specified number of user-defined material properties.
cmname: User-specified 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, you should not use “ABQ_” as the leading string for cmname.
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.
coordMp(nblock,*): Material point coordinates.
density(nblock): Current density at the material points in the midstep configuration.
props(nprops): User-supplied material properties.
tempOld(nblock): Temperatures at each material point at the beginning of the increment.
fieldOld(nblock, nfieldv): Values of the user-defined field variables at each material point at the beginning of the increment.
stateOld(nblock, nstatev): State variables at each material point at the beginning of the increment.
enerThermOld(nblock): Internal thermal energy per unit mass at each material point at the beginning of the increment.
tempNew(nblock): Temperatures at each material point at the end of the increment.
tempgradNew(nblock, ntgrad): Current values of the spatial gradients of temperature.
fieldNew(nblock, nfieldv): Values of the user-defined field variables at each material point at the end of the increment.

Example: Using more than one user-defined material model

To use more than one user-defined material model, the variable cmname can be tested for different material names inside user subroutine VUMATHT, as illustrated below:

if (cmname(1:4) .eq. 'MAT1') then
   call VUMATHT_MAT1(argument_list)
else if (cmname(1:4) .eq. 'MAT2') then
   call VUMATHT_MAT2(argument_list)
end if

VUMATHT_MAT1 and VUMATHT_MAT2 are the actual user material subroutines containing the thermal constitutive material models for each material MAT1 and MAT2, respectively. User subroutine VUMATHT merely acts as a directory here. The argument list can be the same as that used in user subroutine VUMATHT. The material names must be in uppercase characters since cmname is passed in as an uppercase character string.

Example: A simple thermal material

As a simple example of the coding of user subroutine VUMATHT, consider a material with isotropic conductivity and a constant specific heat. The heat transfer equations are developed here, and the example VUMATHT is given. The problem can also be solved by directly specifying thermal conductivity, specific heat, density, and internal heat generation.

First, the heat transfer equations are outlined.

The basic energy balance is

\int_{V} ρ \dot{U} d V = \int_{S} q d S + \int_{V} r d V,

where V is the volume of solid material with surface area S, $ρ$ is the density of the material, $\dot{U}$ is the material time rate of the internal thermal energy, q is the heat flux per unit area of the body flowing into the body, and r is the heat supplied externally into the body per unit volume. q and r can depend on the strains or displacements of the body.

In Abaqus/Explicit the heat transfer and mechanical solutions are obtained simultaneously by an explicit coupling. Thus, the material time rate of the internal thermal energy can be written as

\dot{U} = \frac{\partial U}{\partial θ} \dot{θ} .

A heat flux vector $f$ is defined such that

q = - f \cdot n,

where $n$ is the unit outward normal to the surface S. Introducing the above relation into the energy balance equation and using the divergence theorem, the following relation is obtained:

\int_{V} ρ \frac{\partial U}{\partial θ} \dot{θ} d V = - \int_{V} \frac{\partial}{\partial x} \cdot f d V + \int_{V} r d V .

The corresponding weak form is given by

\int_{V} δ θ ρ \frac{\partial U}{\partial θ} \dot{θ} d V = \int_{V} δ g \cdot f d V + \int_{V} δ θ r d V + \int_{S_{q}} δ θ q d S,

where

g = \frac{\partial θ}{\partial x}

is the temperature gradient and $δ θ$ is an arbitrary variational field satisfying the essential boundary conditions.

In Abaqus/Explicit this nonlinear system is solved using the explicit central-difference time integration rule with a lumped capacitance matrix as described in Fully coupled thermal-stress analysis.

The thermal constitutive behavior for this example is now defined. Isotropic conductivity and a constant specific heat for the material are assumed. The heat conduction in the material is assumed to be governed by Fourier's law.

The internal thermal energy per unit mass is defined as

U = \int c d θ = c (θ - θ^{Z}),

and its variation with respect to temperature is

\frac{\partial U}{\partial θ} = c,

where c is the specific heat of the material and $θ^{Z}$ is the value of absolute zero on the temperature scale being used.

Fourier's law for heat conduction is given as

f = - k \cdot g,

where $k$ is the thermal conductivity matrix.

An approximation to the stability limit for the forward-difference operator in the thermal solution response is given by

Δ t \approx L_{m i n}^{2} \frac{ρ c_{e}}{2 k_{e}},

where $L_{m i n}$ is the smallest element dimension in the mesh and the parameters $c_{e}$ and $k_{e}$ represent the material's effective specific heat and effective thermal conductivity, respectively. For a conservative estimate of the time increment size of the thermal solution response, we can specify

c_{e} = c,

k_{e} = \max (k_{i j}),

where $k_{i j}$ are components of the thermal conductivity matrix $k$ .

No state variables are needed for this material, so the allocation of space for them is not necessary.

A thermal user material definition can be used to read in the two constants for our simple case, namely the specific heat, c, and the coefficient of isotropic thermal conductivity, k, so that

PROPS(1) = k,

PROPS(2) = c .

      subroutine vumatht (
C Read only (unmodifiable) variables -
     *     nblock, nElem, nIntPt, nLayer, nSectPt, 
     *     ntgrad, nstatev, nfieldv, nprops,  
     *     cmname, stepTime, totalTime, dt,  
     *     coordMp, density, props,  
     *     tempOld, fieldOld, stateOld, enerThermOld, 
     *     tempNew, tempgradNew, fieldNew, 
C Write only (modifiable) variables -
     *     stateNew, fluxNew, enerThermNew, dEnerThDTemp, condEff )
C
      include 'vaba_param.inc'
C
      dimension nElem(nblock)
C
      dimension coordMp(nblock,*), density(nblock), props(nprops),
     *    tempOld(nblock), fieldOld(nblock, nfieldv), 
     *    stateOld(nblock, nstatev), enerThermOld(nblock),
     *    tempNew(nblock), tempgradNew(nblock, ntgrad), 
     *    fieldNew(nblock, nfieldv), stateNew(nblock, nstatev),
     *    fluxNew(nblock, ntgrad), enerThermNew(nblock), 
     *    dEnerThDTemp(nblock,2), condEff(nblock)
C
      character*80 cmname
*
      parameter ( zero = 0.d0 )
      parameter ( tempZero = 0.d0)
* properties
* 1. conductivity
* 2. specific heat
*
      cond = props(1)
      specHeat = props(2)
*
      do km = 1, nblock
*
*        update heat flux vector
         do i = 1, ntgrad
            fluxNew(km, i) = -cond*tempgradNew(km,i)
         end do
         condEff(km) = cond
*
*        update internal thermal energy and its derivatives
         if (totalTime .eq. zero) then
            dEnerTh = specHeat*(tempNew(km) - tempZero)
         else
            dEnerTh = specHeat*(tempNew(km) - tempOld(km))
         end if
         enerThermNew(km) = enerThermOld(km) + dEnerTh
         dEnerThDTemp(km,1) = specHeat
         dEnerThDTemp(km,2) = zero
      end do
*
      return
      end