Open Educational Resources

Linear Elastic Materials: Matrix of Material Properties of Linear Elastic Materials

Learning Outcomes

  • Identify that in general, after considering that the stress matrix, the strain matrix, and the coefficients matrix are symmetric, there are 21 constants to describe the relationship between the stress matrix and strain matrix.
  • Explain why the coefficients matrix needs to be symmetric.
  • Identify that for orthotropic materials, the coefficients matrix has 9 constants and the relationship is particular to a specific coordinate system. Explain the meaning of the 9 physical constants for orthotropic materials.
  • Identify that for transversely isotropic materials, the coefficients matrix has 5 constants and that the relationship is particular to a specific coordinate system. Explain the meaning of the 5 physical constants for transversely isotropic materials.
  • Identify that for isotropic materials, the coefficients matrix has 2 constants and that the relationship is valid in any coordinate system. Explain the meaning of the 2 physical constants for isotropic materials.
  • Given 6 stress and/or strain components, calculate the remaining components.
  • Describe the physical meaning of and differentiate between Young’s modulus, Poisson’s ratio, the shear modulus, and the bulk modulus for isotropic materials. Given any two of the  constants above, compute the remaining ones.

Introduction

An elastic material is defined as a material whose Cauchy stress tensor is a function of the current state of deformation independent of the history or path, I.e., \sigma:F\in\mathbb{M}^3_+\mapsto\sigma(F)\in\mathbb{M}^3. For a non-homogeneous material, this function would depend on the position X as well but this dependence is dropped at the moment for brevity. The principle of material frame-indifference restricts the possible forms for \sigma to:

    \[\sigma(F) = R\sigma(U)R^T\]

Equivalently, the constitutive law for the second Piola-Kirchhoff stress is restricted to have the form:

    \[S(F) = S(U)\]

The equivalence can be shown as follows, given S(F)=S(U), therefore:

    \[\sigma(U) = \frac{1}{det (U)} US(U)U\]

and

    \[\sigma(F) = \frac{1}{det (F)}FS(F)F^T\]

Therefore:

    \[S(U) = S(F) \Rightarrow \sigma(F) = R\sigma(U)R^T\]

The opposite statement can be shown similarly. In other words, for an elastic material, the second Piola-Kirchhoff stress tensor can be viewed as a function of the right stretch tensor U independent of the rotation R while the Cauchy stress matrix, as it is a spacial measure of stress depends on both U and R such that \sigma(RU)=R\sigma(U)R^T. For example, if a material is rotating in space, then, the second Piola-Kirchhoff stress tensor is expected to be constant, while the Cauchy stress matrix will change as it is tied to the spacial coordinates.

Small Deformations

Small deformations can be characterized by a deformed configuration function given by:

    \[x = X + u(X) = X +ev(X)\]

with \epsilon being a small parameter and v(X) is a vector field independent of the parameter \epsilon such that:

    \[\lim_{\epsilon \to 0} ||u(X)|| = 0\]

Following the notation adopted in Taylor’s Theorem, the expression \mathcal O (\epsilon^n) indicates an expression (scalar or tensor) that is directly proportional to \epsilon^n.

Measures of Deformation

In the case of small deformations, the deformation gradient can be written as:

    \[F = I + \nabla u = I + \mathcal O (\epsilon)\]

As a reminder, the displacement gradient can be decomposed into the small strain tensor \varepsilon and the infinitesimal rotation tensor W_{inf} resulting in:

    \[F = I + \nabla u = I + \varepsilon + W_{inf}\]

The right Cauchy-Green deformation tensor adopts the following form:

    \[C = F^TF = I + \nabla u + \nabla u^T + \mathcal O(\epsilon^2) = I +2\epsilon + \mathcal O(\epsilon^2)\]

The right stretch tensor U can be approximated as:

    \[U = C^{\frac{1}{2}} = I + \frac{1}{2}(\nabla u + \nabla u^T) + \mathcal O (\epsilon^2) = I + \varepsilon + \mathcal O (\epsilon^2)\]

With the small strain tensor \varepsilon appearing naturally in the approximation. The rotation R can be assumed to have the form:

    \[R = I + W_{inf} + \mathcal O (\epsilon^2)\]

such that:

    \[RU = I + \varepsilon + W_{inf} + \mathcal O (\epsilon^2)\]

It is important to note that the small deformations as shown by the above approximations imply both small strains and small rotations. In this case, the small strain \varepsilon and the Green Strain \varepsilon_{Green} matrices can be assumed equal as they differ by a small term:

    \[\varepsilon_{Green} = \frac{1}{2}(\nabla u + \nabla u^T) + \mathcal O (\epsilon^2) = \varepsilon_{small} + \mathcal O (\epsilon^2)\]

The volumetric strain can be evaluated as follows:

    \[\frac{V -V_0}{V_0} = J - 1 = det(I+\nabla u)-1 = 1 +Trace(\nabla u)-1+\mathcal O (\epsilon^2) = Trace(\varepsilon)+\mathcal O (\epsilon^2)\]

which is based on the following approximation (utilizing the relationship between the determinant and the trace as shown in the invariants section):

    \[J = det(F)\]

    \[= det(I+\nabla u)\]

    \[= \frac{1}{6}((I_1(I+\nabla u))^3-3I_1(I+\nabla u)^2I_1(I+\nabla u)+I_1(I+\nabla u)^3)\]

    \[= 1 + Trace(\nabla u) + \mathcal O (\epsilon^2)\]

Stress Tensors

Under the assumption of small deformations, the stress tensors are all equal by ignoring the terms \mathcal O (\epsilon):

    \[(1+ \mathcal O (\epsilon))\sigma = (I + \nabla u).S.(I+\nabla u^T) = S + \mathcal O (\epsilon)\]

Similarly

    \[(1+ \mathcal O (\epsilon))\sigma = P.(I + \nabla u) = P + \mathcal O (\epsilon)\]

Constitutive Relationship

As shown above, for elastic materials, the second Piola-Kirchhoff stress tensor is a function of U. Utilizing the assumption of small deformation and using a Taylor expansion around the identity matrix we get:

    \[S(U) = S(I+ \varepsilon + \mathcal O (\epsilon^2)) = S(I) + \frac{\partial S}{\partial U}\Big|_{U=I} \varepsilon + \mathcal O (\epsilon^2)\]

