# Introduction and Basic Equations and $$ a _ {\alpha \beta} = \left(\mathbf {e} _ {\alpha}\right) ^ {T} \cdot \mathbf {a} \cdot \mathbf {e} _ {\beta} = \mathbf {e} _ {\alpha} \cdot \mathbf {a} \cdot \mathbf {e} _ {\beta}. $$ These component definitions are particularly convenient for calculating the dot product of two vectors, for we can write $\mathbf { a } \cdot \mathbf { b } = ( a _ { \alpha } \mathbf { e } ^ { \alpha } ) \cdot ( b ^ { \beta } \mathbf { e } _ { \beta } ) = a _ { \alpha } b ^ { \beta } \mathbf { e } ^ { \alpha } \cdot \mathbf { e } _ { \beta } ,$ $$ = a _ {\alpha} b ^ {\beta} \delta_ {\beta} ^ {\alpha}, $$ $$ = a _ {\alpha} b ^ {\alpha}, $$ which is $$ \mathbf {a} \cdot \mathbf {b} = a _ {1} b ^ {1} + a _ {2} b ^ {2} + a _ {3} b ^ {3}. $$ Similarly, the scalar product of two matrices is $$ \mathbf {a}: \mathbf {b} = a _ {\alpha \beta} b ^ {\beta \alpha}; $$ that $\mathbf { i s } ,$ we simply multiply corresponding entries in the $a _ { \alpha \beta }$ and $b ^ { \beta \alpha }$ arrays, arranged as matrices, and then sum the products. Finally, on the computer we need to store only one form of component: $a _ { \alpha } , a _ { \alpha \beta }$ or $a ^ { \alpha } , a ^ { \alpha \beta }$ . We can always go from one to the other using the "metric tensor," $g _ { \alpha \beta }$ , and its inverse, $g ^ { \alpha \beta }$ , which are defined as $$ g _ {\alpha \beta} = \mathbf {e} _ {\alpha} \cdot \mathbf {e} _ {\beta}, $$ and $$ g ^ {\alpha \beta} = \mathbf {e} ^ {\alpha} \cdot \mathbf {e} ^ {\beta}. $$ For $a ^ { \alpha } = \mathbf { a } \cdot \mathbf { e } ^ { \alpha }$ (from above), $= a _ { \beta } \mathbf { e } ^ { \beta } \cdot \mathbf { e } ^ { \alpha }$ (expressing a in component form) ; $= a _ { \beta } g ^ { \beta \alpha }$ (by the de¯nition of $g ^ { \beta \alpha } )$ Thus, $a ^ { \alpha } = g ^ { \beta \alpha } a _ { \beta }$ ; similarly $a _ { \alpha } = g _ { \beta \alpha } a ^ { \beta }$ , and, by extension, for matrices, $$ a ^ {\alpha \beta} = g ^ {\gamma \alpha} g ^ {\delta \beta} a _ {\gamma \delta} $$ and $$ a _ {\alpha \beta} = g _ {\gamma \alpha} g _ {\delta \beta} a ^ {\gamma \delta}. $$ The metric tensor and its inverse are symmetric: $$ g _ {\alpha \beta} = \mathbf {e} _ {\alpha} \cdot \mathbf {e} _ {\beta} = \mathbf {e} _ {\beta} \cdot \mathbf {e} _ {\alpha} = g _ {\beta \alpha}. $$ The two sets of base vectors and components of vectors or matrices associated with them are named as follows:

$\mathbf{e}_{\alpha}$	are covariant base vectors,
$\mathbf{e}^{\alpha}$	are contravariant base vectors,
$a_{\alpha}$ (or $a_{\alpha\beta}$ )	are covariant components of a vector (or matrix),
$a^{\alpha}$ (or $a^{\alpha\beta}$ )	are contravariant components of a vector (or matrix).

