# DUDG(NTGRD) Variation of internal thermal energy per unit mass with respect to the spatial gradients of temperature, $\partial U / \partial ( \partial \theta / \partial \mathbf { x } )$ , at the end of the increment. The size of this array depends on the value of NTGRD as defined below. This term is typically zero in classical heat transfer analysis. # FLUX(NTGRD) Heat flux vector, , at the end of the increment. This variable is passed in with the values at the beginning of the increment and must be updated to the values at the end of the increment. # DFDT(NTGRD) Variation of the heat flux vector with respect to temperature, , evaluated at the end of the increment. # DFDG(NTGRD,NTGRD) Variation of the heat flux vector with respect to the spatial gradients of temperature, , at the end of the increment. The size of this array depends on the value of NTGRD as defined below. # Variables that can be updated # STATEV(NSTATV) An array containing the solution-dependent state variables. In an uncoupled heat transfer analysis STATEV is passed into UMATHT with the values of these variables at the beginning of the increment. However, any changes in STATEV made in user subroutine USDFLD will be included in the values passed into UMATHT, since USDFLD is called before UMATHT. In addition, if UMATHT is being used in a fully coupled temperature-displacement or coupled thermalelectrical-structural analysis and user subroutine CREEP, user subroutine UEXPAN, user subroutine UMAT, or user subroutine UTRS is used to define the mechanical behavior of the material, those routines are called before this routine; therefore, any updating of STATEV done in CREEP, UEXPAN, UMAT, or UTRS will be included in the values passed into UMATHT. In all cases STATEV should be passed back from UMATHT as the values of the state variables at the end of the current increment. # PNEWDT Ratio of suggested new time increment to the time increment being used (DTIME, see below). This variable allows you to provide input to the automatic time incrementation algorithms in Abaqus/Standard (if automatic time incrementation is chosen). PNEWDT is set to a large value before each call to UMATHT. If PNEWDT is redefined to be less than 1.0, Abaqus/Standard must abandon the time increment and attempt it again with a smaller time increment. The suggested new time increment provided to the automatic time integration algorithms is PNEWDT × DTIME, where the PNEWDT used is the minimum value for all calls to user subroutines that allow redefinition of PNEWDT for this iteration. If PNEWDT is given a value that is greater than 1.0 for all calls to user subroutines for this iteration and the increment converges in this iteration, Abaqus/Standard may increase the time increment. The suggested new time increment provided to the automatic time integration algorithms is PNEWDT × DTIME, where the PNEWDT used is the minimum value for all calls to user subroutines for this iteration. If automatic time incrementation is not selected in the analysis procedure, values of PNEWDT that are greater than 1.0 will be ignored and values of PNEWDT that are less than 1.0 will cause the job to terminate. # Variables passed in for information # TEMP Temperature at the start of the increment. # DTEMP Increment of temperature. # DTEMDX(NTGRD) Current values of the spatial gradients of temperature, # TIME(1) Value of step time at the beginning of the current increment. # TIME(2) Value of total time at the beginning of the current increment. # DTIME Time increment. # PREDEF Array of interpolated values of predefined field variables at this point at the start of the increment, based on the values read in at the nodes. # DPRED Array of increments of predefined field variables. # CMNAME User-defined material name, left justified. # NTGRD Number of spatial gradients of temperature. # NSTATV Number of solution-dependent state variables associated with this material type (defined as described in “Allocating space” in “User subroutines: overview,” Section 18.1.1 of the Abaqus Analysis User’s Guide). # PROPS(NPROPS) User-specified array of material constants associated with this user material. # NPROPS User-defined number of material constants associated with this user material. # COORDS An array containing the coordinates of this point. These are the current coordinates in a fully coupled temperature-displacement or coupled thermal-electrical-structural analysis if geometric nonlinearity is accounted for during the step (see “Defining an analysis,” Section 6.1.2 of the Abaqus Analysis User’s Guide); otherwise, the array contains the original coordinates of the point. # NOEL Element number. # NPT Integration point number. # LAYER Layer number (for composite shells and layered solids). # KSPT Section point number within the current layer. # KSTEP Step number. # KINC Increment number. # Example: Using more than one user-defined thermal material model To use more than one user-defined thermal material model, the variable CMNAME can be tested for different material names inside user subroutine UMATHT, as illustrated below: ```sql IF (CMNAME(1:4) .EQ. 'MAT1') THEN CALL UMATHT_MAT1 (argument_list) ELSE IF(CMNAME(1:4) .EQ. 'MAT2') THEN CALL UMATHT_MAT2 (argument_list) END IF ``` UMATHT\_MAT1 and UMATHT\_MAT2 are the actual user material subroutines containing the constitutive material models for each material MAT1 and MAT2, respectively. Subroutine UMATHT merely acts as a directory here. The argument list can be the same as that used in subroutine UMATHT. # Example: Uncoupled heat transfer As a simple example of the coding of user subroutine UMATHT, consider uncoupled heat transfer analysis in a material. The equations for this case are developed here, and the corresponding UMATHT is given. This problem can also be solved by specifying thermal conductivity, specific heat, density, and internal heat generation directly. First, the equations for an uncoupled heat transfer analysis are outlined. The basic energy balance is $$ \int_ {V} \rho \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 , \rho$ 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. A heat flux vector is defined such that $$ q = - \mathbf {f} \cdot \mathbf {n}, $$ where 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} \rho \dot {U} d V = - \int_ {V} \frac {\partial}{\partial \mathbf {x}} \cdot \mathbf {f} d V + \int_ {V} r d V. $$ The corresponding weak form is given by $$ \int_ {V} \delta \theta \rho \dot {U} d V - \int_ {V} \delta \mathbf {g} \cdot \mathbf {f} d V = \int_ {V} \delta \theta r d V + \int_ {S} \delta \theta q d S, $$ where $$ \mathbf {g} = \frac {\partial \theta}{\partial \mathbf {x}} $$ is the temperature gradient and is an arbitrary variational field satisfying the essential boundary conditions. Introducing the backward difference integration algorithm: $$ \dot {U} _ {t + \Delta t} = (U _ {t + \Delta t} - U _ {t}) (1 / \Delta t), $$ the weak form of the energy balance equation becomes $$ \frac {1}{\Delta t} \int_ {v} \delta \theta \rho (U _ {t + \Delta t} - U _ {t}) d V = \int_ {V} \delta \mathbf {g} \cdot \mathbf {f} d V + \int_ {V} \delta \theta r d V + \int_ {S} \delta \theta q d S. $$ This nonlinear system is solved using Newton’s method. In the above equations the thermal constitutive behavior of the material is given by $$ U = U (\theta , t, \partial \theta / \partial \mathbf {x}, s ^ {i}, \dots) U = U (\theta , \mathbf {g}, t, s ^ {i}, \dots) $$ and $$ \mathbf {f} = \mathbf {f} (\theta , t, \partial \theta / \partial \mathbf {x}, s ^ {i}, \dots), \mathbf {f} = \mathbf {f} (\theta , \mathbf {g}, t, s ^ {i}, \dots), $$ where $s ^ { i }$ are state variables. The Jacobian for Newton’s method is given by (after dropping the subscripts $t + \Delta t$ on U) $$ \begin{array}{l} \frac {1}{\Delta t} \int_ {V} \delta \theta \rho \frac {\partial U}{\partial \theta} d \theta d V + \frac {1}{\Delta t} \int_ {V} \delta \theta \rho \frac {\partial U}{\partial \mathbf {g}} \cdot d \mathbf {g} d V \\ - \int_ {V} \delta \mathbf {g} \cdot \frac {\partial \mathbf {f}}{\partial \theta} d \theta d V - \int_ {V} \delta \mathbf {g} \cdot \frac {\partial \mathbf {f}}{\partial \mathbf {g}} \cdot d \mathbf {g} d V \\ - \int_ {V} \delta \theta \frac {\partial r}{\partial \theta} d \theta d V - \int_ {S} \delta \theta \frac {\partial q}{\partial \theta} d \theta d S. \\ \end{array} $$ The thermal constitutive behavior for this example is now defined. We assume a constant specific heat for the material. 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 = U (\theta), $$ with $$ \frac {\partial U}{\partial \theta} = c, $$ where c is the specific heat of the material and $$ \frac {\partial U}{\partial \mathbf {g}} = 0. $$ Fourier’s law for heat conduction is given as $$ \mathbf {f} = - \mathbf {k} \cdot \mathbf {g}, $$ where is the thermal conductivity matrix and is position, so that $$ \frac {\partial \mathbf {f}}{\partial \mathbf {g}} = - \mathbf {k} $$ and $$ \frac {\partial \mathbf {f}}{\partial \theta} = - \frac {\partial \mathbf {k}}{\partial \theta} \cdot \mathbf {g}. $$ The assumption of conductivity without any temperature dependence implies that $$ \frac {\partial \mathbf {f}}{\partial \theta} = 0. $$ 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 thermal conductivity, k, so that $$ \operatorname{PROPS} (1) = k, $$ $$ \operatorname{PROPS} (2) = c. $$```fortran SUBROUTINE UMATHT(U, DUDT, DUDG, FLUX, DFT, DFDG, 1 STATEV, TEMP, DTEMP, DTEMDX, TIME, DTIME, PREDEF, DPRED, 2 CMNAME, NTGRD, NSTATV, PROPS, NPROPS, COORDS, PNEWDT, 3 NOEL, NPT, LAYER, KSPT, KSTEP, KINC) C INCLUDE 'ABA_PARAM.INC' C CHARACTER*80 CMNAME DIMENSION DUDG(NTGRD), FLUX(NTGRD), DFT(NTGRD), 1 DFDG(NTGRD, NTGRD), STATEV(NSTATV), DTEMDX(NTGRD), 2 TIME(2), PREDEF(1), DPRED(1), PROPS(NPROPS), COORDS(3) C COND = PROPS(1) SPECHT = PROPS(2) C DUDT = SPECHT DU = DUDT*DTEMP U = U+DU C DO I=1, NTGRD FLUX(I) = -COND*DTEMDX(I) DFDG(I, I) = -COND END DO C RETURN END ``` # 1.1.46 UMESHMOTION: User subroutine to specify mesh motion constraints during adaptive meshing. Product: Abaqus/Standard # References • “Defining ALE adaptive mesh domains in Abaqus/Standard,” Section 12.2.6 of the Abaqus Analysis User’s Guide • \*ADAPTIVE MESH • \*ADAPTIVE MESH CONSTRAINT # Overview User subroutine UMESHMOTION: • is called at the end of any increment where adaptive meshing is performed (as specified by the frequency in increments); • can be used to define the motion of nodes in an adaptive mesh constraint node set; and • can call utility routines GETVRN, GETNODETOELEMCONN, and GETVRMAVGATNODE to access results data at the node. # Accessing node point data You are provided with access to the values of the node point quantities at the end of the increment through the utility routine GETVRN described in “Obtaining node point information,” Section 2.1.9. You can also access values of material point quantities extrapolated to, and averaged, at nodes at the end of the increment through the utility routine GETVRMAVGATNODE described in “Obtaining material point information averaged at a node,” Section 2.1.8. GETVRMAVGATNODE requires the list of elements attached to the node, which is obtained by calling the utility routine GETNODETOELEMCONN described in “Obtaining node to element connectivity,” Section 2.1.10. # User subroutine interface ```txt SUBROUTINE UMESHMOTION(UREF, ULOCAL, NODE, NNDOF, * LNODETYPE, ALOCAL, NDIM, TIME, DTIME, PNEWDT, * KSTEP, KINC, KMESHSWEEP, JMATYP, JGVBLOCK, LSMOOTH) C INCLUDE 'ABA_PARAM.INC' C DIMENSION ULOCAL(NDIM), JELEMLIST(*) DIMENSION ALOCAL(NDIM, *), TIME(2) ``` DIMENSION JMATYP(\*),JGVBLOCK(\*) C user coding to define ULOCAL and, optionally PNEWDT RETURN END # Variable to be defined # ULOCAL Components of the mesh displacement or velocity of the adaptive mesh constraint node, described in the coordinate system ALOCAL. ULOCAL will be passed into the routine as values determined by the mesh smoothing algorithm. All components of the mesh displacement or velocity will be applied; i.e., you do not have the ability to select the directions in which the mesh displacement should be applied. # Variables that can be updated # PNEWDT Ratio of suggested new time increment to the time increment currently being used (DTIME, see below). This variable allows you to provide input to the automatic time incrementation algorithms in Abaqus/Standard (if automatic time incrementation is chosen). PNEWDT is set to a large value before each call to UMESHMOTION. The suggested new time increment provided to the automatic time integration algorithms is PNEWDT × DTIME, where the PNEWDT used is the minimum value for all calls to user subroutines that allow redefinition of PNEWDT for this increment. If automatic time incrementation is not selected in the analysis procedure, values of PNEWDT greater than 1.0 will be ignored and values of PNEWDT less than 1.0 will cause the job to terminate. # LSMOOTH Flag specifying that surface smoothing be applied after application of the mesh motion constraint. Set LSMOOTH to 1 to enable surface smoothing. When this flag is set, the constraint defined in ULOCAL will be modified by the smoothing algorithm. In cases where ULOCAL describes mesh motion normal to a surface, the smoothing will have a minor impact on this normal component of mesh motion. # Variables passed in for information # UREF The value of the user-specified displacement or velocity provided as part of the adaptive mesh constraint definition. This value is updated based on any amplitude definitions used with the adaptive mesh constraint or default ramp amplitude variations associated with the current step. # NODE Node number. # NNDOF Number of degrees of freedom at the node. # LNODETYPE Node type flag. LNODETYPE=1 indicates that the node is on the interior of the adaptive mesh region. LNODETYPE=2 indicates that the node is involved in a tied constraint. LNODETYPE=3 indicates that the node is at the corner of the boundary of an adaptive mesh region. LNODETYPE=4 indicates that the node lies on the edge of a boundary of an adaptive mesh region. LNODETYPE=5 indicates that the node lies on a flat surface on a boundary of the adaptive mesh region. LNODETYPE=6 indicates that the node participates in a constraint (other than a tied constraint) as a master node. LNODETYPE=7 indicates that the node participates in a constraint (other than a tied constraint) as a slave node. LNODETYPE=10 indicates that a concentrated load is applied to the node. # ALOCAL Local coordinate system aligned with the tangent to the adaptive mesh domain at the node. If the node is on the interior of the adaptive mesh domain, ALOCAL will be set to the identity matrix. In other cases the 1-direction is along an edge or in the plane of a flat surface. When NDIM=2, the 2-direction is normal to the surface. When NDIM=3, the 2-direction also lies in the plane of a flat surface or is arbitrary if the node is on an edge. When NDIM=3 the 3-direction is normal to the surface or is arbitrary if the node is on an edge. # NDIM Number of coordinate dimensions. # TIME(1) Current value of step time. # TIME(2) Current value of total time. # DTIME Time increment. # KSTEP Step number. # KINC Increment number. # KMESHSWEEP Mesh sweep number. # JMATYP Variable that must be passed into the GETVRMAVGATNODE utility routine to access local results at the node. # JGVBLOCK Variable that must be passed into the GETVRN, GETNODETOELEMCONN, and GETVRMAVGATNODE utility routines to access local results at the node.