Assuming that the reference configuration is stress free, i.e., S(I)=0 and that the Cauchy and Piola-Kirchhoff stress tensors are equal:

    \[\sigma = \frac{\partial S}{\partial U}\Big|_{U=I} \varepsilon\]

where \frac{\partial S}{\partial U} is a matrix of 81 constants.

Matrix of Material Properties of Linear Elastic Materials

A linear elastic material is a material that exhibits a linear relationship between the components of the stress tensor and the components of the strain tensor. A linear elastic material constitutive law, under the assumption of small deformation, is fully represented by a linear map C_s between the stress matrix \sigma\in\mathbb{M}^3 and the infinitesimal strain matrix \varepsilon\in\mathbb{M}^3 such that C_s:\sigma\mapsto\varepsilon. In general, the stress matrix has 3\times 3=9 components and the strain matrix has 3\times 3=9 components. A linear relationship between these components would be represented by a fourth-order tensor with 3\times 3\times 3\times 3=81 components. These components are unique for a particular linear elastic material, and thus, they are considered to be material constants. By adopting the Einstein summation convention the linear map between the stress and the strain tensors can be written in component form as

    \[\sigma_{ij} = \sum_{k,l=1}^3 (C_{ijkl}\varepsilon_{kl}) = C_{ijkl}\varepsilon_{kl}\]

The components of the fourth-order tensor C_s can be visualized in the following form:

    \[C_s = \begin{pmatrix} \begin{pmatrix} C_{1111} & C_{1112} & C_{1113} \\ C_{1121} & C_{1122} & C_{1123} \\ C_{1131} & C_{1132} & C_{1133} \\ \end{pmatrix} & \begin{pmatrix} C_{1211} & C_{1212} & C_{1213} \\ C_{1221} & C_{1222} & C_{1223} \\ C_{1231} & C_{1232} & C_{1233} \\ \end{pmatrix} & \begin{pmatrix} C_{1311} & C_{1312} & C_{1313} \\ C_{1321} & C_{1322} & C_{1323} \\ C_{1331} & C_{1332} & C_{1333} \\ \end{pmatrix} \\ \begin{pmatrix} C_{2111} & C_{2112} & C_{2113} \\ C_{2121} & C_{2122} & C_{2123} \\ C_{2131} & C_{2132} & C_{2133} \\ \end{pmatrix} & \begin{pmatrix} C_{2211} & C_{2212} & C_{2213} \\ C_{2221} & C_{2222} & C_{2223} \\ C_{2231} & C_{2232} & C_{2233} \\ \end{pmatrix} & \begin{pmatrix} C_{2311} & C_{2312} & C_{2313} \\ C_{2321} & C_{2322} & C_{2323} \\ C_{2331} & C_{2332} & C_{2333} \\ \end{pmatrix} \\ \begin{pmatrix} C_{3111} & C_{3112} & C_{3113} \\ C_{3121} & C_{3122} & C_{3123} \\ C_{3131} & C_{3132} & C_{3133} \\ \end{pmatrix} & \begin{pmatrix} C_{3211} & C_{3212} & C_{3213} \\ C_{3221} & C_{3222} & C_{3223} \\ C_{3231} & C_{3232} & C_{3233} \\ \end{pmatrix} & \begin{pmatrix} C_{3311} & C_{3312} & C_{3313} \\ C_{3321} & C_{3322} & C_{3323} \\ C_{3331} & C_{3332} & C_{3333} \\ \end{pmatrix} \\ \end{pmatrix}\]

For example, the elements of the matrix on the top left of the above equation can be viewed as the matrix of coefficients that are used to calculate the stress component \sigma_{11} as follows:

    \[\sigma_{11} = C_{1111}\varepsilon_{11} + C_{1112}\varepsilon_{12} + C_{1113}\varepsilon_{13} + C_{1121}\varepsilon_{21} + C_{1122}\varepsilon_{22} + C_{1123}\varepsilon_{23} + C_{1131}\varepsilon_{31} + C_{1132}\varepsilon_{32} + C_{1133}\varepsilon_{33}\]

Reductions Due to Symmetries of the Stress and the Strain

Because the stress matrix is symmetric, i.e., \sigma_{ij}=\sigma_{ji}, then this leads to the following relationship among the components of C_s

    \[C_{ijkl} = C_{jikl}\]

This can be further explained using an example. Assume a strain state defined by a strain matrix with zero values except for the strain \varepsilon_{22}. Then, \sigma_{12}=C_{1222}\varepsilon_{22} and \sigma_{21}=C_{2122}\varepsilon_{22}. Since \sigma_{21}=\sigma_{12}, then, C_{1222}=C_{2122}. This applies to all the coefficients in the form of C_{ijkl} and C_{jikl}.

Similarly, because the strain matrix is symmetric, i.e., \varepsilon_{ij}=\varepsilon_{ji}, then this leads to the following relationship among the components of C_s

    \[C_{ijkl} = C_{ijlk}\]

This can be further explained using an example. Assume a strain state defined by a strain matrix with zero values except for the strain components \varepsilon_{12}=\varepsilon_{21}. Then, \sigma_{23}=C_{2312}\varepsilon_{12}+C_{2321}\varepsilon_{21}. Since \varepsilon_{12}=\varepsilon_{21} we have: \sigma_{23}=\frac{C_{2312}+C_{2321}}{2}\varepsilon_{12}+\frac{C_{2312}+C_{2321}}{2}\varepsilon_{21}, i.e., the coefficients C_{2312} and C_{2321} can be chosen to be equal. This applies to all the coefficients in the form of C_{ijkl} and C_{ijlk}.

Therefore, because of the symmetries of the stress and the strain matrices, the number of independent coefficients in the matrix C_s can be reduced to 6\times 6=36 independent coefficients. This simplifies the relationship which can now be represented in matrix form by assuming that the stress and the strain matrices are vectors in \mathbb{R}^6 with the following forms:

    \[\sigma = \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} \hspace{5mm} \varepsilon = \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} = \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{23} \\ \gamma_{12} \\ \gamma_{13} \\ \gamma_{23} \\ \end{pmatrix}\]