Thus, the contravariant components are those associated with the covariant base vectors, $\mathbf { a } = a ^ { \alpha } \mathbf { e } _ { \alpha }$ , and vice versa. The simplest case is when the basis is a set of orthogonal unit vectors (a rectangular Cartesian system) because then--from the definition $\mathbf { e } ^ { \alpha } \cdot \mathbf { e } _ { \beta } = \delta _ { \beta } ^ { \alpha }$ --we see that ${ \mathbf e } ^ { \alpha } = { \mathbf e } _ { \alpha }$ , and so $a _ { \alpha } = a ^ { \alpha }$ and we need not distinguish the type of component. Whenever possible a rectangular Cartesian system is chosen, so the type of component need not be distinguished. This system is discussed in more detail in the sections on beam elements and shell elements. # Components of a derivative Consider a vector-valued function, b, which is expressed in component form on a basis system, $\mathbf { e } _ { \alpha }$ . Let the vector-valued function a depend on b: a(b). Then $$ d \mathbf {a} = \left(\frac {\partial \mathbf {a}}{\partial \mathbf {b}}\right) \cdot \mathbf {e} _ {\alpha} d b ^ {\alpha} $$ so that the component of da associated with a change $d b ^ { \alpha }$ is $$ \left(\frac {\partial \mathbf {a}}{\partial \mathbf {b}}\right) \cdot \mathbf {e} _ {\alpha}, $$ which we write, for convenience, as $$ \left(\frac {\partial \mathbf {a}}{\partial b ^ {\alpha}}\right), $$ meaning $$ d \mathbf {a} = \left(\frac {\partial \mathbf {a}}{\partial b ^ {\alpha}}\right) d b ^ {\alpha}. $$ # Introduction and Basic Equations Now suppose a is written on a different basis-- $\mathbf { E } _ { \alpha }$ , say--so that we store a as the components $$ d a _ {\alpha} = d \mathbf {a} \cdot \mathbf {E} _ {\alpha}. $$ Then $$ d a _ {\alpha} = \mathbf {E} _ {\alpha} \cdot \frac {\partial \mathbf {a}}{\partial b ^ {\beta}} d b ^ {\beta}. $$ Typically we would then write $$ d a _ {\alpha} = H _ {\alpha \beta} d b ^ {\beta}, $$ where $$ H _ {\alpha \beta} = \mathbf {E} _ {\alpha} \cdot \frac {\partial \mathbf {a}}{\partial b ^ {\beta}} = \mathbf {E} _ {\alpha} \cdot \frac {\partial \mathbf {a}}{\partial \mathbf {b}} \cdot \mathbf {e} _ {\beta}. $$ Readers who are familiar with general curvilinear tensor analysis will recognize $H _ { \alpha \beta }$ as the covariant derivative of $a _ { \alpha }$ with respect to $b ^ { \beta }$ , often written as $a _ { \alpha \vert \beta }$ . The advantage of the direct matrix notation is clear: because we can imagine a and b as vectors in space, we have a physical understanding of what we mean by ${ \partial \mathbf { a } } / { \partial \mathbf { b } } ;$ ; it is the change in the vector-valued function a as a function of another vector-valued function b. For computations we must express a and b in component form. Then $$ H _ {\alpha \beta} = a _ {\alpha | \beta} = \mathbf {E} _ {\alpha} \cdot \frac {\partial \mathbf {a}}{\partial \mathbf {b}} \cdot \mathbf {e} _ {\beta} $$ provides the necessary components once we have chosen convenient basis systems: $\mathbf { e } _ { \alpha }$ for b and $\mathbf { E } _ { \alpha }$ for a. Typically ${ \bf e } _ { \alpha }$ and ${ \bf E } _ { \alpha }$ will both be the simple rectangular Cartesian bases $$ \mathbf {e} _ {1} = (1, 0, 0) $$ $$ \mathbf {e} _ {2} = (0, 1, 0) $$ $$ \mathbf {e} _ {3} = (0, 0, 1) $$ everywhere. But sometimes we must use more complicated basis systems--examples are when we need quantities associated with the surface of a general shell and when the symmetry of the geometry and, possibly, of the deformation makes it convenient to work in an axisymmetric system. The careful projection of the general results written in direct matrix notation onto the chosen basis system allows us to implement the theory for computation. As an example, consider the usual expression for strain rate, $$ \dot {\pmb {\varepsilon}} = \frac {1}{2} \left(\frac {\partial \dot {\mathbf {u}}}{\partial \mathbf {x}} + \left[ \frac {\partial \dot {\mathbf {u}}}{\partial \mathbf {x}} \right] ^ {T}\right), $$ # Introduction and Basic Equations which requires the matrix $\partial \dot { \bf u } / \partial { \bf x }$ to be evaluated, where u\_ is the velocity of the material currently flowing through the point x in space. Let us now derive the components of "\_ when the basis system for both u\_ and x is the cylindrical system that we usually choose for axisymmetric problems, with the basis vectors ${ \bf e } _ { 1 } \ ( \mathrm { r a d i a l } ) \ = ( \cos \theta , \sin \theta , 0 )$ ${ \bf e } _ { 2 } { \bf \Psi } ( \mathrm { a x i a l } ) { \bf \Psi } = ( 0 , 0 , 1 )$ $\mathbf { e } _ { 3 } ~ ( \mathrm { c i r c u m f e r e n t i a l } ) ~ = ( - \sin \theta , \cos \theta , 0 )$ (in ABAQUS for axisymmetric cases we always take the components in this order--radial, axial, circumferential). These basis vectors are orthogonal and of unit length, so that ${ \mathbf e } ^ { \alpha } = { \mathbf e } _ { \alpha }$ : We consider position to be defined by the coordinates $( r , z , \theta )$ , with $$ d \mathbf {x} = d r \mathbf {e} _ {1} + d z \mathbf {e} _ {2} + r d \theta \mathbf {e} _ {3}, $$ so that $$ d x ^ {1} = d r, \quad d x ^ {2} = d z, \quad \text { and } \quad d x ^ {3} = r d \theta . $$ Thus, $$ \frac {\partial \dot {\mathbf {u}}}{\partial x ^ {\alpha}} = \left(\frac {\partial \dot {\mathbf {u}}}{\partial r}, \frac {\partial \dot {\mathbf {u}}}{\partial z}, \frac {1}{r} \frac {\partial \dot {\mathbf {u}}}{\partial \theta}\right), $$ where $$ \dot {\mathbf {u}} = \dot {u} _ {r} \mathbf {e} _ {1} + \dot {u} _ {z} \mathbf {e} _ {2} + \dot {u} _ {\theta} \mathbf {e} _ {3}, $$ so that $$ \frac {\partial \dot {\mathbf {u}}}{\partial r} = \frac {\partial \dot {u} _ {r}}{\partial r} \mathbf {e} _ {1} + \frac {\partial \dot {u} _ {z}}{\partial r} \mathbf {e} _ {2} + \frac {\partial \dot {u} _ {\theta}}{\partial r} \mathbf {e} _ {3} $$ $$ \frac {\partial \dot {\mathbf {u}}}{\partial z} = \frac {\partial \dot {u} _ {r}}{\partial z} \mathbf {e} _ {1} + \frac {\partial \dot {u} _ {z}}{\partial z} \mathbf {e} _ {2} + \frac {\partial \dot {u} _ {\theta}}{\partial z} \mathbf {e} _ {3} $$ $$ \frac {1}{r} \frac {\partial \dot {\mathbf {u}}}{\partial \theta} = \frac {1}{r} \frac {\partial \dot {u} _ {r}}{\partial \theta} \mathbf {e} _ {1} + \frac {1}{r} \frac {\partial \dot {u} _ {z}}{\partial \theta} \mathbf {e} _ {2} + \frac {1}{r} \frac {\partial \dot {u} _ {\theta}}{\partial \theta} \mathbf {e} _ {3} + \frac {1}{r} \dot {u} _ {r} \frac {\partial \mathbf {e} _ {1}}{\partial \theta} + \frac {1}{r} \dot {u} _ {\theta} \frac {\partial \mathbf {e} _ {3}}{\partial \theta}. $$ We know that $$ \frac {\partial \mathbf {e} _ {1}}{\partial \theta} = (- \sin \theta , \cos \theta , 0) = \mathbf {e} _ {3}, \quad \text {and} \quad \frac {\partial \mathbf {e} _ {3}}{\partial \theta} = (- \cos \theta , - \sin \theta , 0) = - \mathbf {e} _ {1}, $$ so that $$ \frac {1}{r} \frac {\partial \dot {\mathbf {u}}}{\partial \theta} = \frac {1}{r} \left(\frac {\partial \dot {u} _ {r}}{\partial \theta} - \dot {u} _ {\theta}\right) \mathbf {e} _ {1} + \frac {1}{r} \frac {\partial \dot {u} _ {z}}{\partial \theta} \mathbf {e} _ {2} + \frac {1}{r} \left(\frac {\partial \dot {u} _ {\theta}}{\partial \theta} + \dot {u} _ {r}\right) \mathbf {e} _ {3}, $$ and thus, $$ \frac {\partial \dot {u} _ {\alpha}}{\partial x _ {\beta}} = \left[ \begin{array}{c c c} \partial \dot {u} _ {r} / \partial r & \partial \dot {u} _ {r} / \partial z & (1 / r) (\partial \dot {u} _ {r} / \partial \theta - \dot {u} _ {\theta}) \\ \partial \dot {u} _ {z} / \partial r & \partial \dot {u} _ {z} / \partial z & (1 / r) \partial \dot {u} _ {z} / \partial \theta \\ \partial \dot {u} _ {\theta} / \partial r & \partial \dot {u} _ {\theta} / \partial z & (1 / r) (\partial \dot {u} _ {\theta} / \partial \theta + \dot {u} _ {r}) \end{array} \right]. $$ The components of the strain rate are thus $$ \dot {\varepsilon} _ {r r} = \frac {\partial \dot {u} _ {r}}{\partial r}, \quad \dot {\varepsilon} _ {z z} = \frac {\partial \dot {u} _ {z}}{\partial z}, \quad \dot {\varepsilon} _ {\theta \theta} = \frac {1}{r} \left(\frac {\partial \dot {u} _ {\theta}}{\partial \theta} + \dot {u} _ {r}\right), $$ $$ \dot {\gamma} _ {r z} = 2 \dot {\varepsilon} _ {r z} = 2 \dot {\varepsilon} _ {z r} = \frac {\partial \dot {u} _ {r}}{\partial z} + \frac {\partial \dot {u} _ {z}}{\partial r}, \quad \dot {\gamma} _ {r \theta} = 2 \dot {\varepsilon} _ {r \theta} = 2 \dot {\varepsilon} _ {\theta r} = \frac {1}{r} \left(\frac {\partial \dot {u} _ {r}}{\partial \theta} - \dot {u} _ {\theta}\right) + \frac {\partial \dot {u} _ {\theta}}{\partial r}, $$ and $$ \dot {\gamma} _ {z \theta} = 2 \dot {\varepsilon} _ {z \theta} = 2 \dot {\varepsilon} _ {\theta z} = \frac {1}{r} \frac {\partial \dot {u} _ {z}}{\partial \theta} + \frac {\partial \dot {u} _ {\theta}}{\partial z}. $$ For the case of purely axisymmetric deformation, $\dot { u } _ { \theta } = 0$ and $\partial \dot { u } _ { r } / \partial \theta = \partial \dot { u } _ { z } / \partial \theta = 0$ , so these results simplify to the familiar expressions $$ \dot {\varepsilon} _ {r r} = \frac {\partial \dot {u} _ {r}}{\partial r}, \quad \dot {\varepsilon} _ {z z} = \frac {\partial \dot {u} _ {z}}{\partial z}, \quad \dot {\varepsilon} _ {\theta \theta} = \frac {\dot {u} _ {r}}{r} $$ $$ \dot {\gamma} _ {r z} = \frac {\partial \dot {u} _ {r}}{\partial z} + \frac {\partial \dot {u} _ {z}}{\partial r}, \quad \dot {\gamma} _ {r \theta} = \dot {\gamma} _ {z \theta} = 0. $$ In summary, direct matrix notation allows us to obtain all our fundamental results without reference to any particular choice of coordinate system. Careful application of the concept of the covariant derivative then allows these general results to be projected into component form for computation. # Virtual quantities The concepts of virtual displacements and virtual work are fundamental to the development. Virtual quantities are infinitesimally small variations of physical measures, such as displacement, strain, velocity, and so on. The virtual variation of a scalar quantity a is indicated by $\delta \boldsymbol { a } ;$ of a vector or matrix a by ±a. We extend this notation to such expressions as $$ \delta \pmb {\varepsilon} = \mathrm{sym} \left(\frac {\partial \delta \mathbf {v}}{\partial \mathbf {x}}\right), $$ which is the symmetric part of the spatial gradient of a virtual vector field ±v. This notation corresponds to the virtual rate of deformation (a measure of strain rate) if ±v is a virtual velocity field. # Initial and current positions Most structural problems concern the description of the way a structure behaves as it is loaded and moves from its reference configuration. Thus, we often compare positions of a point in the current (deformed) configuration and a reference configuration that is usually chosen as the configuration when the structure is unloaded or, in the case of geotechnical problems, when the model is subject only to geostatic stresses. To distinguish these configurations, we use lowercase type ( x) to indicate the current position and uppercase type (X) to indicate the initial position of the same material point in the same spatial coordinate frame. In ABAQUS we almost always store the rectangular Cartesian components of X and x. The exception is in axisymmetric structures, where radial (r) and axial (z) components are stored. # Nodal variables So far we have discussed quantities that are considered to be associated with all points in a model. The finite element approximation is based on assuming interpolations, by which displacement, position, and--often--other variables at any material point are defined by a finite number of nodal variables. In this manual we use uppercase superscripts to refer to individual nodal variables or nodal vectors and adopt the summation convention for these indices. Hence, the interpolation can be written quite generally as $$ \mathbf {a} = \mathbf {N} ^ {N} a ^ {N}, $$ where a is some vector-valued function at any point in the structure; ${ \bf N } ^ { N } ( g ^ { \alpha } ) , N = 1 , 2 . . .$ up to the total number of variables in the problem, is a set of N vector interpolation functions (these are functions of the material coordinates, $g ^ { \alpha } ) _ { ; }$ ; and $a ^ { N } , N = 1 , 2 . .$ : is a set of nodal variables. In some sections in this manual we need to describe operations on nodal variables for the complete system of finite element equations. In these sections we use the classical matrix-vector notation. In this notation $\{ a \}$ represents a column vector containing nodal variables, bac represents a row vector, and a matrix is written as [A]. Common operations are the scalar product between two vectors, $$ c = \left\lfloor a \right\rfloor \{b \} $$ (which is equivalent to $c = a ^ { N } b ^ { N }$ in index notation) and the matrix-vector product $$ \{c \} = [ A ] \{b \} $$ (which is equivalent to $c ^ { M } = A ^ { M N } b ^ { N }$ in index notation). # 1.3 Finite rotations # 1.3.1 Rotation variables Since ABAQUS contains such capabilities as structural elements (beams and shells) for which it is necessary to define arbitrarily large magnitudes of rotation, a convenient method for storing the rotation at a node is required. The components of a rotation vector $\phi$ are stored as the degrees of freedom 4, 5, and 6 at any node where a rotation is required. This method of storing rotation is described in \`\`Conventions,'' Section 1.2.2 of the ABAQUS/Standard User's Manual and the ABAQUS/Explicit User's Manual. The finite rotation vector $\phi$ consists of a rotation magnitude $\phi = \| \phi \|$ and a rotation axis or direction in space, $\mathbf { p } = \pmb { \phi } / | | \pmb { \phi } | |$ . Physically, the rotation $\phi$ is interpreted as a rotation by $\phi$ radians around the axis $\mathbf { p } .$ . To characterize this finite rotation mathematically, the rotation vector $\phi$ is used to define an orthogonal transformation or rotation matrix. To do so, first define the skew-symmetric matrix $\hat { \pmb { \phi } }$ associated with $\phi$ by the relationships $$ \hat {\pmb {\phi}} \cdot \pmb {\phi} = \mathbf {0} \quad \mathrm{and} \quad \hat {\pmb {\phi}} \cdot \mathbf {v} = \pmb {\phi} \times \mathbf {v} \quad \mathrm{forallvectors} \quad \mathbf {v}. $$ $\phi$ is called the axial vector of the skew-symmetric matrix $\hat { \pmb { \phi } } .$ In matrix components relative to the standard Euclidean basis, if ${ \pmb { \phi } } = \{ \phi ^ { 1 } \phi ^ { 2 } \phi ^ { 3 } \} ^ { T }$ , then $$ [ \hat {\pmb {\phi}} ] = \left[ \begin{array}{c c c} 0 & - \phi^ {3} & \phi^ {2} \\ \phi^ {3} & 0 & - \phi^ {1} \\ - \phi^ {2} & \phi^ {1} & 0 \end{array} \right]. $$ In what follows, a^ will be used to denote the skew-symmetric matrix with axial vector a. A well-known result from linear algebra is that the exponential of a skew-symmetric matrix $\hat { \pmb { \phi } }$ is an orthogonal (rotation) matrix that produces the finite rotation $\phi .$ Let the rotation matrix be $\mathbf { C } ,$ such that $\mathbf { C } ^ { - 1 } = \mathbf { C } ^ { T }$ . Then by definition, $$ \mathbf {C} = \exp [ \hat {\pmb {\phi}} ] = \mathbf {I} + \hat {\pmb {\phi}} + \frac {1}{2 !} \hat {\pmb {\phi}} ^ {2} + \dots . $$ However, the above infinite series has the following closed-form expression Equation 1.3.1-1 $$ \mathbf {C} = \exp [ \hat {\phi} ] = \cos \| \phi \| \mathbf {I} + \frac {\sin \| \phi \|}{\| \phi \|} \hat {\phi} + \frac {(1 - \cos \| \phi \|)}{\| \phi \| ^ {2}} \phi \phi . $$ In components, $$ C _ {i j} = \cos \phi \delta_ {i j} + (1 - \cos \phi) p _ {i} p _ {j} + \sin \phi \epsilon_ {i k j} p _ {k}, $$ where $\mathbf { p } = \{ p _ { 1 } \quad p _ { 2 } \quad p _ { 3 } \} ^ { T }$ and $\epsilon _ { i j k }$ is the alternator tensor, defined by $$ \epsilon_ {1 2 3} = \epsilon_ {2 3 1} = \epsilon_ {3 1 2} = 1; \quad \epsilon_ {1 3 2} = \epsilon_ {2 1 3} = \epsilon_ {3 2 1} = - 1; \text { all other } \epsilon_ {i j k} = 0. $$ It is this closed-form expression that allows the exact and numerically efficient geometric representation of finite rotations. # Quaternion parametrization Even though ABAQUS stores and outputs the rotation vector, quaternion parameters prove to be an efficient and convenient way to treat finite rotations computationally. Let $q _ { 0 } \in R$ be a scalar, and let $\mathbf { q } \in R ^ { 3 }$ be a vector field. The quaternion q is simply the pairing $$ \pmb {q} = (q _ {0}, \mathbf {q}). $$ To associate q with the finite rotation vector $\phi ,$ define the following: Equation 1.3.1-2 $$ q _ {0} = \cos \left\| \phi / 2 \right\| \quad \mathrm{and} \quad \mathbf {q} = \frac {\sin \left\| \phi / 2 \right\|}{\left\| \phi \right\|} \phi . $$ By trigonometric identities it follows that the orthogonal matrix $\exp [ \hat { \phi } ]$ in Equation 1.3.1-1 is given in terms of q as Equation 1.3.1-3 $$ \exp [ \hat {\pmb {\phi}} ] = (2 q _ {0} ^ {2} - 1) \mathbf {I} + 2 q _ {0} \hat {\mathbf {q}} + 2 \mathbf {q} \mathbf {q}. $$ By the convention introduced above, $\hat { \mathbf { q } }$ is the skew-symmetric matrix with axial vector q. For a more detailed discussion of quaternion algebra and its relation to other representations of finite rotations, see the discussion by Spring (1986). # Compound rotations A compound rotation is the successive application of two or more rotation fields. In geometrically linear problems compound rotations are obtained simply as the linear superposition of the individual (linearized) rotation vectors. This fact follows directly from the series expansion for $\exp [ \hat { \phi } ]$ . Let $\phi _ { 1 }$ and $\pmb { \phi } _ { 2 }$ be infinitesimal rotations. Thus, $\exp [ \hat { \pmb { \phi } } _ { 1 } ] \approx \mathbf { I } + \hat { \pmb { \phi } } _ { 1 } , \exp [ \hat { \pmb { \phi } } _ { 2 } ] \approx \mathbf { I } + \hat { \pmb { \phi } } _ { 2 }$ , and $$ \exp [ \hat {\pmb {\phi}} _ {1} ] \cdot \exp [ \hat {\pmb {\phi}} _ {2} ] \approx \exp [ \hat {\pmb {\phi}} _ {2} ] \cdot \exp [ \hat {\pmb {\phi}} _ {1} ] \approx \mathbf {I} + \hat {\pmb {\phi}} _ {1} + \hat {\pmb {\phi}} _ {2}. $$ In geometrically nonlinear analysis compound rotations are no longer additive. Furthermore, they are not commutative; that is, the order of application is important. A significant exception occurs when the multiple rotations share the same rotation axis. This special case is investigated further below. A detailed example of a finite compound rotation is given in \`\`Conventions,'' Section 1.2.2 of the ABAQUS/Standard User's Manual and the ABAQUS/Explicit User's Manual. Let $\mathbf { C } _ { i }$ be the orthogonal transformation representing the compound rotation defined as the product of a set of individual or incremental rotations $\Delta \mathbf { C } _ { p }$ , for $p = 1 , 2 , \ldots , i .$ . (For the case of specified boundary conditions $\mathbf { C } _ { i }$ is the final product after i steps of all the specified rotations $\Delta \mathbf { C } _ { p } ;$ for the iterative numerical solution procedure $\mathbf { C } _ { i }$ is the total rotation after i increments, where $\Delta \mathbf { C } _ { p }$ , for $p = 1 , \ldots , i ,$ is the converged rotation field solution at each increment.) By definition, the compound rotation is the product $$ \mathbf {C} _ {i} = \Delta \mathbf {C} _ {i} \cdot \Delta \mathbf {C} _ {i - 1} \dots \Delta \mathbf {C} _ {1}, $$ or equivalently by the recursion relation, $$ \mathbf {C} _ {i} = \Delta \mathbf {C} _ {i} \cdot \mathbf {C} _ {i - 1}. $$ It is important to note that $\Delta { \bf C } \cdot { \bf C }$ , which is interpreted as the finite rotation $\Delta \mathbf { C }$ superposed on the finite rotation C, is different from $\mathbf { C } \cdot \Delta \mathbf { C }$ , which is interpreted as the finite rotation C superposed on the finite rotation $\Delta \mathbf { C }$ . Although compound rotations are defined in terms of orthogonal matrices, in a numerical context the rotation vectors (or equivalently the quaternion parameters) associated with the rotation matrices are the degrees of freedom. Compound rotations are performed as follows: Given a quaternion parametrization ${ \pmb q } = ( q _ { 0 } , { \pmb q } )$ and an incremental (finite) rotation $\Delta \pmb q = \left( \Delta q _ { 0 } , \Delta \mathbf q \right)$ , where $\Delta \pmb q$ is defined in terms of an incremental rotation vector $\Delta \phi$ by Equation 1.3.1-2, the total or compound rotation is given by the quaternion ${ \pmb r } = ( r _ { 0 } , { \bf r } )$ , which is calculated as $$ \mathbf {r} = \Delta \mathbf {q} \circ \mathbf {q}. $$ Here ± denotes the quaternion product and is defined as Equation 1.3.1-4 $$ \Delta \mathbf {q} \circ \mathbf {q} \stackrel {\text {def}} {=} \left(\Delta q _ {0} q _ {0} - \Delta \mathbf {q} \cdot \mathbf {q}, \Delta q _ {0} \mathbf {q} + q _ {0} \Delta \mathbf {q} + \Delta \mathbf {q} \times \mathbf {q}\right). $$ Equation 1.3.1-4 allows for the update of rotation fields without ever calculating the orthogonal matrix from the quaternion and without performing a matrix multiplication. Furthermore, all operations are singularity free regardless of the magnitude of the incremental rotation field $\Delta \phi .$ . The final (total) rotation vector can be calculated from the quaternion r by inverting Equation 1.3.1-2. For the special case when compound rotations share the same rotation axis, the compound rotation reduces to an additive form. Let $\Delta \pmb q$ and q have the same rotation axis p. Then $$ \pmb {q} = (\cos \| \pmb {\phi} / 2 \|, \sin \| \pmb {\phi} / 2 \| \mathbf {p}), \Delta \pmb {q} = (\cos \| \Delta \pmb {\phi} / 2 \|, \sin \| \Delta \pmb {\phi} / 2 \| \mathbf {p}), \mathrm{and} $$ $$ \Delta \boldsymbol {q} \circ \boldsymbol {q} = (\cos \| \Delta \boldsymbol {\phi} / 2 \| \cos \| \boldsymbol {\phi} / 2 \| - \sin \| \Delta \boldsymbol {\phi} / 2 \| \sin \| \boldsymbol {\phi} / 2 \|, $$ $$ \cos \left\| \Delta \phi / 2 \right\| \sin \left\| \phi / 2 \right\| \mathbf {p} + \cos \left\| \phi / 2 \right\| \sin \left\| \Delta \phi / 2 \right\| \mathbf {p}), $$ which reduces to $$ \Delta \pmb {q} \circ \pmb {q} = (\cos \left\| (\pmb {\phi} + \Delta \pmb {\phi}) / 2 \right\|, \sin \left\| (\pmb {\phi} + \Delta \pmb {\phi}) / 2 \right\| \mathbf {p}). $$ # Rotation vector extraction For output purposes it is necessary to extract the rotation vector corresponding to a given quaternion. The extraction procedure is as follows: Let ${ \pmb r } = ( r _ { 0 } , { \bf r } )$ be the quaternion, and let $\phi$ be the rotation vector. Thus, Equation 1.3.1-5 $$ \| \boldsymbol {\phi} \| = 2 \tan^ {- 1} \left[ \frac {\| \mathbf {r} \|}{r _ {0}} \right] \quad \mathrm{and} \quad \boldsymbol {\phi} = \| \boldsymbol {\phi} \| \frac {\mathbf {r}}{\| \mathbf {r} \|}. $$ It is important to note that the extraction of the rotation vector from the quaternion is not unique. The magnitude $\| \phi \|$ is determined only up to the addition of 2n¼, $n = \ldots , - 1 , 0 , 1 , \ldots$ : ABAQUS will always choose that rotation vector such that $\| \phi \| < 2 \pi$ . # Director and rotation field updates As an example of the utility of the quaternion parameters, consider the incremental update of a director field for either a beam or shell analysis. At some stage of the solution the director field $\mathbf { t } ^ { i }$ , the quaternion parametrization of the rotation field ${ \pmb q } ^ { i }$ , and the incremental rotation field $\Delta \phi ^ { i }$ are known at increment i. To update the director field by the incremental rotation to increment $i + 1$ , proceed as follows: First calculate the quaternion parametrization of the incremental rotation: $$ \Delta q _ {0} ^ {i} = \cos \left\| \Delta \phi^ {i} / 2 \right\| \quad \mathrm{and} \quad \Delta \mathbf {q} ^ {i} = \frac {\sin \left\| \Delta \phi^ {i} / 2 \right\|}{\left\| \Delta \phi^ {i} \right\|} \Delta \phi^ {i}. $$ The director field at $i + 1$ is then defined as $\mathbf { t } ^ { i + 1 } = \exp [ \widehat { \Delta \phi } ^ { i } ] \cdot \mathbf { t } ^ { i }$ , where $\exp [ \widehat { \Delta \phi } ^ { i } ]$ is calculated with Equation 1.3.1-3. Thus, the director is calculated directly from the quaternion as $$ \mathbf {t} ^ {i + 1} = (2 (\Delta q _ {0} ^ {i}) ^ {2} - 1) \mathbf {t} ^ {i} + 2 \Delta q _ {0} ^ {i} \Delta \mathbf {q} ^ {i} \times \mathbf {t} ^ {i} + 2 (\Delta \mathbf {q} ^ {i} \cdot \mathbf {t} ^ {i}) \Delta \mathbf {q} ^ {i}. $$ Furthermore, the update of the rotation field is obtained by quaternion multiplication $\pmb q ^ { i + 1 } = \Delta \pmb q ^ { i } \circ \pmb q ^ { i }$ and is defined by $$ \boldsymbol {q} ^ {i + 1} = \left(q _ {0} ^ {i + 1}, \mathbf {q} ^ {i + 1}\right) = \left(\Delta q _ {0} ^ {i} q _ {0} ^ {i} - \Delta \mathbf {q} ^ {i} \cdot \mathbf {q} ^ {i}, \Delta q _ {0} ^ {i} \mathbf {q} ^ {i} + q _ {0} ^ {i} \Delta \mathbf {q} ^ {i} + \Delta \mathbf {q} ^ {i} \times \mathbf {q} ^ {i}\right). $$ # Variations of the rotation field In the development of the balance equations, it is necessary to calculate the variation of the rotation field. Consider the vector field a, which is obtained by rotation of the reference vector field A: