SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD, 1 RPL,DDSDDT,DRPLDE,DRPLDT, 2 STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME, 3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT, 4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC) ! INCLUDE 'ABA_PARAM.INC' ! WARNING - the aba_param.inc file declares ! Implicit real*8(a-h,o-z) ! This means that, by default, any variables with ! first letter between a-h or o-z are double precision. ! The rest are integers. ! Note that this also means that if you type a variable ! name incorrectly, the compiler won't catch your typo. ! CHARACTER*80 CMNAME DIMENSION STRESS(NTENS),STATEV(NSTATV), 1 DDSDDE(NTENS,NTENS),DDSDDT(NTENS),DRPLDE(NTENS), 2 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1), 3 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3) ! ! DDSDDE(NTENS,NTENS) ! Jacobian matrix of the constitutive model. ! DDSDDE(I,J) defines the change in the Ith stress component ! at the end of the time increment caused by an infinitesimal ! perturbation of the Jth component of the strain increment array. ! Unless you invoke the unsymmetric equation solution capability ! for the user-defined material, ABAQUS/Standard will use only ! the symmetric part of DDSDDE. The symmetric part of the matrix ! is calculated by taking one half the sum of the matrix and its transpose. ! STRESS(NTENS) ! This array is passed in as the stress tensor at the beginning ! of the increment and must be updated in this routine to be the ! stress tensor at the end of the increment. If you specified ! initial stresses (“Initial conditions,” Section 19.2.1), this ! array will contain the initial stresses at the start of the ! analysis. The size of this array depends on the value of NTENS ! as defined below. In finite-strain problems the stress tensor ! has already been rotated to account for rigid body motion in ! the increment before UMAT is called, so that only the corotational ! part of the stress integration should be done in UMAT. The ! measure of stress used is “true” (Cauchy) stress. ! ! STATEV(NSTATV) ! An array containing the solution-dependent state variables. ! These are passed in as the values at the beginning of the ! increment unless they are updated in user subroutines USDFLD ! (“USDFLD,” Section 25.2.39) or UEXPAN (“UEXPAN,” Section 25.2.20), ! in which case the updated values are passed in. In all cases ! STATEV must be returned as the values at the end of the increment. ! The size of the array is defined as described in ! “Allocating space” in “User subroutines: overview,” Section 25.1.1. ! ! In finite-strain problems any vector-valued or tensor-valued ! state variables must be rotated to account for rigid body ! motion of the material, in addition to any update in the ! values associated with constitutive behavior. The rotation ! increment matrix, DROT, is provided for this purpose. ! ! SSE, SPD, SCD ! Specific elastic strain energy, plastic dissipation, and ! “creep” dissipation, respectively. These are passed in as ! the values at the start of the increment and should be ! updated to the corresponding specific energy values at ! the end of the increment. They have no effect on the solution, ! except that they are used for energy output. ! ! Only in a fully coupled thermal-stress analysis ! RPL ! Volumetric heat generation per unit time at the end of the increment ! caused by mechanical working of the material. ! ! DDSDDT(NTENS) ! Variation of the stress increments with respect to the temperature. ! ! DRPLDE(NTENS) ! Variation of RPL with respect to the strain increments. ! ! DRPLDT ! Variation of RPL with respect to the temperature. ! ! Variables that can be updated ! ! PNEWDT ! Ratio of suggested new time increment to the time increment being ! used (DTIME, see discussion later in this section). This variable ! allows you to provide input to the automatic time incrementation ! algorithms in ABAQUS/Standard (if automatic time incrementation is chosen). ! For a quasi-static procedure the automatic time stepping that ABAQUS/Standard ! uses, which is based on techniques for integrating standard creep laws ! (see “Quasi-static analysis,” Section 6.2.5), cannot be controlled from within ! the UMAT subroutine. ! PNEWDT is set to a large value before each call to UMAT. ! 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 ! ! STRAN(NTENS) ! An array containing the total strains at the beginning of the increment. ! If thermal expansion is included in the same material definition, the ! strains passed into UMAT are the mechanical strains only (that is, the ! thermal strains computed based upon the thermal expansion coefficient have ! been subtracted from the total strains). These strains are available for output ! as the “elastic” strains. ! ! In finite-strain problems the strain components have been rotated to account for ! rigid body motion in the increment before UMAT is called and are approximations ! to logarithmic strain. ! DSTRAN(NTENS) ! Array of strain increments. If thermal expansion is included in the same ! material definition, these are the mechanical strain increments (the total ! strain increments minus the thermal strain increments). ! ! 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. ! ! TEMP ! Temperature at the start of the increment. ! ! DTEMP ! Increment of temperature. ! ! 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. 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. ! ! NDI ! Number of direct stress components at this point. ! ! NSHR ! Number of engineering shear stress components at this point. ! ! NTENS ! Size of the stress or strain component array (NDI + NSHR). ! ! NSTATV ! Number of solution-dependent state variables that are associated with ! this material type (defined as described in “Allocating space” in “User ! subroutines: overview,” Section 25.1.1). ! ! 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 if geometric nonlinearity is accounted for during the step ! (see “Procedures: overview,” Section 6.1.1); otherwise, the array contains ! the original coordinates of the point. ! ! DROT(3,3) ! Rotation increment matrix. This matrix represents the increment of rigid ! body rotation of the basis system in which the components of stress ! (STRESS) and strain (STRAN) are stored. It is provided so that vector- or ! tensor-valued state variables can be rotated appropriately in this subroutine: ! stress and strain components are already rotated by this amount before UMAT ! is called. This matrix is passed in as a unit matrix for small-displacement ! analysis and for large-displacement analysis if the basis system for the ! material point rotates with the material (as in a shell element or when a ! local orientation is used). ! ! CELENT ! Characteristic element length, which is a typical length of a line across ! an element for a first-order element; it is half of the same typical length ! for a second-order element. For beams and trusses it is a characteristic length ! along the element axis. For membranes and shells it is a characteristic length ! in the reference surface. For axisymmetric elements it is a characteristic length ! in the plane only. For cohesive elements it is equal to the constitutive ! thickness. ! ! DFGRD0(3,3) ! Array containing the deformation gradient at the beginning of the increment. ! See the discussion regarding the availability of the deformation gradient for ! various element types. ! ! DFGRD1(3,3) ! Array containing the deformation gradient at the end of the increment. ! The components of this array are set to zero if nonlinear geometric effects ! are not included in the step definition associated with this increment. See ! the discussion regarding the availability of the deformation gradient for ! various element types. ! ! 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. ! user coding to define DDSDDE, STRESS, STATEV, SSE, SPD, SCD ! and, if necessary, RPL, DDSDDT, DRPLDE, DRPLDT, PNEWDT ! ! Local variables double precision E,xnu integer i,j ddsdde = 0.d0 E = props(1) xnu = props(2) ! for debugging, you can use ! write(6,*) ' Hello ' ! Output is then written to the .dat file If (ndi==3 .and. nshr==1) then ! Plane strain or axisymmetry ddsdde(1,1) = 1.d0-xnu ddsdde(1,2) = xnu ddsdde(1,3) = xnu ddsdde(2,1) = xnu ddsdde(2,2) = 1.d0-xnu ddsdde(2,3) = xnu ddsdde(3,1) = xnu ddsdde(3,2) = xnu ddsdde(3,3) = 1.d0-xnu ddsdde(4,4) = 0.5d0*(1.d0-2.d0*xnu) ddsdde = ddsdde*E/( (1.d0+xnu)*(1.d0-2.d0*xnu) ) else if (ndi==2 .and. nshr==1) then ! Plane stress ddsdde(1,1) = 1.d0 ddsdde(1,2) = xnu ddsdde(2,1) = xnu ddsdde(2,2) = 1.d0 ddsdde(3,3) = 0.5d0*(1.d0-xnu) ddsdde = ddsdde*E/( (1.d0+xnu*xnu) ) else ! 3D ! if (kinc==1) then ! write(*,*) 'hello' ! endif ddsdde(1,1) = 1.d0-xnu ddsdde(1,2) = xnu ddsdde(1,3) = xnu ddsdde(2,1) = xnu ddsdde(2,2) = 1.d0-xnu ddsdde(2,3) = xnu ddsdde(3,1) = xnu ddsdde(3,2) = xnu ddsdde(3,3) = 1.d0-xnu ddsdde(4,4) = 0.5d0*(1.d0-2.d0*xnu) ddsdde(5,5) = ddsdde(4,4) ddsdde(6,6) = ddsdde(4,4) ddsdde = ddsdde*E/( (1.d0+xnu)*(1.d0-2.d0*xnu) ) endif ! ! NOTE: ABAQUS uses engineering shear strains, ! i.e. stran(ndi+1) = 2*e_12, etc... do i = 1,ntens do j = 1,ntens stress(i) = stress(i) + ddsdde(i,j)*dstran(j) end do end do RETURN END