Notice that traditionally, the engineering shear strains \gamma_{12}, \gamma_{13}, and \gamma_{23} are used in the constitutive relationships and these are equal to double the corresponding shear strains. The linear map between the stress and the strain can be represented by a matrix C\in\mathbb{M}^6. The relationship between the vector representation of the stress and the vector representation of the strain has the following form:

    \[\sigma = C\varepsilon \hspace{5mm} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \begin{pmatrix} C_{1111} & C_{1122} & C_{1133} & C_{1112} & C_{1113} & C_{1123} \\ C_{2211} & C_{2222} & C_{2233} & C_{2212} & C_{2213} & C_{2223} \\ C_{3311} & C_{3322} & C_{3333} & C_{3312} & C_{3313} & C_{3323} \\ C_{1211} & C_{1222} & C_{1233} & C_{1212} & C_{1213} & C_{1223} \\ C_{1311} & C_{1322} & C_{1333} & C_{1312} & C_{1313} & C_{1323} \\ C_{2311} & C_{2322} & C_{2333} & C_{2312} & C_{2313} & C_{2323} \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix}\]

The above relationship is described by 36 independent coefficients. So far, there are no restrictions imposed on C for it to be symmetric. However, in the next section we will show that the restriction imposed by the preservation of energy leads to the symmetry of C.

Reductions Due to Preservation of Energy

For a linear elastic material, the energy stored during loading and/or the energy released during unloading have to be independent of the path of loading and/or unloading. Thus, for such a material, a potential energy function \overline{U} has to exist such that the stress is derived from this potential energy function as follows:

    \[\sigma_{ij} = \frac{\partial \overline{U}}{\partial \varepsilon_{ij}}\]

See chapter 6 in the book “Elasticity: Theory, Applications, and Numerics” by Martin H. Sadd for a detailed description of this concept.

The stress component \sigma_{ij} can be obtained using the following relationship:

    \[\sigma_{ij} = \sum_{k,l=1}^3 C_{ijkl}\varepsilon_{kl} \Rightarrow C_{ijkl} = \frac{\partial \sigma_{ij}}{\partial \varepsilon_{kl}} = \frac{\partial^2 \overline{U}}{\partial \varepsilon_{ij}\partial \varepsilon_{kl}}\]

The equality of the second partial derivatives of \overline{U} leads to the symmetry of C as follows:

    \[\frac{\partial^2 \overline{U}}{\partial \varepsilon_{ij} \partial \varepsilon_{kl}} = \frac{\partial^2 \overline{U}}{\partial \varepsilon_{kl} \partial \varepsilon_{ij}} \Rightarrow C_{ijkl} = C_{klij}\]

In matrix form, the relationship between the vector representation of the stress and the vector representation of the strain has the following form with 21 independent coefficients:

    \[\begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \begin{pmatrix} C_{1111} & C_{1122} & C_{1133} & C_{1112} & C_{1113} & C_{1123} \\ C_{2211} & C_{2222} & C_{2233} & C_{2212} & C_{2213} & C_{2223} \\ C_{3311} & C_{3322} & C_{3333} & C_{3312} & C_{3313} & C_{3323} \\ C_{1211} & C_{1222} & C_{1233} & C_{1212} & C_{1213} & C_{1223} \\ C_{1311} & C_{1322} & C_{1333} & C_{1312} & C_{1313} & C_{1323} \\ C_{2311} & C_{2322} & C_{2333} & C_{2312} & C_{2313} & C_{2323} \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix}\]

The symmetry of the coefficient matrix is crucial so that the material model would not predict energy creation or dissipation during loading cycles. This can be illustrated with an example in which the relationship between the in-plane strain and in-plane stress of a material is given by the following non-symmetric relationship:

    \[\begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \end{pmatrix} = \frac{1}{E} \begin{pmatrix} 1 & -0.3 \\ -0.4 & 1 \\ \end{pmatrix} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \end{pmatrix}\]

Assume that this material is loaded such that \sigma_{11} is increased gradually from 0 to 100MPa. Then, \sigma_{11} is kept constant and \sigma_{22} is increased gradually from 0 to 100MPa. Then, assume that both stresses are decreased simultaneously from 100MPa to 0. Since this is an elastic material, no energy should be created or dissipated through the entire loading process since the final state of stress is identical to the initial state of stress even though the paths of loading and unloading are different. The energy in a path from state a to state b can be calculated using the equation:

    \[\overline{U} = \int_a^b \sigma_{11} d\varepsilon_{11} + \int_a^b \sigma_{22} d\varepsilon_{22}\]

To simplify the integrals, we will use the following relationships:

    \[d\varepsilon_{11} = \frac{d \sigma_{11}}{E} -0.3\frac{d \sigma_{22}}{E}\]

    \[d\varepsilon_{22} = -0.4\frac{d \sigma_{11}}{E} + 1\frac{d \sigma_{22}}{E}\]

Therefore, the integrals have the following form:

    \[\overline{U} = \int_a^b \frac{\sigma_{11}}{E} d\sigma_{11} + \int_a^b -0.3\frac{\sigma_{11}}{E} d \sigma_{22} + \int_a^b -0.4\frac{\sigma_{22}}{E} d\sigma_{11} + \int_a^b \frac{\sigma_{22}}{E} d\sigma_{22}\]

The loading can be divided into three paths. In the first path, \sigma_{11} increases from 0 to 100MPa while \sigma_{22} stays constant at zero. Therefore, the energy can be calculated as follows:

    \[\overline{U}_{path 1} = \int_0^{100} \frac{\sigma_{11}}{E} d\sigma_{11} + 0 + 0 + 0 = \frac{100^2}{2E}\]

In the second path, \sigma_{11} stays constant at 100MPa while \sigma_{22} increases from 0 to 100MPa. Therefore, the energy in the second path can be calculated as follows:

    \[\overline{U}_{path 2} = 0 + \int_0^{100} -0.3\frac{\sigma_{11}}{E} d\sigma_{22} + 0 + \int_0^{100} \frac{\sigma_{22}}{E} d\sigma_{22} = (-0.3)\frac{100}{E}100 + \frac{100^2}{2E}\]

In the third path, \sigma_{11}=\sigma_{22}, \mathrm{d}\sigma_{11}=\mathrm{d}\sigma_{22} and both change from 100MPa to 0. Therefore, the energy in the third path can be calculated as follows:

    \[\overline{U}_{path 3} = \int_{100}^0 \frac{\sigma_{11}}{E} d\sigma_{11} + \int_{100}^0 -0.3\frac{\sigma_{11}}{E}d\sigma_{22} + \int_{100}^0 -0.4\frac{\sigma_{22}}{E} d\sigma_{11} + \int_{100}^0 \frac{\sigma_{22}}{E} d\sigma_{22}\]

    \[= -\frac{100^2}{2E} + 0.3\frac{100^2}{2E} + 0.4\frac{100^2}{2E} -\frac{100^2}{2E}\]

