```txt bi1 = sInvariant(k,1) term = bi1-three Udev(k) = Udev(k) * + a1*term + a2*term**2 + a3*term**3 duDi(k,1) = a1 + two*a2*term + three*a3*term**2 d2uDiDi(k,indx(1,1)) = two*a2 + three*two*a3*term C - I2 bi2 = sInvariant(k,2) term = bi2-three Udev(k) = Udev(k) * + b1*term + b2*term**2 + b3*term**3 duDi(k,2) = b1 + two*b2*term + three*b3*term**2 d2uDiDi(k,indx(2,2)) = two*b2 + three*two*b3*term C - I3 (=J) bi3 = sInvariant(k,3) term = bi3-one duDi(k,3) = two*dinv*term d2uDiDi(k,indx(3,3)) = two*dinv C - I4(11) nI411 = indxInv4(1,1) bi411 = sInvariant(k,nI411) term = bi411-one Udev(k) = Udev(k) * + c2*term**2 + c3*term**3 + c4*term**4 * + c5*term**5 + c6*term**6 duDi(k,nI411) = * two*c2*term * + three*c3*term**2 * + four*c4*term**3 * + five*c5*term**4 * + six*c6*term**5 d2uDiDi(k,indx(nI411,nI411)) = * two*c2 * + three*two*c3*term * + four*three*c4*term**2 * + five*four*c5*term**3 * + six*five*c6*term**4 C - I5(11) nI511 = indxInv5(1,1) bi511 = sInvariant(k,nI511) term = bi511-one Udev(k) = Udev(k) ``` ```txt * + d2*term**2 + d3*term**3 + d4*term**4 * + d5*term**5 + d6*term**6 duDi(k,nI511) = * two*d2*term * + three*d3*term**2 * + four*d4*term**3 * + five*d5*term**4 * + six*d6*term**5 d2uDiDi(k,indx(nI511,nI511)) = * two*d2 * + three*two*d3*term * + four*three*d4*term**2 * + five*four*d5*term**3 * + six*five*d6*term**4 C - I4(22) nI422 = indxInv4(2,2) bi422 = sInvariant(k,nI422) term = bi422-one Udev(k) = Udev(k) * + e2*term**2 + e3*term**3 + e4*term**4 * + e5*term**5 + e6*term**6 duDi(k,nI422) = * two*e2*term * + three*e3*term**2 * + four*e4*term**3 * + five*e5*term**4 * + six*e6*term**5 d2uDiDi(k,indx(nI422,nI422)) = * two*e2 * + three*two*e3*term * + four*three*e4*term**2 * + five*four*e5*term**3 * + six*five*e6*term**4 C - I5(22) nI522 = indxInv5(2,2) bi522 = sInvariant(k,nI522) term = bi522-one Udev(k) = Udev(k) * + f2*term**2 + f3*term**3 + f4*term**4 * + f5*term**5 + f6*term**6 duDi(k,nI522) = * two*f2*term ``` ```txt * + three*f3*term**2 * + four*f4*term**3 * + five*f5*term**4 * + six*f6*term**5 d2uDiDi(k, index(nI522, nI522)) = * two*f2 * + three*two*f3*term * + four*three*f4*term**2 * + five*four*f5*term**3 * + six*five*f6*term**4 C - I4(12) nI412 = indexInv4(1, 2) bi412 = sInvariant(k, nI412) term = zeta(k, 1) * (bi412 - zeta(k, 1)) Udev(k) = Udev(k) * + g2*term**2 + g3*term**3 * + g4*term**4 + g5*term**5 * + g6*term**6 duDi(k, nI412) = zeta(k, 1) * ( * two*g2*term * + three*g3*term**2 * + four*g4*term**3 * + five*g5*term**4 * + six*g6*term**5 ) d2uDiDi(k, index(nI412, nI412)) = zeta(k, 1) ** 2 * ( * two*g2 * + three*two*g3*term * + four*three*g4*term**2 * + five*four*g5*term**3 * + six*five*g6*term**4 ) C end do C return end C Function to map index from Square to Triangular storage C of symmetric matrix C integer function index(i, j) include 'vaba_param.inc' ii = min(i, j) ``` ```matlab jj = max(i, j) index = ii + jj*(jj-1)/2 return end C C Function to generate enumeration of scalar C Pseudo-Invariants of type 4 C integer function indexInv4(i, j) include 'vaba_param.inc' ii = min(i, j) jj = max(i, j) indexInv4 = 4 + jj*(jj-1) + 2*(ii-1) return end C C Function to generate enumeration of scalar C Pseudo-Invariants of type 5 C integer function indexInv5(i, j) include 'vaba_param.inc' ii = min(i, j) jj = max(i, j) indexInv5 = 5 + jj*(jj-1) + 2*(ii-1) return end ``` # Additional reference • Kaliske, M., and J. Schmidt, “Formulation of Finite Nonlinear Anisotropic Elasticity,” CADFEM GmbH Infoplaner 2/2005, vol. 2, pp. 22–23, 2005. # 1.2.11 VUANISOHYPER\_STRAIN: User subroutine to define anisotropic hyperelastic material behavior based on Green strain. Product: Abaqus/Explicit # References • “Anisotropic hyperelastic behavior,” Section 22.5.3 of the Abaqus Analysis User’s Guide • \*ANISOTROPIC HYPERELASTIC • “UANISOHYPER\_INV and VUANISOHYPER\_INV,” Section 4.1.13 of the Abaqus Verification Guide # Overview User subroutine VUANISOHYPER\_STRAIN: • can be used to define the strain energy potential of anisotropic hyperelastic materials as a function of the components of the Green strain tensor; • will be called for blocks of material calculation points for which the material definition contains user-defined anisotropic hyperelastic behavior with Green strain-based formulation (“Anisotropic hyperelastic behavior,” Section 22.5.3 of the Abaqus Analysis User’s Guide); • can use and update solution-dependent state variables; • can use any field variables that are passed in; and • requires that the values of the derivatives of the strain energy density function be defined with respect to the components of the modified Green strain tensor and the volume ratio. # Component ordering in tensors The component ordering depends upon whether the tensor is second or fourth order. # Symmetric second-order tensors For symmetric second-order tensors, such as the modified Green strain tensor, there are ndir+nshr components; the component order is given as a natural permutation of the indices of the tensor. The direct components are first and then the indirect components, beginning with the 12-component. For example, a stress tensor contains ndir direct stress components and nshr shear stress components, which are passed in as
Component2D Case3D Case
1 $\overline{\varepsilon}_{11}^{G}$ $\overline{\varepsilon}_{11}^{G}$
2 $\overline{\varepsilon}_{22}^{G}$ $\overline{\varepsilon}_{22}^{G}$
Component2D Case3D Case
3 $\overline{\varepsilon}_{33}^{G}$ $\overline{\varepsilon}_{33}^{G}$
4 $\overline{\varepsilon}_{12}^{G}$ $\overline{\varepsilon}_{12}^{G}$
5 $\overline{\varepsilon}_{23}^{G}$
6 $\overline{\varepsilon}_{31}^{G}$
The shear strain components are stored as tensor components and not as engineering components. # Symmetric fourth-order tensors For symmetric fourth-order tensors, such as the deviatoric elasticity tensor $\partial ^ { 2 } U / \partial \overline { { \varepsilon } } _ { i j } ^ { G } \partial \overline { { \varepsilon } } _ { k l } ^ { G }$ , there are (ndir+nshr)\*(ndir+nshr+1)/2 independent components. These components are ordered using the following triangular storage scheme:
Component2D Case3D Case
1 $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{11}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{11}^{G}$
2 $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{22}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{22}^{G}$
3 $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{22}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{22}^{G}$
4 $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{33}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{33}^{G}$
5 $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{33}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{33}^{G}$
6 $\partial^{2}U/\partial\overline{\varepsilon}_{33}^{G}\partial\overline{\varepsilon}_{33}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{33}^{G}\partial\overline{\varepsilon}_{33}^{G}$
7 $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{12}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{12}^{G}$
8 $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{12}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{12}^{G}$
9 $\partial^{2}U/\partial\overline{\varepsilon}_{33}^{G}\partial\overline{\varepsilon}_{12}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{33}^{G}\partial\overline{\varepsilon}_{12}^{G}$
10 $\partial^{2}U/\partial\overline{\varepsilon}_{12}^{G}\partial\overline{\varepsilon}_{12}^{G}$ $\partial^{2}U/\partial\overline{\varepsilon}_{12}^{G}\partial\overline{\varepsilon}_{12}^{G}$
11 $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{23}^{G}$
12 $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{23}^{G}$
13 $\partial^{2}U/\partial\overline{\varepsilon}_{33}^{G}\partial\overline{\varepsilon}_{23}^{G}$
14 $\partial^{2}U/\partial\overline{\varepsilon}_{12}^{G}\partial\overline{\varepsilon}_{23}^{G}$
15 $\partial^{2}U/\partial\overline{\varepsilon}_{23}^{G}\partial\overline{\varepsilon}_{23}^{G}$
16 $\partial^{2}U/\partial\overline{\varepsilon}_{11}^{G}\partial\overline{\varepsilon}_{31}^{G}$
Component2D Case3D Case
17 $\partial^{2}U/\partial\overline{\varepsilon}_{22}^{G}\partial\overline{\varepsilon}_{31}^{G}$
18 $\partial^{2}U/\partial\overline{\varepsilon}_{33}^{G}\partial\overline{\varepsilon}_{31}^{G}$
19 $\partial^{2}U/\partial\overline{\varepsilon}_{12}^{G}\partial\overline{\varepsilon}_{31}^{G}$
20 $\partial^{2}U/\partial\overline{\varepsilon}_{23}^{G}\partial\overline{\varepsilon}_{31}^{G}$
21 $\partial^{2}U/\partial\overline{\varepsilon}_{31}^{G}\partial\overline{\varepsilon}_{31}^{G}$
If Q denotes the component number of term $\partial ^ { 2 } U / \partial \overline { { \varepsilon } } _ { i j } ^ { G } \partial \overline { { \varepsilon } } _ { k l } ^ { G }$ in the above table and M and N (with $M \ \leq \ N )$ denote the component numbers of $\overline { { \varepsilon } } _ { i j } ^ { G }$ and $\overline { { \varepsilon } } _ { k l } ^ { G }$ , respectively, in the table for second-order tensors, Q is given by the relationship $Q = M ^ { ' } + N \times ( N - 1 ) / 2$ . For example, consider the term $\partial ^ { 2 } U / \partial \overline { { \varepsilon } } _ { 1 1 } ^ { G } \partial \overline { { \varepsilon } } _ { 2 3 } ^ { G }$ . The component numbers for $\overline { { \varepsilon } } _ { 1 1 } ^ { G }$ and $\overline { { \varepsilon } } _ { 2 3 } ^ { G }$ are $M = 1$ and $N = 5$ , respectively, giving $Q = 1 + ( 5 \times 4 ) / 2 = 1 1$ . # Special consideration for shell elements When VUANISOHYPER\_STRAIN is used to define the material response of shell elements, Abaqus/Explicit cannot calculate a default value for the transverse shear stiffness of the element. Hence, you must define the element’s transverse shear stiffness. See “Shell section behavior,” Section 29.6.4 of the Abaqus Analysis User’s Guide, for guidelines on choosing this stiffness. # Material point deletion Material points that satisfy a user-defined failure criterion can be deleted from the model (see “User-defined mechanical material behavior,” Section 26.7.1 of the Abaqus Analysis User’s Guide). You must specify the state variable number controlling the element deletion flag when you allocate space for the solution-dependent state variables, as explained in “User-defined mechanical material behavior,” Section 26.7.1 of the Abaqus Analysis User’s Guide. The deletion state variable should be set to a value of one or zero in VUANISOHYPER\_STRAIN. A value of one indicates that the material point is active, and a value of zero indicates that Abaqus/Explicit should delete the material point from the model by setting the stresses to zero. The structure of the block of material points passed to user subroutine VUANISOHYPER\_STRAIN remains unchanged during the analysis; deleted material points are not removed from the block. Abaqus/Explicit will “freeze” the values of the strains passed to VUANISOHYPER\_STRAIN for all deleted material points; that is, the strain values remain constant after deletion is triggered. Once a material point has been flagged as deleted, it cannot be reactivated. # User subroutine interface subroutine vuanisohyper\_strain( C Read only (unmodifiable) variables – 1 nblock,jElem,kIntPt,kLayer,kSecPt,cmname, ```csv 2 ndir,nshr,nstatev,nfieldv,nprops, 3 props,tempOld,tempNew,fieldOld,fieldNew, 4 stateOld, ebar,detu, C Write only (modifiable) variables - 4 udev,duDe,duDj, 5 d2uDeDe,d2uDjDj,d2uDeDj, 6 stateNew) C include 'vaba_param.inc' C dimension jElem(nblock), 1 props(nprops), 2 tempOld(nblock), 3 fieldOld(nblock,nfieldv), 4 stateOld(nblock,nstatev), 5 tempNew(nblock), 6 fieldNew(nblock,nfieldv), 7 ebar(nblock,ndir+nshr), detu(nblock), 8 uDev(nblock), 9 duDe(nblock,ndir+nshr), duDj(nblock), * d2uDeDe(nblock,(ndir+nshr)*(ndir+nshr+1)/2), 1 d2uDjDj(nblock), 2 d2uDeDj(nblock,ndir+nshr), 3 stateNew(nblock,nstatev) C character*80 cmname C do 100 km = 1,nblock user coding 100 continue return end ``` # Variables to be defined udev(nblock) $\tilde { U } _ { d e v ; }$ , the deviatoric part of the strain energy density of the primary material response. This quantity is needed only if the current material definition also includes Mullins effect (see “Mullins effect,” Section 22.6.1 of the Abaqus Analysis User’s Guide). duDe(nblock,ndir+nshr) Derivatives of strain energy potential with respect to the components of the modified Green strain tensor, ${ \partial U } / { \partial \overline { { \varepsilon } } _ { i j } ^ { G } }$ . duDj(nblock,ndir+nshr) Derivatives of strain energy potential with respect to volume ratio, $\partial U / \partial J$ d2uDeDe(nblock,(ndir+nshr)\*(ndir+nshr+1)/2) Second derivatives of strain energy potential with respect to the components of the modified Green strain tensor (using triangular storage), $\partial ^ { 2 } U / \partial \overline { { \varepsilon } } _ { i j } ^ { G } \partial \overline { { \varepsilon } } _ { k l } ^ { G }$ . d2uDjDj(nblock) Second derivatives of strain energy potential with respect to volume ratio, $\partial ^ { 2 } U / \partial J ^ { 2 }$ . d2uDeDj(nblock,ndir+nshr) Cross derivatives of strain energy potential with respect to components of the modified Green strain tensor and volume ratio, $\partial ^ { 2 } U / \bar { \partial \varepsilon } _ { i j } ^ { G } \partial J ^ { 2 }$ . 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 “User subroutines: overview,” Section 18.1.1 of the Abaqus Analysis User’s Guide, for more information). # Variables passed in for information nblock Number of material points to be processed in this call to VUANISOHYPER\_STRAIN. 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 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. ndir Number of direct components in a symmetric tensor. nshr Number of indirect components in a symmetric tensor. nstatev Number of user-defined state variables that are associated with this material type (you define this as described in “Allocating space” in “User subroutines: overview,” Section 18.1.1 of the Abaqus Analysis User’s Guide). nfieldv Number of user-defined external field variables. nprops User-specified number of user-defined material properties. props(nprops) User-supplied material properties. tempOld(nblock) Temperatures at each material point at the beginning of the increment. tempNew(nblock) Temperatures at each material point at the end of the increment. fieldOld(nblock,nfieldv) Values of the user-defined field variables at each material point at the beginning of the increment. fieldNew(nblock,nfieldv) Values of the user-defined field variables at each material point at the end of the increment. stateOld(nblock,nstatev) State variables at each material point at the beginning of the increment. ebar(nblock,ndir+nshr) Modified Green strain tensor, $\overline { { \varepsilon } } ^ { G }$ , at each material point at the end of the increment. detu(nblock) J, determinant of deformation gradient (volume ratio) at the end of the increment. # Example: Using more than one user-defined anisotropic hyperelastic material model To use more than one user-defined anisotropic hyperelastic material model, the variable cmname can be tested for different material names inside user subroutine VUANISOHYPER\_STRAIN, as illustrated below: