14 KiB
* tim .gt. tAccelerateMotor ) then
svars(2) = svars(2) + 1.0
end if
nrCrossings = int(svars(2))
c stop the motor if sensor crosses zero the second
c time
if (nrCrossings.eq.2) then
ampValueNew = 0.0
lFlagsDefine(iConcludeStep)=1
end if
c store sensor value
svars(1) = vTrans_CU1
end if
end if
return
end
1.2.10 VUANISOHYPER_INV: User subroutine to define anisotropic hyperelastic material behavior using the invariant formulation.
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_INV:
• can be used to define the strain energy potential of anisotropic hyperelastic materials as a function of an irreducible set of scalar invariants;
• will be called for blocks of material calculation points for which the material definition contains user-defined anisotropic hyperelastic behavior with invariant-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 scalar invariants.
Enumeration of invariants
To facilitate coding and provide easy access to the array of invariants passed to user subroutine VUANISOHYPER_INV, an enumerated representation of each invariant is introduced. Any scalar invariant can, therefore, be represented uniquely by an enumerated invariant, I _ { n } ^ { * } , where the subscript n denotes the order of the invariant according to the enumeration scheme in the following table:
| Invariant | Enumeration, n |
| $\overline{I}_{1}$ | 1 |
| $\overline{I}_{2}$ | 2 |
| J | 3 |
| $\overline{I}_{4(\alpha\beta)}$ | $4 + 2(\alpha - 1) + \beta(\beta - 1)$ ; $\alpha \leq \beta$ |
| $\overline{I}_{5(\alpha\beta)}$ | $5 + 2(\alpha - 1) + \beta(\beta - 1)$ ; $\alpha \leq \beta$ |
For example, in the case of three families of fibers there are a total of 15 invariants: \overline { { I } } _ { 1 } , \overline { { I } } _ { 2 } , J , six invariants of type \overline { { I } } _ { 4 ( \alpha \beta ) } , and six invariants of type \overline { { I } } _ { 5 ( \alpha \beta ) } , with \alpha , \beta = 1 , 2 , 3 \left( \alpha \leq \beta \right) . The following correspondence exists between each of these invariants and their enumerated counterpart:
| Enumerated invariant | Invariant |
| $I_1^*$ | $\overline{I}_1$ |
| $I_2^*$ | $\overline{I}_2$ |
| $I_3^*$ | $J$ |
| $I_4^*$ | $\overline{I}_{4(11)}$ |
| $I_5^*$ | $\overline{I}_{5(11)}$ |
| $I_6^*$ | $\overline{I}_{4(12)}$ |
| $I_7^*$ | $\overline{I}_{5(12)}$ |
| $I_8^*$ | $\overline{I}_{4(22)}$ |
| $I_9^*$ | $\overline{I}_{5(22)}$ |
| $I_{10}^*$ | $\overline{I}_{4(13)}$ |
| $I_{11}^*$ | $\overline{I}_{5(13)}$ |
| $I_{12}^*$ | $\overline{I}_{4(23)}$ |
| $I_{13}^*$ | $\overline{I}_{5(23)}$ |
| $I_{14}^*$ | $\overline{I}_{4(33)}$ |
| $I_{15}^*$ | $\overline{I}_{5(33)}$ |
A similar scheme is used for the array zeta of terms \zeta _ { \alpha \beta } = { \bf A } _ { \alpha } \cdot { \bf A } _ { \beta } . Each term can be represented uniquely by an enumerated counterpart \zeta _ { m } ^ { * } , as shown below:
| Dot product | Enumeration, m |
| $\zeta_{\alpha\beta}$ | $\alpha + \frac{1}{2}(\beta - 2)(\beta - 1)$ ; $\alpha < \beta$ |
As an example, for the case of three families of fibers there are three \zeta _ { \alpha \beta } terms: \zeta _ { 1 2 } , \zeta _ { 1 3 } , and \zeta _ { 2 3 } . These are stored in the zeta array as \left( \zeta _ { 1 } ^ { * } , \zeta _ { 2 } ^ { * } , \zeta _ { 3 } ^ { * } \right) .
Storage of arrays of derivatives of energy function
The components of the array duDi of first derivatives of the strain energy potential with respect to the scalar invariants, { \partial U } / { \partial I _ { i } ^ { * } } , are stored using the enumeration scheme discussed above for the scalar invariants.
The elements of the array d2uDiDi of second derivatives of the strain energy function, \partial ^ { 2 } U / \partial I _ { i } ^ { * } \partial I _ { j } ^ { * } , are laid out in memory using triangular storage: if denotes the component in this array corresponding to the term \partial ^ { 2 } U / \partial I _ { i } ^ { * } \partial I _ { i } ^ { * } , then k = i + j \times ( j - 1 ) / 2 ; ( i \leq j ) . For example, the term \partial ^ { 2 } U / \partial I _ { 2 } ^ { * } \partial I _ { 5 } ^ { * } is stored in component \overset { \cdot } { k } = 2 + ( 5 \times 4 ) / 2 = 1 2 in the d2uDiDi array.
Special considerations for shell elements
When VUANISOHYPER_INV 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_INV. 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_INV remains unchanged during the analysis; deleted material points are not removed from the block. Abaqus/Explicit will “freeze” the values of the invariants passed to VUANISOHYPER_INV for all deleted material points; that is, the 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_inv (
C Read only (unmodifiable) variables -
1 nblock, nFiber, nInv,
2 jElem, kIntPt, kLayer, kSecPt,
3 cmname,
4 nstatev, nfieldv, nprops,
5 props, tempOld, tempNew, fieldOld, fieldNew,
6 stateOld, sInvariant, zeta,
C Write only (modifiable) variables -
7 uDev, duDi, d2uDiDi,
8 stateNew)
C
include 'vaba_param.inc'
C
dimension props(nprops),
1 tempOld(nblock),
2 fieldOld(nblock,nfieldv),
3 stateOld(nblock,nstatev),
4 tempNew(nblock),
5 fieldNew(nblock,nfieldv),
6 sInvariant(nblock,nInv),
7 zeta(nblock,nFiber*(nFiber-1)/2),
8 uDev(nblock), duDi(nblock,nInv),
9 d2uDiDi(nblock,nInv*(nInv+1)/2),
* 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).
duDi(nblock,nInv)
Array of derivatives of strain energy potential with respect to the scalar invariants, { \partial U } / { \partial I _ { i } ^ { * } } , ordered using the enumeration scheme discussed above.
d2uDiDi(nblock,nInv*(nInv+1)/2)
Arrays of second derivatives of strain energy potential with respect to the scalar invariants (using triangular storage), \partial ^ { 2 } U / \partial I _ { i } ^ { * } \partial I _ { j } ^ { * } .
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.
nFiber
Number of families of fibers defined for this material.
nInv
Number of scalar invariants.
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
internal material models are given names starting with the “ABQ_” character string. To
you should not use “ABQ_” as the leading string for cmname.
nstatev
Number of user-defined state variables that are associated with this material type (as described in “Allocating space” in “User subroutines: overview,” Section 18.1.1
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.
sInvariant(nblock,nInv)
Array of scalar invariants, , at each material point at the end of the increment. The invariants are ordered using the enumeration scheme discussed above.
zeta(nblock,nFiber*(nFiber-1)/2) )
Array of dot product between the directions of different families of fiber in the reference configuration, . The array contains the enumerated values using the scheme discussed above.
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_INV, as illustrated below:
if (cmname(1:4) .eq. 'MAT1') then
call VUANISOHYPER_INV1(argument_list)
else if (cmname(1:4) .eq. 'MAT2') then
call VUANISOHYPER_INV2(argument_list)
end if
VUANISOHYPER_INV1 and VUANISOHYPER_INV2 are the actual subroutines containing the anisotropic hyperelastic models for each material MAT1 and MAT2, respectively. Subroutine VUANISOHYPER_INV merely acts as a directory here. The argument list can be the same as that used in subroutine VUANISOHYPER_INV. The material names must be in uppercase characters since cmname is passed in as an uppercase character string.
Example: Anisotropic hyperelastic model of Kaliske and Schmidt
As an example of the coding of subroutine VUANISOHYPER_INV, consider the model proposed by Kaliske and Schmidt (2005) for nonlinear anisotropic elasticity with two families of fibers. The strain energy function is given by a polynomial series expansion in the form
\begin{array}{l} U = \frac {1}{D} (J - 1) ^ {2} + \sum_ {i = 1} ^ {3} a _ {i} (\overline {{I}} _ {1} - 3) ^ {i} + \sum_ {j = 1} ^ {3} b _ {j} (\overline {{I}} _ {2} - 3) ^ {j} + \sum_ {k = 2} ^ {6} c _ {k} (\overline {{I}} _ {4 (1 1)} - 1) ^ {k} + \sum_ {l = 2} ^ {6} d _ {l} (\overline {{I}} _ {5 (1 1)} - 1) ^ {l} \\ + \sum_ {m = 2} ^ {6} e _ {m} (\overline {{I}} _ {4 (2 2)} - 1) ^ {m} + \sum_ {n = 2} ^ {6} f _ {n} (\overline {{I}} _ {5 (2 2)} - 1) ^ {n} + \sum_ {p = 2} ^ {6} g _ {p} (\zeta_ {1 2} \overline {{I}} _ {4 (1 2)} - \zeta_ {1 2} ^ {2}) ^ {p}. \\ \end{array}
The code in subroutine VUANISOHYPER_INV must return the derivatives of the strain energy function with respect to the scalar invariants, which are readily computed from the above expression. In this example auxiliary functions are used to facilitate enumeration of pseudo-invarinats of type \overline { { I } } _ { 4 ( \alpha \beta ) } and \overline { { I } } _ { 5 ( \alpha \beta ) } , as well as for indexing into the array of second derivatives using symmetric storage. The subroutine would be coded as follows:
subroutine vuanisohyper_inv (
C Read only -
* nblock, nFiber, nInv,
* jElem, kIntPt, kLayer, kSecPt,
* cmname,
* nstatev, nfieldv, nprops,
* props, tempOld, tempNew, fieldOld, fieldNew,
* stateOld, sInvariant, zeta,
C Write only -
* uDev, duDi, d2uDiDi,
* stateNew )
C
include 'vaba_param.inc'
C
dimension props(nprops),
* tempOld(nblock),
* fieldOld(nblock,nfieldv),
* stateOld(nblock,nstatev),
* tempNew(nblock),
* fieldNew(nblock,nfieldv),
* sInvariant(nblock,nInv),
* zeta(nblock,nFiber*(nFiber-1)/2),
* uDev(nblock), duDi(nblock,*),
* d2uDiDi(nblock,*),
* stateNew(nblock,nstatev)
C
character*80 cmname
C
parameter ( zero = 0.d0, one = 1.d0, two = 2.d0,
* three = 3.d0, four = 4.d0, five = 5.d0, six = 6.d0 )
C Kaliske energy function (3D)
C
C Read material properties
d=props(1)
dinv = one / d
a1=props(2)
a2=props(3)
a3=props(4)
b1=props(5)
b2=props(6)
b3=props(7)
c2=props(8)
c3=props(9)
c4=props(10)
c5=props(11)
c6=props(12)
d2=props(13)
d3=props(14)
d4=props(15)
d5=props(16)
d6=props(17)
e2=props(18)
e3=props(19)
e4=props(20)
e5=props(21)
e6=props(22)
f2=props(23)
f3=props(24)
f4=props(25)
f5=props(26)
f6=props(27)
g2=props(28)
g3=props(29)
g4=props(30)
g5=props(31)
g6=props(32)
C
do k = 1, nblock
Udev(k) = zero
C Compute Udev and 1st and 2nd derivatives w.r.t invariants
C - I1