Therefore, the total energy after the three paths is equal to:

    \[\overline{U} = \overline{U}_{path 1} + \overline{U}_{path 2} + \overline{U}_{path 3} = (0.4-0.3)\frac{100^2}{2E}\]

which is contrary to the expectation that the material model preserves energy. The above loading cycle shows a dissipation of energy. The amount of energy dissipated is related to the difference between the off-diagonal components in the matrix of material coefficients. Had the matrix of coefficients been symmetric, the total energy stored or generated during the loading cycle would have been zero. Can you identify a different loading cycle in which the material would generate energy?

Reductions Due to Material Symmetry

The number of constants can be further reduced, depending on the symmetries of the material. If a material possesses a certain symmetry, this implies that the components of C_s do not change upon coordinate transformations within the specified symmetry. Since C_s is a fourth-order tensor, upon changing the coordinate system using an orthogonal matrix Q\in\mathbb{M}^3, the new material constants have the following form:

    \[C_{ijkl}' = Q_{im}Q_{jn}Q_{ko}Q_{lp}C_{mnop}\]

If a material’s behaviour is independent of a certain change of the coordinate system described by a particular orthogonal matrix Q, then the components of C'_s should be equal to those of C_s. Such orthogonal matrix Q is termed a material symmetry. For example, if a material response is symmetric around a plane, then a reflection matrix Q_1\in\mathbb{M}^3 across this plane would produce the same components for the coefficients matrix. If a material response is the same within any rotation inside a plane, then a rotation Q_2\in\mathbb{M}^3 within that plane would produce the same components for C'_s as those for C_s. The number of the independent constants of C_s can be further reduced from 21, depending on the level of anisotropy of the material. In the following section, the different levels of anisotropy are described, and the physical interpretation of some of the material constants is given.

Orthotropic Linear Elastic Materials

In an orthotropic linear elastic material, the 21 material constants are reduced to 9 material constants fully describing the relationship between the stress and the small strain. A proof of this argument is illustrated using the mathematica code found here. Assuming that the basis vectors e_1, e_2, and e_3 are the axis of orthotropy, then the relationship between the stress and the strain is given by:

(1)   \begin{equation*} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \begin{pmatrix} C_{1111} & C_{1122} & C_{1133} & 0 & 0 & 0 \\ C_{1122} & C_{2222} & C_{2233} & 0 & 0 & 0 \\ C_{1133} & C_{2233} & C_{3333} & 0 & 0 & 0 \\ 0 & 0 & 0 & C_{1212} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{1313} & 0 \\ 0 & 0 & 0 & 0 & 0 & C_{2323} \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} \end{equation*}

Notice that symmetry around each plane leads to the normal stresses being independent of the shear strains and vice versa. In other words, if an orthotropic material is pulled along a direction perpendicular to a plane of symmetry, then the associated shear strains are zero. For orthotropic materials, the relationship between the small strains and the stresses is traditionally given as function of nine physical material constants E_{11}, E_{22}, E_{33}, \nu_{12}, \nu_{13}, \nu_{23}, G_{12}, G_{13}, and G_{23} where: E_{ii}: The ratio of the uniaxial stress to uniaxial strain in the direction of the basis vector e_i. \nu_{ij}: Negative the ratio between the axial strain in the direction of the basis vector e_j to the axial strain in the direction of the basis vector e_i when the material is uniaxially stressed in the direction of the basis vector e_i. G_{ij}: The ratio between the shear stress \sigma_{ij} to the engineering shear strain \gamma_{ij}=2\varepsilon_{ij} in the plane defined by the directions of e_i and e_j.

The relationship can be written in terms of the above physical material constants as follows:

    \[\varepsilon_{11} = \frac{\sigma_{11}}{E_{11}} - \nu_{21}\frac{\sigma_{22}}{E_{22}} - \nu_{31}\frac{\sigma_{33}}{E_{33}}\]

    \[\varepsilon_{22} = -\nu_{12}\frac{\sigma_{11}}{E_{11}} + \frac{\sigma_{22}}{E_{22}} - \nu_{32}\frac{\sigma_{33}}{E_{33}}\]

    \[\varepsilon_{33} = -\nu_{13}\frac{\sigma_{11}}{E_{11}} - \nu_{23}\frac{\sigma_{22}}{E_{22}} + \frac{\sigma_{33}}{E_{33}}\]

    \[2\varepsilon_{12} = \gamma_{12} = \frac{\sigma_{12}}{G_{12}}\]

    \[2\varepsilon_{13} = \gamma_{13} = \frac{\sigma_{13}}{G_{13}}\]

    \[2\varepsilon_{23} = \gamma_{23} = \frac{\sigma_{23}}{G_{23}}\]

The relationship has the following matrix form:

(2)   \begin{equation*} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} = \begin{pmatrix} \frac{1}{E_{11}} & -\frac{\nu_{21}}{E_{22}} & -\frac{\nu_{31}}{E_{33}} & 0 & 0 & 0 \\ -\frac{\nu_{12}}{E_{11}} & \frac{1}{E_{22}} & -\frac{\nu_{32}}{E_{33}} & 0 & 0 & 0 \\ -\frac{\nu_{13}}{E_{11}} & -\frac{\nu_{23}}{E_{22}} & \frac{1}{E_{33}} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{G_{12}} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G_{13}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{G_{23}} \\ \end{pmatrix} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} \end{equation*}

Notice that symmetry of the matrix C implies that

    \[\frac{\nu_{ij}}{E_{ii}} = \frac{\nu_{ji}}{E_{jj}}\]

Transversely Isotropic Linear Elastic Materials

In a transversely isotropic linear elastic material, five material constants describe the relationship between the stress and the small strain. A proof of this argument is illustrated using the mathematica code found here. Assuming full isotropy in the plane formed by the basis vectors e_1 and e_2. i.e., the material response is identical given any rotation around e_3. Also, assuming that the plane formed by the vectors e_1 and e_2 is a plane of symmetry, then, the relationship between the stress and the strain is given by five material constants:

