FRIC

User subroutine to define frictional behavior for contact surfaces.

User subroutine FRIC:

can be used to define the frictional behavior between contacting surfaces;
can be used when the extended versions of the classical Coulomb friction model provided in Abaqus are too restrictive and a more complex definition of shear transmission between contacting surfaces is required;
will be called at points on the slave surface of a contact pair and at the integration points in a contact element (only when the contact point is closed) for which the contact interaction property model contains user-subroutine-defined friction;
must provide the entire definition of shear interaction between the contacting surfaces; and
can use and update solution-dependent state variables.

The following topics are discussed:

User subroutine interface
Variables to be defined
Variables that can be updated
Variables passed in for information

User subroutine interface

      SUBROUTINE FRIC(LM,TAU,DDTDDG,DDTDDP,DSLIP,SED,SFD,
     1 DDTDDT,PNEWDT,STATEV,DGAM,TAULM,PRESS,DPRESS,DDPDDH,SLIP,
     2 KSTEP,KINC,TIME,DTIME,NOEL,CINAME,SLNAME,MSNAME,NPT,NODE,
     3 NPATCH,COORDS,RCOORD,DROT,TEMP,PREDEF,NFDIR,MCRD,NPRED,
     4 NSTATV,CHRLNGTH,PROPS,NPROPS)
C
      INCLUDE 'ABA_PARAM.INC'
C
      CHARACTER*80 CINAME,SLNAME,MSNAME
C
      DIMENSION TAU(NFDIR),DDTDDG(NFDIR,NFDIR),DDTDDP(NFDIR),
     1 DSLIP(NFDIR),DDTDDT(NFDIR,2),STATEV(*),DGAM(NFDIR),
     2 TAULM(NFDIR),SLIP(NFDIR),TIME(2),COORDS(MCRD),
     3 RCOORD(MCRD),DROT(2,2),TEMP(2),PREDEF(2,*),PROPS(NPROPS)


      user coding to define LM, TAU, DDTDDG, DDTDDP,
      and, optionally, DSLIP, SED, SFD, DDTDDT, PNEWDT, STATEV


      RETURN
      END

Variables to be defined

In all cases

LM

Relative motion flag. User subroutine FRIC is called only if the contact point is determined to be closed; that is, if the contact pressure is positive (the contact point was closed in the previous iteration) or if the contact point is overclosed (the contact point was open in the previous iteration).

During iterations LM is passed into the subroutine as the value defined during the previous iteration. At the start of an increment or if the contact point opened during the previous iteration, this variable will be passed into the routine depending on the contact condition in the previous increment. If the contact point was slipping, LM is equal to 0; if the contact point was sticking, LM is equal to 1; and if the contact point was open, LM is equal to 2.

Set LM equal to 0 if relative motion is allowed (either due to slip or elastic stick). In this case the subroutine must specify the frictional stress $τ_{1}$ (and $τ_{2}$ for three-dimensional analysis) as a function of the relative sliding motion $γ_{1}$ (and $γ_{2}$ ), the interface contact pressure p, and other predefined or user-defined state variables. In addition, the subroutine must define the derivatives of the frictional stress with respect to $γ_{1}$ , ( $γ_{2}$ ), and p. For instance, in the case of isotropic elastic sticking, $\partial τ_{1} / \partial γ_{1} = \partial τ_{2} / \partial γ_{2} = k_{e l a s}$ , $\partial τ_{1} / \partial γ_{2} = \partial τ_{2} / \partial γ_{1} = 0$ , where $k_{e l a s}$ is the elastic stiffness of the interface.

Set LM equal to 1 if no relative motion is allowed; a rigid sticking condition at the interface is enforced by a Lagrange multiplier method. In this case no further variables need to be updated. If LM is always set to 1, a “perfectly rough” interface is created. It is not advisable to set LM to 1 when the finite-sliding, surface-to-surface contact formulation is used.

Set LM equal to 2 if friction is ignored (frictionless sliding is assumed). In this case no further variables need to be updated. If LM is always set to 2, a “perfectly smooth” interface is created.

You can make decisions about the stick/slip condition based on incremental slip information and calculated frictional stresses. These quantities are passed in by Abaqus/Standard, as discussed below.

To avoid convergence problems for the general class of frictional contact problems, set LM to 2 and exit this routine if the contact point was open at the end of the previous increment; that is, if Abaqus/Standard sets LM=2 when it calls this routine, simply exit the routine.

If the return value of LM is 0
TAU(NFDIR): These values are passed in as the values of the frictional stress components, $τ_{α}$ , at the beginning of the increment and must be updated to the values at the end of the increment. Here, and in the rest of this description, Greek subscripts ( $α$ , $β$ ) refer to frictional shear directions. The orientation of these directions on contact surfaces is defined in Contact formulations in Abaqus/Standard.
DDTDDG(NFDIR,NFDIR): $\partial Δ τ_{α}$ / $\partial Δ γ_{β}$ , partial derivative of the frictional stress in direction $α$ with respect to the relative motion in direction $β$ .
DDTDDP(NFDIR): $\partial Δ τ_{α}$ / $\partial Δ p$ , partial derivative of the frictional stress in direction $α$ with respect to the contact pressure. Since these terms yield an unsymmetric contribution to the stiffness matrix, they are used only if the unsymmetric equation solver is used (see Defining an analysis).

Variables that can be updated

DSLIP(NFDIR)

$Δ γ_{α}^{s l}$ , increment in nonrecoverable sliding motion (slip). If LM was 0 in the previous iteration, this array is passed in as the user-defined values during the previous iteration; otherwise, it will be zero. The array should be updated only if the return value of LM is 0.

