ProductsAbaqus/Explicit Component ordering in tensorsThe component ordering depends upon whether the tensor is second or fourth order. Symmetric second-order tensorsFor 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
The shear strain components are stored as tensor components and not as engineering components. Symmetric fourth-order tensorsFor symmetric fourth-order tensors, such as the deviatoric elasticity tensor , there are (ndir+nshr)*(ndir+nshr+1)/2 independent components. These components are ordered using the following triangular storage scheme:
If Q denotes the component number of term in the above table and M and N (with ) denote the component numbers of and , respectively, in the table for second-order tensors, Q is given by the relationship . For example, consider the term . The component numbers for and are and , respectively, giving . Special consideration for shell elementsWhen 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 for guidelines on choosing this stiffness. Material point deletionMaterial points that satisfy a user-defined failure criterion can be deleted from the model (see User-defined mechanical material behavior). 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. 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,
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
Variables passed in for information
Example: Using more than one user-defined anisotropic hyperelastic material modelTo 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: if (cmname(1:4) .eq. 'MAT1') then call VUANISOHYPER_STRAIN1(argument_list) else if (cmname(1:4) .eq. 'MAT2') then call VUANISOHYPER_STRAIN2(argument_list) end if VUANISOHYPER_STRAIN1 and VUANISOHYPER_STRAIN2 are the actual subroutines containing the anisotropic hyperelastic models for each material MAT1 and MAT2, respectively. Subroutine VUANISOHYPER_STRAIN merely acts as a directory here. The argument list can be the same as that used in subroutine VUANISOHYPER_STRAIN. The material names must be in uppercase characters since cmname is passed in as an uppercase character string. Example: Orthotropic Saint-Venant Kirchhoff modelAs a simple example of the coding of subroutine VUANISOHYPER_STRAIN, consider the generalization to anisotropic hyperelasticity of the Saint-Venant Kirchhoff model. The strain energy function of the Saint-Venant Kirchhoff model can be expressed as a quadratic function of the Green strain tensor, , as where is the fourth-order elasticity tensor. The derivatives of the strain energy function with respect to the Green strain are given as However, subroutine VUANISOHYPER_STRAIN must return the derivatives of the strain energy function with respect to the modified Green strain tensor, , and the volume ratio, J, which can be accomplished easily using the following relationship between , , and : where is the second-order identity tensor. Thus, using the chain rule we find where and In this example an auxiliary function is used to facilitate indexing into a fourth-order symmetric tensor. The subroutine would be coded as follows: subroutine vuanisohyper_strain ( C Read only - * nblock, * jElem, kIntPt, kLayer, kSecPt, * cmname, * ndir, nshr, nstatev, nfieldv, nprops, * props, tempOld, tempNew, fieldOld, fieldNew, * stateOld, ebar, detu, C Write only - * uDev, duDe, duDj, * d2uDeDe, d2uDjDj, d2uDeDj, * stateNew ) C include 'vaba_param.inc' C dimension props(nprops), * tempOld(nblock), * fieldOld(nblock,nfieldv), * stateOld(nblock,nstatev), * tempNew(nblock), * fieldNew(nblock,nfieldv), * ebar(nblock,ndir+nshr), detu(nblock), * uDev(nblock), duDe(nblock,ndir+nshr), duDj(nblock), * d2uDeDe(nblock,*), d2uDjDj(nblock), * d2uDeDj(nblock,ndir+nshr), * stateNew(nblock,nstatev) C character*80 cmname C parameter( half = 0.5d0, one = 1.d0, two = 2.d0, * third = 1.d0/3.d0, twothds = 2.d0/3.d0, four = 4.d0, * dinv = 0.d0 ) C C Orthotropic Saint-Venant Kirchhoff strain energy function C (3D) C D1111 = props(1) D1122 = props(2) D2222 = props(3) D1133 = props(4) D2233 = props(5) D3333 = props(6) D1212 = props(7) D1313 = props(8) D2323 = props(9) C do k = 1, nblock C d2UdE11dE11 = D1111 d2UdE11dE22 = D1122 d2UdE11dE33 = D1133 d2UdE22dE11 = d2UdE11dE22 d2UdE22dE22 = D2222 d2UdE22dE33 = D2233 d2UdE33dE11 = d2UdE11dE33 d2UdE33dE22 = d2UdE22dE33 d2UdE33dE33 = D3333 d2UdE12dE12 = D1212 d2UdE13dE13 = D1313 d2UdE23dE23 = D2323 C xpow = exp ( log(detu(k)) * twothds ) detuInv = one / detu(k) C E11 = xpow * ebar(k,1) + half * ( xpow - one ) E22 = xpow * ebar(k,2) + half * ( xpow - one ) E33 = xpow * ebar(k,3) + half * ( xpow - one ) E12 = xpow * ebar(k,4) E23 = xpow * ebar(k,5) E13 = xpow * ebar(k,6) C term1 = twothds * detuInv dE11Dj = term1 * ( E11 + half ) dE22Dj = term1 * ( E22 + half ) dE33Dj = term1 * ( E33 + half ) dE12Dj = term1 * E12 dE23Dj = term1 * E23 dE13Dj = term1 * E13 term2 = - third * detuInv d2E11DjDj = term2 * dE11Dj d2E22DjDj = term2 * dE22Dj d2E33DjDj = term2 * dE33Dj d2E12DjDj = term2 * dE12Dj d2E23DjDj = term2 * dE23Dj d2E13DjDj = term2 * dE13Dj C dUdE11 = d2UdE11dE11 * E11 * + d2UdE11dE22 * E22 * + d2UdE11dE33 * E33 dUdE22 = d2UdE22dE11 * E11 * + d2UdE22dE22 * E22 * + d2UdE22dE33 * E33 dUdE33 = d2UdE33dE11 * E11 * + d2UdE33dE22 * E22 * + d2UdE33dE33 * E33 dUdE12 = two * d2UdE12dE12 * E12 dUdE23 = two * d2UdE23dE23 * E23 dUdE13 = two * d2UdE13dE13 * E13 C U = half * ( E11*dUdE11 + E22*dUdE22 + E33*dUdE33 ) * + E12*dUdE12 + E13*dUdE13 + E23*dUdE23 uDev(k) = U C duDe(k,1) = xpow * dUdE11 duDe(k,2) = xpow * dUdE22 duDe(k,3) = xpow * dUdE33 duDe(k,4) = xpow * dUdE12 duDe(k,5) = xpow * dUdE23 duDe(k,6) = xpow * dUdE13 C xpow2 = xpow * xpow C Only update nonzero components d2uDeDe(k,indx(1,1)) = xpow2 * d2UdE11dE11 d2uDeDe(k,indx(1,2)) = xpow2 * d2UdE11dE22 d2uDeDe(k,indx(2,2)) = xpow2 * d2UdE22dE22 d2uDeDe(k,indx(1,3)) = xpow2 * d2UdE11dE33 d2uDeDe(k,indx(2,3)) = xpow2 * d2UdE22dE33 d2uDeDe(k,indx(3,3)) = xpow2 * d2UdE33dE33 d2uDeDe(k,indx(4,4)) = xpow2 * d2UdE12dE12 d2uDeDe(k,indx(5,5)) = xpow2 * d2UdE23dE23 d2uDeDe(k,indx(6,6)) = xpow2 * d2UdE13dE13 C duDj(k) = dUdE11*dE11Dj + dUdE22*dE22Dj + dUdE33*dE33Dj * + two * ( dUdE12*dE12Dj + dUdE13*dE13Dj * + dUdE23*dE23Dj ) d2uDjDj(k)= dUdE11*d2E11DjDj+dUdE22*d2E22DjDj * +dUdE33*d2E33DjDj * + two*(dUdE12*d2E12DjDj+dUdE13*d2E13DjDj * +dUdE23*d2E23DjDj) * + d2UdE11dE11 * dE11Dj * dE11Dj * + d2UdE22dE22 * dE22Dj * dE22Dj * + d2UdE33dE33 * dE33Dj * dE33Dj * + two * ( d2UdE11dE22 * dE11Dj * dE22Dj * + d2UdE11dE33 * dE11Dj * dE33Dj * + d2UdE22dE33 * dE22Dj * dE33Dj ) * + four * ( d2UdE12dE12 * dE12Dj * dE12Dj * d2UdE13dE13 * dE13Dj * dE13Dj * d2UdE23dE23 * dE23Dj * dE23Dj ) C d2uDeDj(k,1) = xpow * ( term1 * dUdE11 * + d2UdE11dE11 * dE11Dj * + d2UdE11dE22 * dE22Dj * + d2UdE11dE33 * dE33Dj ) d2uDeDj(k,2) = xpow * ( term1 * dUdE22 * + d2UdE22dE11 * dE11Dj * + d2UdE22dE22 * dE22Dj * + d2UdE22dE33 * dE33Dj ) d2uDeDj(k,3) = xpow * ( term1 * dUdE33 * + d2UdE33dE11 * dE11Dj * + d2UdE33dE22 * dE22Dj * + d2UdE33dE33 * dE33Dj ) d2uDeDj(k,4) = xpow * ( term1 * dUdE12 * + two * d2UdE12dE12 * dE12Dj ) d2uDeDj(k,5) = xpow * ( term1 * dUdE23 * + two * d2UdE23dE23 * dE23Dj ) d2uDeDj(k,6) = xpow * ( term1 * dUdE13 * + two * d2UdE13dE13 * dE13Dj ) end do C return end C integer function indx( i, j ) C include 'vaba_param.inc' C C Function to map index from Square to Triangular storage C of symmetric matrix C ii = min(i,j) jj = max(i,j) C indx = ii + jj*(jj-1)/2 C return end |