(3)   \begin{equation*} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \begin{pmatrix} C_{1111} & C_{1122} & C_{1133} & 0 & 0 & 0 \\ C_{1122} & C_{2222} & C_{2233} & 0 & 0 & 0 \\ C_{1133} & C_{1133} & C_{3333} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{C_{1111} - C_{1122}}{2} & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{1313} & 0 \\ 0 & 0 & 0 & 0 & 0 & C_{1313} \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} \end{equation*}

For transversely isotropic materials, the relationship between the small strains and the stresses are traditionally given as function of five physical material constants E, \nu, E_{33}, \nu_{13}, and G_{13} where:

E: The ratio of the uniaxial stress to uniaxial strain in any direction in the isotropy plane.
\nu: Negative the ratio between the transverse strain to the axial strain when the material is uniaxially stressed in any direction in the isotropy plane.
E_{33}: The ratio of the uniaxial stress to uniaxial strain in the direction perpendicular to the isotropy plane.
\nu_{13}: Negative the ratio between the axial strain \varepsilon_{33} to the axial strain \varepsilon_{11} when the material is stressed uniaxially in the direction of e_1
G_{13}: The ratio between the shear stress \sigma_{13} to the engineering shear strain \gamma_{13}=2\varepsilon_{13} in the plane defined by the directions of e_1 and e_3.

In this case, the relationship between the strains and the stresses can be written as:

    \[\varepsilon_{11} = \frac{\sigma_{11}}{E} - \nu\frac{\sigma_{22}}{E} - \nu_{31}\frac{\sigma_{33}}{E_{33}}\]

    \[\varepsilon_{22} = -\nu\frac{\sigma_{11}}{E} + \frac{\sigma_{22}}{E} - \nu_{31}\frac{\sigma_{33}}{E_{33}}\]

    \[\varepsilon_{33} = -\nu_{13}\frac{\sigma_{11}}{E} - \nu_{13}\frac{\sigma_{22}}{E} + \frac{\sigma_{33}}{E_{33}}\]

    \[2\varepsilon_{12} = \gamma_{12} = \frac{2(1+\nu)}{E}\sigma_{12}\]

    \[2\varepsilon_{13} = \gamma_{13} = \frac{\sigma_{13}}{G_{13}}\]

    \[2\varepsilon_{23} = \gamma_{23} = \frac{\sigma_{23}}{G_{13}}\]

The relationship has the following matrix form:

(4)   \begin{equation*} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} = \begin{pmatrix} \frac{1}{E} & -\frac{\nu}{E} & -\frac{\nu_{31}}{E_{33}} & 0 & 0 & 0 \\ -\frac{\nu}{E} & \frac{1}{E} & -\frac{\nu_{31}}{E_{33}} & 0 & 0 & 0 \\ -\frac{\nu_{13}}{E} & -\frac{\nu_{13}}{E} & \frac{1}{E_{33}} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{2(1+\nu)}{E} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G_{13}} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{G_{13}} \\ \end{pmatrix} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} \end{equation*}

The symmetry of the matrix C implies that

    \[\frac{\nu_{13}}{E} = \frac{\nu_{31}}{E_{33}}\]

Isotropic Linear Elastic Materials

In an isotropic linear elastic material, the response of the material to stress is independent of the coordinate system chosen. In other words, the relationship between the stress and the strain is independent of the coordinate system. This assumption reduces the number of material constants that describe the behaviour of the material in any direction to only two. A proof of this argument is illustrated using the mathematica code found here. The relationship between the stress and the strain components has the following form independent of the chosen coordinate system:

(5)   \begin{equation*} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \begin{pmatrix} C_{1111} & C_{1122} & C_{1122} & 0 & 0 & 0 \\ C_{1122} & C_{1111} & C_{1111} & 0 & 0 & 0 \\ C_{1122} & C_{1122} & C_{1111} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{C_{1111}-C_{1122}}{2} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{C_{1111}-C_{1122}}{2} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{C_{1111}-C_{1122}}{2} \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} \end{equation*}

The two constants that are generally used are the Young’s modulus E and Poisson’s ratio \nu, with the following physical definitions: Young’s modulus E: The slope of the stress-strain curve during uniaxial tensile tests of a sample of a linear elastic isotropic material. Poisson’s ratio \nu: Negative the ratio of the transverse strain (perpendicular to the applied load) to the axial strain (in the direction of the applied load) during uniaxial tensile tests of a sample of a linear elastic isotropic material.

For such materials, the strains are given by the following equations

    \[\varepsilon_{11} = \frac{\sigma_{11}}{E} - \nu\frac{\sigma_{22}}{E} - \nu\frac{\sigma_{33}}{E}\]

    \[\varepsilon_{22} = -\nu\frac{\sigma_{11}}{E} + \frac{\sigma_{22}}{E} - \nu\frac{\sigma_{33}}{E}\]

    \[\varepsilon_{33} = -\nu\frac{\sigma_{11}}{E} - \nu\frac{\sigma_{22}}{E} + \frac{\sigma_{33}}{E}\]

    \[2\varepsilon_{12} = \gamma_{12} = \frac{2(1+\nu)}{E}\sigma_{12}\]

    \[2\varepsilon_{13} = \gamma_{13} = \frac{2(1+\nu)}{E}\sigma_{13}\]

    \[2\varepsilon_{23} = \gamma_{23} = \frac{2(1+\nu)}{E}\sigma_{23}\]

The relationship has the following matrix form:

(6)   \begin{equation*} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} = \begin{pmatrix} \frac{1}{E} & -\frac{\nu}{E} & -\frac{\nu}{E} & 0 & 0 & 0 \\ -\frac{\nu}{E} & \frac{1}{E} & -\frac{\nu}{E} & 0 & 0 & 0 \\ -\frac{\nu}{E} & -\frac{\nu}{E} & \frac{1}{E} & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{G} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{G} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{G} \\ \end{pmatrix} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} \end{equation*}

Where:

    \[G = \frac{E}{2(1+\nu)}\]

The inverse of the above relationship has the following form:

(7)   \begin{equation*} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \frac{E}{(1-2\nu)(1+\nu)} \begin{pmatrix} (1-\nu) & \nu & \nu & 0 & 0 & 0 \\ \nu & (1-\nu) & \nu & 0 & 0 & 0 \\ \nu & \nu & (1-\nu) & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1-2\nu}{2} & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1-2\nu}{2} & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1-2\nu}{2} \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} \end{equation*}

Note that equations 1, 2, 3, and 4 are only valid in particular coordinate systems while equations 5, 6, and 7 can be used in any coordinate system.

Alternative Material Constants

The stress-strain relationships (Equations 6 and 7) rely on the definitions of E and \nu. However, other material constants can be used. The wikipedia page “Elastic Modulus” offers an extensive list of such constants. In particular, the shear modulus G, the bulk modulus K, and the Lamé’s constants \lambda and \mu are sometimes used. The following are their respective definitions:

Shear Modulus

The shear modulus is the slope of the shear stress versus the engineering shear strain of a linear elastic isotropic material. It is related to E and \nu as follows:

    \[G = \frac{E}{2(1+\nu)}\]

The above relationship can be deduced for an isotropic material as follows. We first use Equation 6 but assume that G, E, and \nu are independent of each other. Assume a state of pure shear stress with \sigma_{12}=\sigma_{21} being the only non-zero stress component. Using the isotropic linear elastic constitutive law equation (Equation 6), then the only non-zero strain components are 2\varepsilon_{12}=2\varepsilon_{21}=\frac{\sigma_{12}}{G}. Therefore, the stress and strain matrices are:

    \[\sigma = \begin{pmatrix} 0 & \sigma_{12} & 0 \\ \sigma_{12} & 0 & 0 \\ 0 & 0 & 0 \\ \end{pmatrix} \hspace{5mm} \varepsilon = \begin{pmatrix} 0 & \varepsilon_{12} & 0 \\ \varepsilon_{12} & 0 & 0 \\ 0 & 0 & 0 \\ \end{pmatrix}\]

We now apply a 45^{\circ} coordinate transformation in the plane described by e_1 and e_2 using the rotation matrix Q:

    \[Q = \begin{pmatrix} \frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} & 0 \\ -\frac{1}{\sqrt{2}} & \frac{1}{\sqrt{2}} & 0 \\ 0 & 0 & 1 \\ \end{pmatrix}\]