This array is useful to detect slip reversals between iterations. It is used by the output options to indicate whether this point is sticking or slipping. Upon convergence of an increment, the values in DSLIP(NFDIR) are accumulated in SLIP(NFDIR), which are stored as the plastic strains.

SED

This variable is passed in as the value of the elastic energy density at the start of the increment and should be updated to the elastic energy density at the end of the increment. This variable is used for output only and has no effect on other solution variables.

SFD

This variable should be defined as the incremental frictional dissipation. The units are energy per unit area if the contact element or contact pair calling FRIC uses stresses as opposed to forces. For regular stress analysis this variable is used for output only and has no effect on other solution variables. In coupled temperature-displacement and coupled thermal-electrical-structural analyses the dissipation is converted into heat if the gap heat generation model is used. If SFD is not defined, the heat generation is calculated based on the dissipation obtained as the product of the slip increment, DSLIP, and the frictional stress, TAU.

DDTDDT(NFDIR,2)

$\partial Δ τ_{α}$ / $\partial Δ θ_{1}$ , $\partial Δ τ_{α}$ / $\partial Δ θ_{2}$ partial derivatives of the frictional stress in direction $α$ with respect to the temperatures of the two surfaces. This is required only for coupled temperature-displacement and coupled thermal-electrical-structural elements, in which the frictional stress is a function of the surface temperatures.

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 FRIC.

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 greater than 1.0 will be ignored and values of PNEWDT less than 1.0 will cause the job to terminate.

STATEV(NSTATV)

An array containing the user-defined solution-dependent state variables. You specify the number of available state variables; see Allocating space for details. This array will be passed in containing the values of these variables at the start of the increment. If any of the solution-dependent state variables is being used in conjunction with the friction behavior, they must be updated in this subroutine to their values at the end of the increment.

Variables passed in for information

DGAM(NFDIR): If LM was set to 0 in the previous iteration, this value is the increment of sliding motion in the current increment, $Δ γ_{α}$ . Otherwise, it will be zero. Comparison with DSLIP(NFDIR) makes it possible to determine whether slip changes to stick at this point and/or if there is a slip direction reversal occurring at this point.
TAULM(NFDIR): If LM was set to 1 in the previous iteration, this value is the current value of the constraint stress at the end of the increment, $τ_{α}^{L M}$ . Otherwise, it will be zero. Comparison with the critical shear stress makes it possible to determine whether stick changes to slip at this point.
PRESS: p, contact pressure at end of increment.
DPRESS: $Δ p$ , increment in contact pressure.
DDPDDH: $\partial Δ p$ / $\partial Δ h$ , current contact stiffness, in the case of soft contact (Contact pressure-overclosure relationships).
SLIP(NFDIR): Total nonrecoverable sliding motion (slip) at the beginning of the increment, $γ_{α}^{s l}$ . This value is the accumulated value of DSLIP(NFDIR) from previous increments.
KSTEP: Step number.
KINC: Increment number.
TIME(1): Value of step time at the end of the increment.
TIME(2): Value of total time at the end of the increment.
DTIME: Current increment in time.
NOEL: Element label for contact elements. Passed in as zero if contact surfaces are defined.
CINAME: User-specified surface interaction name associated with the friction definition, left justified. For contact elements it is the element set name given for the interface definition associated with the friction definition; if an optional name is assigned to the interface definition, CINAME is passed in as this name, left justified.
SLNAME: Slave surface name. Passed in as blank if contact elements are used.
MSNAME: Master surface name. Passed in as blank if contact elements are used.
NPT: Integration point number for contact elements. Passed in as zero if contact surfaces are defined.
NODE: User-defined global slave node number (or internal node number for models defined in terms of an assembly of part instances) involved with this contact point. Corresponds to the predominant slave node of the constraint if the surface-to-surface contact formulation is used. Passed in as zero if called from a contact element.
NPATCH: Not used.
COORDS(MCRD): An array containing the current coordinates of this point.
RCOORD(MCRD): If the master surface is defined as a rigid surface, this array is passed in containing the coordinates of the opposing point on the rigid surface in its current position and orientation.
DROT(2,2): Rotation increment matrix. For contact with a three-dimensional rigid surface, this matrix represents the incremental rotation of the surface directions relative to the rigid surface. It is provided so that vector- or tensor-valued state variables can be rotated appropriately in this subroutine. Stress and slip components are already rotated by this amount before FRIC is called. This matrix is passed in as a unit matrix for two-dimensional and axisymmetric contact problems.
TEMP(2): Current temperature at the slave node and the opposing master surface, respectively.
PREDEF(2,NPRED): An array containing pairs of values of all the user-specified field variables at the end of the current increment (initial values at the beginning of the analysis and current values during the analysis). If FRIC is called from a contact pair, the first value in a pair corresponds to the slave node and the second value corresponds to the nearest point on the master surface. If FRIC is called from a large-sliding contact element, PREDEF(1,NPRED) corresponds to the value at the integration point of the element and PFREDEF(2,NPRED) corresponds to the nearest point on the opposing surface. If FRIC is called from a small-sliding contact element, PREDEF(1,NPRED) corresponds to the value at the integration point of the first side and PFREDEF(2,NPRED) corresponds to the value at the integration point on the opposite face of the element.
NFDIR: Number of friction directions.
MCRD: Number of coordinate directions at the contact point.
NPRED: Number of predefined field variables.
NSTATV: Number of user-defined state variables.
CHRLNGTH: Characteristic contact surface face dimension, which can be used to define the maximum allowable elastic slip.
PROPS(NPROPS): Array of user-specified property values that are used to define the frictional behavior between the contacting surfaces.
NPROPS: User-specified number of property values associated with this friction model.