In the new coordinate system, \sigma' and \varepsilon' have the following forms:

    \[\sigma ' = Q\sigma Q^T = \begin{pmatrix} \sigma_{12} & 0 & 0 \\ 0 & -\sigma_{12} & 0 \\ 0 & 0 & 0 \\ \end{pmatrix} = \begin{pmatrix} \sigma_{11}' & 0 & 0 \\ 0 & \sigma_{22}' & 0 \\ 0 & 0 & 0 \\ \end{pmatrix}\]

    \[\varepsilon ' = Q\varepsilon Q^T = \begin{pmatrix} \varepsilon_{12} & 0 & 0 \\ 0 & -\varepsilon_{12} & 0 \\ 0 & 0 & 0 \\ \end{pmatrix} = \begin{pmatrix} \varepsilon_{11}' & 0 & 0 \\ 0 & \varepsilon_{22}' & 0 \\ 0 & 0 & 0 \\ \end{pmatrix}\]

Since the material is isotropic, the relationship shown in Equation 6 applies in the new coordinate system as well. Therefore:

    \[\varepsilon_{11}' = \frac{\sigma_{11}'}{E} - \nu\frac{\sigma_{22}'}{E} \Rightarrow \varepsilon_{12} = \frac{\sigma_{12}}{E}(1+\nu)\]

However, we already know from applying the constitutive equation in the original coordinate system that:

    \[\varepsilon_{12} =\frac{\sigma_{12}}{2G}\]

However, since \varepsilon'_{11}=\varepsilon_{12} we get the required relationship: G=\frac{E}{2(1+\nu)}.

Bulk Modulus

The bulk modulus is the slope of the applied hydrostatic stress versus the volumetric strain (\varepsilon_{11}+\varepsilon_{22}+\varepsilon_{33}) of a linear elastic isotropic material. It is related to E and \nu as follows:

    \[K = \frac{E}{3(1-2\nu)}\]

The relationship for the bulk modulus can be obtained as follows. First we assume a hydrostatic state of stress \sigma=pI where I is the identity matrix, and p\in\mathbb{R} is a real number representing the applied pressure. In this case, when adopting Equation 6, the developed strains are \varepsilon_{11}=\varepsilon_{22}=\varepsilon_{33}=p\frac{1-2\nu}{E}. Therefore,

    \[\varepsilon_{11} + \varepsilon_{22} + \varepsilon_{33} = p\frac{3(1-2\nu)}{E}\]

Therefore, the slope of the relationship between p and the volumetric strain \varepsilon_{11}+\varepsilon_{22}+\varepsilon_{33} is equal to \frac{E}{3(1-2\nu)}.

LAMÉ’S CONSTANTS:

These constants are defined as follows:

    \[\lambda = \frac{E\nu}{(1+\nu)(1-2\nu)} \hspace{5mm} \mu = G = \frac{E}{2(1+\nu)}\]

These two constants provide a simplified way of writing 7 as follows.

(8)   \begin{equation*} \begin{pmatrix} \sigma_{11} \\ \sigma_{22} \\ \sigma_{33} \\ \sigma_{12} \\ \sigma_{13} \\ \sigma_{23} \\ \end{pmatrix} = \begin{pmatrix} \lambda + 2\mu & \lambda & \lambda & 0 & 0 & 0 \\ \lambda & \lambda + 2\mu & \lambda & 0 & 0 & 0 \\ \lambda & \lambda & \lambda + 2\mu & 0 & 0 & 0 \\ 0 & 0 & 0 & \mu & 0 & 0 \\ 0 & 0 & 0 & 0 & \mu & 0 \\ 0 & 0 & 0 & 0 & 0 & \mu \\ \end{pmatrix} \begin{pmatrix} \varepsilon_{11} \\ \varepsilon_{22} \\ \varepsilon_{33} \\ 2\varepsilon_{12} \\ 2\varepsilon_{13} \\ 2\varepsilon_{23} \\ \end{pmatrix} \end{equation*}

or in simple component form

    \[\sigma_{ij} = 2\mu \varepsilon_{ij} + \sum_{k=1}^3 \lambda \varepsilon_{kk} \delta_{ij}\]

Videos

Derivation of the Reduced Number of Independent Matrix Coefficients due to Material Symmetries

The following Mathematica code is used to reduce the number of independent coefficients in the fourth-order tensor. In the first part, the symmetry of the stress tensor and the symmetry of the strain tensor are utilized to reduce the number of independent coefficients from 81 to 36. Then, the preservation of energy concept is utilized to reduce the number of independent coefficients from 36 to 21. Next, reflections along the directions of e_1, e_2, and e_3 are applied to find the number of coefficients required for an orthotropic material, with the planes of symmetry being the planes perpendicular to the three basis vectors. In this case, the number of independent constants is shown to be 9. Then, the number of independent constants in a transversely isotropic material is found by utilizing the orthotropic material fourth-order tensor and assuming that the behaviour of the material is independent of a 45^\circ rotation in the plane of e_1 and e_2. The resulting number of independent material constants is shown to be 5. Finally, another 45^\circ rotation around e_1 is used to show that the number of independent material constants for an isotropic material is just 2. Copy and paste the code into your Mathematica browser, and follow it step by step. You should also show that the results for the transversely isotropic and isotropic materials are independent of the angle of rotation used in the derivation, as long as the angle is not a multiple of \frac{\pi}{2}.

View Mathematica Code

(*General Ctensor (fourth order tensor)*)
Ctensor=Table[Row[{c,i,j,k,l}],{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["A general C matrix relating 9 stress tensor components to 9 strain tensor components:"] 
Ctensor//MatrixForm
(*Applying stress symmetry*)
RuleStressSymmetry={};
Do[AppendTo[RuleStressSymmetry,Ctensor[[i,j,k,l]]->Ctensor[[j,i,k,l]]],{i,2,3},{j,1,i-1},{k,1,3},{l,1,3}];
Print["Stress tensor symmetry leads to the following relations:",RuleStressSymmetry] 
CtensorSt=Ctensor/.RuleStressSymmetry;
Print["The resulting C matrix from stress tensor symmetry is:"] 
CtensorSt//MatrixForm
(*Applying Strain symmetry*)
RuleStrainSymmetry={};
Do[AppendTo[RuleStrainSymmetry,CtensorSt[[i,j,k,l]]->CtensorSt[[i,j,l,k]]],{i,1,3},{j,1,3},{k,2,3},{l,1,k-1}];
Print["Strain tensor symmetry leads to the following relations:",RuleStrainSymmetry] 
CtensorStStr=CtensorSt/.RuleStrainSymmetry;
Print["The resulting C matrix from both stress and strain tensors symmetry is:"] 
CtensorStStr//MatrixForm
(*Using symmetry due to Energy Preservation*)
RuleEnergySymmetry={};
Do[AppendTo[RuleEnergySymmetry,CtensorStStr[[i,j,k,l]]->CtensorStStr[[k,l,i,j]]],{i,1,3},{j,i,3},{k,1,i},{l,k,If[k<i,3,j]}];
Print["Symmetry due to the existence of an energy function leads to the following relations:",RuleEnergySymmetry] 
CtensorStStrEn=CtensorStStr/.RuleEnergySymmetry;
Print["The resulting C matrix after using stress, strain and energy symmetries is:"] 
CtensorStStrEn//MatrixForm
Print["The original C tensor and the fully symmetric Ctensor:"] 
{Ctensor//MatrixForm,CtensorStStrEn//MatrixForm} 
Print["The difference (minus) between C tensor and the fully symmetric Ctensor:"] 
Ctensor-CtensorStStrEn//MatrixForm
(*Adopting the new Ctensor after symmetries*)
Ctensor=CtensorStStrEn;
(*Order of the stress vector*)
order={{1,1},{2,2},{3,3},{1,2},{1,3},{2,3}};
(*Converting Ctensor to Cmatrix*)
Cmatrix=Table[ToExpression[ToString[Row[{"Ctensor","[[",Row[Table[If[OddQ[i],Flatten[{order[[j]],order[[k]]}][[(i+1)/2]],","],{i,1,7}]],"]]"}]]],{j,1,6},{k,1,6}];
Print["The Following is the resulting C matrix after adopting all the symmetries:"] 
Cmatrix//MatrixForm
(*Orthotropic material*)
Print["Deriving C tensor and C matrix for Orthotropic Materials"] 
(*The following are three reflections around the three planes of symmetry perpendicular to e1=ex and e2=ey and e3=ez*)
Qx=ReflectionMatrix[{1,0,0}];
Qy=ReflectionMatrix[{0,1,0}];
Qz=ReflectionMatrix[{0,0,1}];
(*Using plane perpendicular to x as a symmetry plane*)
Cdashtensor=Table[Sum[Qx[[i,m]]*Qx[[j,n]]*Qx[[k,o]]*Qx[[l,p]]*Ctensor[[m,n,o,p]],{m,1,3},{n,1,3},{o,1,3},{p,1,3}],{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["The following is the C tensor after reflection around the plane perpendicular to e1=ex:"] 
Cdashtensor//MatrixForm
Print["Following is the difference between the C tensor after reflection around the plane perpendicular to e1=ex and the original symmetric C tensor:"] 
(Ctensor-Cdashtensor)//MatrixForm
zero=Table[0,{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["Following are the relationships after equating C tensor after reflection to the Ctensor before reflection around the plane perpendicular to e1=ex:"] 
ax=Solve[(Ctensor-Cdashtensor)==zero,Flatten[Ctensor]] 
Ctensorx=Ctensor/.ax[[1]];
Print["Following is the resulting C tensor after adopting the previous relations:"] 
Ctensorx//MatrixForm
(*Using plane perpendicular to y as a symmetry plane*)
Cdashtensor=Table[Sum[Qy[[i,m]]*Qy[[j,n]]*Qy[[k,o]]*Qy[[l,p]]*Ctensor[[m,n,o,p]],{m,1,3},{n,1,3},{o,1,3},{p,1,3}],{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["Following is the C tensor after reflection around the plane perpendicular to e2=ey:"] 
Cdashtensor//MatrixForm
Print["Following is the difference between the C tensor after reflection around the plane perpendicular to e2=ey and the original symmetric C tensor:"] 
(Ctensor-Cdashtensor)//MatrixForm
Print["Following are the relationships after equating C tensor after reflection to the Ctensor before reflection around the plane perpendicular to e2=ey:"] 
ay=Solve[(Ctensor-Cdashtensor)==zero,Flatten[Ctensor]] 
Ctensory=Ctensor/.ay[[1]];
Print["Following is the resulting C tensor after adopting the previous relations:"] 
Ctensory//MatrixForm
(*Using plane perpendicular to z as a symmetry plane*)
Cdashtensor=Table[Sum[Qz[[i,m]]*Qz[[j,n]]*Qz[[k,o]]*Qz[[l,p]]*Ctensor[[m,n,o,p]],{m,1,3},{n,1,3},{o,1,3},{p,1,3}],{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["Following is the C tensor after reflection around the plane perpendicular to e3=ez:"] 
Cdashtensor//MatrixForm
Print["Following is the difference between the C tensor after reflection around the plane perpendicular to e3=ez and the original symmetric C tensor:"] 
(Ctensor-Cdashtensor)//MatrixForm
Print["Following are the relationships after equating C tensor after reflection to the Ctensor before reflection around the plane perpendicular to e3=ez:"] 
az=Solve[(Ctensor-Cdashtensor)==zero,Flatten[Ctensor]] 
Ctensorz=Ctensor/.az[[1]];
Print["Following is the resulting C tensor after adopting the previous relations:"] 
Ctensorz//MatrixForm
CtensorOrthotropic=Ctensor/.Flatten[{az[[1]],ay[[1]],ax[[1]]}];
Print["This is the resulting C tensor after adopting the orthotropic material symmetries:"] 
CtensorOrthotropic//MatrixForm
(*Order of the stress vector*)
order={{1,1},{2,2},{3,3},{1,2},{1,3},{2,3}};
(*Converting Ctensor to Cmatrix*)
CmatrixOrthotropic=Table[ToExpression[ToString[Row[{"CtensorOrthotropic","[[",Row[Table[If[OddQ[i],Flatten[{order[[j]],order[[k]]}][[(i+1)/2]],","],{i,1,7}]],"]]"}]]],{j,1,6},{k,1,6}];
Print["Following is the resulting C matrix after adopting the orthotropic material symmetries:"] 
CmatrixOrthotropic//MatrixForm
(*Transverse Isotropy*)
Print["Deriving C tensor of a Transversely Isotropic material from C matrix of Orthotropic Materials assumine that e1=ex and e2=ey form the plane of isotropy"] 
(*Rotation around z with 45 Degrees*)
Q=RotationMatrix[45Degree,{0,0,1}];
Cdashtensor=Table[Sum[Q[[i,m]]*Q[[j,n]]*Q[[k,o]]*Q[[l,p]]*CtensorOrthotropic[[m,n,o,p]],{m,1,3},{n,1,3},{o,1,3},{p,1,3}],{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["Following is the resulting C tensor after rotating the orthotropic material C matrix with some angle in the plane of isotropy:"] 
Cdashtensor=FullSimplify[Cdashtensor];
Cdashtensor//MatrixForm
Print["Following are the differences between the previous C tensor and the Orthotropic material C tensor:"] 
(CtensorOrthotropic-Cdashtensor)//MatrixForm
zero=Table[0,{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Vars={};
Do[If[CtensorOrthotropic[[i,j,k,l]]==0,,,AppendTo[Vars,CtensorOrthotropic[[i,j,k,l]]]],{i,3,1,-1},{j,3,1,-1},{k,3,1,-1},{l,3,1,-1}];
Vars;
Print["Following are the relationships after equating C tensor after rotating in the plane of isotropy to the Ctensor of orthotropic material:"] 
a=Solve[(CtensorOrthotropic-Cdashtensor)==zero,Vars] CtensorTransverseIsotropic=CtensorOrthotropic/.a[[1]];
Print["Following is the C tensor for a Transversely Isotropic Material:"] 
CtensorTransverseIsotropic=FullSimplify[CtensorTransverseIsotropic];
CtensorTransverseIsotropic//MatrixForm
(*Order of the stress vector*)
order={{1,1},{2,2},{3,3},{1,2},{1,3},{2,3}};
(*Converting Ctensor to Cmatrix*)
CmatrixTransverseIsotropic=Table[ToExpression[ToString[Row[{"CtensorTransverseIsotropic","[[",Row[Table[If[OddQ[i],Flatten[{order[[j]],order[[k]]}][[(i+1)/2]],","],{i,1,7}]],"]]"}]]],{j,1,6},{k,1,6}];
CmatrixTransverseIsotropic=FullSimplify[CmatrixTransverseIsotropic];
Print["Following is the C matrix for a Transversely Isotropic Material:"] 
CmatrixTransverseIsotropic//MatrixForm
(*Full Isotropy*)
Print["Deriving C tensor of an Isotropic material from the previous Transversely Isotropic Materials C tensor"] 
(*Rotation around x with 45 Degrees*)
Q=RotationMatrix[45Degree,{1,0,0}];
Cdashtensor=Table[Sum[Q[[i,m]]*Q[[j,n]]*Q[[k,o]]*Q[[l,p]]*CtensorTransverseIsotropic[[m,n,o,p]],{m,1,3},{n,1,3},{o,1,3},{p,1,3}],{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["Following is the resulting C tensor after rotating the Transversely Isotropic material C matrix with some angle around x:"] 
Cdashtensor=FullSimplify[Cdashtensor];
Cdashtensor//MatrixForm
Print["Following are the differences between the previous C tensor and the Transversely Isotropic material C tensor:"] 
(CtensorTransverseIsotropic-Cdashtensor)//MatrixForm
zero=Table[0,{i,1,3},{j,1,3},{k,1,3},{l,1,3}];
Print["Following are the relationships after equating C tensor after rotating to the Ctensor of transversely isotropic material:"] 
a=Solve[(CtensorTransverseIsotropic-Cdashtensor)==zero,Vars] 
CtensorIsotropic=CtensorTransverseIsotropic/.a[[1]];
Print["Following is the C tensor for an Isotropic Material:"] 
CtensorIsotropic=FullSimplify[CtensorIsotropic];
CtensorIsotropic//MatrixForm
(*Order of the stress vector*)
order={{1,1},{2,2},{3,3},{1,2},{1,3},{2,3}};
(*Converting Ctensor to Cmatrix*)
CmatrixIsotropic=Table[ToExpression[ToString[Row[{"CtensorIsotropic","[[",Row[Table[If[OddQ[i],Flatten[{order[[j]],order[[k]]}][[(i+1)/2]],","],{i,1,7}]],"]]"}]]],{j,1,6},{k,1,6}];
CmatrixIsotropic=FullSimplify[CmatrixIsotropic];
Print["Following is the C matrix for an Isotropic Material:"] 
CmatrixIsotropic//MatrixForm

Quiz-Elasticity Tensor

Solution Guide

Leave a Reply

Your email address will not be published.