Distances of Stiffnesses to Symmetry Classes

For a given elastic stiffness tetrad an algorithm is provided to determine the distance of this particular tetrad to all tetrads of a prescribed symmetry class. If the particular tetrad already belongs to this class then the distance is zero and the presentation of this tetrad with respect to the symmetry axes can be obtained. If the distance turns out to be positive, the algorithm provides a measure to see how close it is to this symmetry class. Moreover, the closest element of this class to it is also determined. This applies in cases where the tetrad is not ideal due to scattering of its measurement. The algorithm is entirely algebraic and applies to all symmetry classes, although the isotropic and the cubic class need a different treatment from all other classes.


Introduction
Many materials in engineering applications exhibit an inner structure, e.g., crystal lattices or fibers, which causes anisotropy. Such a structure induces certain material symmetries. Knowledge of these symmetries can be useful in order to apply representation theorems or to find the number of independent elastic constants that are needed to describe the elastic behavior fully.
Based on material symmetries, materials are often categorized into symmetry classes. All material laws within the same symmetry class share the same symmetry group. However, it is a complicated task to determine the symmetry of a particular material. In the literature, see [14], the symmetry class of a stiffness is usually determined by the form of its VOIGT representation. Since the representation is non-unique, this method of identification is only applicable if the representation uses a particular base called symmetry base. Furthermore, measured stiffnesses always show a certain amount of scattering. As a result, such stiffnesses appear to be triclinic (the "most anisotropic" case) even if represented in a special base.
A practical method to determine the symmetry class of a material and the corresponding stiffness is needed. The method presented in this paper is able to identify the (exact) symmetry class of a material if it exists. Furthermore, in the presence of scattering, this method will determine the closest non-triclinic symmetry class of a stiffness. This is achieved by using the distance between the closest element in a non-triclinic symmetry class and the stiffness.
The determination of the symmetry class of a stiffness as well as other closely related topics have been discussed for a long time. In [10] the closest cubic stiffness is calculated by projecting a stiffness onto a cubic base which has some degrees of freedom. These degrees of freedom are determined by minimizing the residual. Furthermore, GAZIS et al. generalize the concept to find closest stiffnesses in others than the cubic class.
An approach to determine the symmetry class of a stiffness based on MAXWELL multipoles is shown in [1]. Here symmetry axes are found by finding MAXWELL multipoles of the parts of the harmonic decomposition of the stiffness. These axes allow for the classification of the symmetry class.
A method for symmetry classification based on eigenvalues of the stiffness is presented by [8]. Furthermore, the material symmetries are characterized by planes of symmetry. As a side result, in an article concerning the harmonic decomposition of stiffnesses, Baerheim [2] provides an algorithm to calculate the rotation which transforms a stiffness from an arbitrary base into the symmetry base where the typical VOIGT representations are observable. However, BAERHEIM points out that the algorithm is not applicable to cubic stiffnesses.
The identification of measured stiffnesses is the topic of [9]. They introduce the concepts of a nearest stiffness in a given symmetry class as well as of the distance between a stiffness and a symmetry class.
Another method to identify ideal stiffnesses using mirror planes is proposed in [16]. This method is derived by application of an irreducible decomposition of the stiffness tensor as well as MAXWELL multipoles. A characteristic function is introduced the roots of which determine the symmetry class of the material.
In [15] the concepts of the nearest stiffness in a given symmetry class and distance of a stiffness to a symmetry class are used. However, the rotation of the stiffness into a special base is performed using a numerical optimization scheme. The resulting method is applied to fiber-reinforced composites.
In our work we present an algorithm for the identification of the symmetry class of a stiffness is proposed which is based on [2,9,15]. Instead of using a numerical optimization scheme we employ the algorithm presented by BAERHEIM to calculate the rotation from an arbitrary into a special base. The cubic case is integrated into the procedure by applying an eigenvalue method. Using the rotated stiffness representation, the closest stiffness in a given symmetry class can be determined. The relative deviation between the rotated and the closest stiffness is calculated and considered as distance of the stiffness to the symmetry class. This approach allows for the classification of ideal as well as measured stiffnesses. Furthermore, our method presents a fully algebraic algorithm for the classification problem which is applicable to all symmetry classes.

Notations and Definitions
We employ a direct tensor notation. Scalars, vectors, dyadics, and nth-order tensors are represented by a, a, A and n A, respectively. The scalar product for vectors is written as a · b.
Furthermore, the dyadic product of two arbitrary vectors a and b is noted as a ⊗ b. The n-fold contraction between two higher-order tensors, (2.1) If m = n = p this coincides with the usual scalar product of higher-order tensors.
The Rayleigh product between a second-order tensor Q and the nth-order tensor n B is defined by with respect to any vector base {e i }. In the present context, this product will be often used for a proper orthogonal Q resulting in a rotation of n B. For n being even, the following identity Consequently, for tensors of even order the distinction between proper orthogonal and orthogonal tensors in the RAYLEIGH product is irrelevant.
To denote a specific rotation around an axis a by an angle α the notation R α a is used. In the sequel, we will consider fourth-order tensors which represent the elastic stiffness of a material. For them we assume the usual major symmetry and the left and right subsymmetry. Such tetradics form the following set K. Such a group can be discrete or continuous. However, even if the symmetry group is continuous as for the isotropic and transversally isotropic case, a discrete subset can be used to completely fix the symmetry of the stiffnesses. In these cases it is possible to generate a finite subset for the isotropic case by adding another rotation to the generators of the cubic symmetry group and consider the resulting symmetry transformations, see [5] or [12].
In the transversally isotropic case, one rotation around a canonic axis by the angle 2π /5, generates a sufficient subset of the transversally isotropic symmetry group. However, the hexagonal symmetry group can be used for brevity.
It may happen that the symmetry groups of two stiffness tetradics differ only by some rotation R ∈ SO(3) In such cases we call the two symmetry groups isomorphic. The isomorphic symmetry group can be considered as the orbit of G 4 K under SO(3). For a given symmetry group, the symmetry class is defined as the set of all stiffnesses with either the same or some isomorphic symmetry group, i.e., (2.5)

Symmetry Bases of Stiffnesses
If an ideal stiffness is represented in the symmetry base, its symmetry class can be identified by the form of the component matrix in VOIGT-notation. The typical forms of the component matrices can be found in various publications like, e.g., [12,14]. The symmetry base is often referred to as the crystal base. However, this notion applies also in cases without crystal lattice but some other structure inducing anisotropy like fiber-reinforced composites.
Usually the symmetry base of a given stiffness is not known. This yields that the symmetry class can not be immediately identified by the form of the VOIGT-representation. Therefore a rotation is needed that transforms the given VOIGT-representation into a representation for a certain symmetry. If we only consider orthonormal vector bases, such a transformation is given by some orthogonal tensor the determination of which is non-trivial.
An algorithm to determine this tensor has been provided by [2], based on [1] and [7], which applies only for the non-cubic cases as we will discuss in the sequel. It makes use of the fact that every tensor can be decomposed into a totally symmetric part and a residual. Consider the tensor of order n in a given base {e i }. For brevity we introduce the index set Next we consider permutations of some given index set I . We denote a permutation of the index set as σ i (I ) like, e.g., σ 2 (I ) = {i 2 , i 1 , . . . , i n }. It is known, that there are n! of such permutations. The totally symmetric part of a tensor of nth-order is defined as Similarly, the residual is defined as To further explain this concept consider as an example the fourth-order tensor The symmetric part can be calculated as However, if the fourth-order tensor shows the major and the two subsymmetries as defined in Eq. (2.3), the calculation simplifies to (3.7)

Decomposition of the Stiffness Tensor
Every stiffness 4 K ∈ S can be decomposed into a totally symmetric and a residual part, i.e., As shown in [2], the residual can be constructed from a symmetric second-order tensor T as Note the cyclic permutations that need to be performed with the bases vectors in order to obtain the correct result. The symmetric second-order tensor can be calculated as It was shown by BAERHEIM that the base of a stiffness 4 K is a symmetry base if and only if it is the spectral base of T. This fact can be used to find the transformation Q ∈ SO(3) that transforms a non-symmetry into a symmetry base.
However, this is not to be confused with finding the symmetry class of 4 K. Finding the rotation Q allows us to rotate the stiffness into a base where the symmetry transformations in the symmetry groups are known and it is therefore a pre-processing step for finding the symmetry class. The decision to which symmetry class a stiffness belongs is based on a distance measure and thus completely independent of any part of the harmonic decomposition.
The transformation Q ∈ SO(3) can be calculated by following the steps below. K to get the stiffness in a symmetry base.
The symmetry base is not unique for two reasons. First, one needs to choose an ordering of eigenvalues which will affect the transformation. Second, while calculating the eigenvectors one can choose to multiply any eigenvector by −1. The result is still an eigenvector.
As we have already mentioned this method is not applicable in the case of a cubic symmetry. This can be shown as follows. Every cubic stiffness can be written as with the three projectors see [3,4]. This yields It is obvious that the residual is an isotropic tensor in case of a cubic stiffness i.e., rotations of cubic stiffnesses can not be detected. In order to calculate the stiffness components in a symmetry base for the cubic case, another method must be employed. One may consider the following scalar products, We substitute the two sets of candidates into the characteristic equation using VOIGTnotation which is only zero for the eigenvalues. After identifying the eigenvalues, the stiff-ness components in the symmetry base can be calculated as (3.20) Both, the BAERHEIM-method and the eigenvalue method for cubic stiffnesses provide algorithms to find the symmetry base for every stiffness. In the case of ideal stiffnesses the symmetry class of a stiffness in a symmetry base can be identified by means of typical forms of the VOIGT-matrices, see, e.g., [14].

Distance Between a Stiffness and a Symmetry Class
The identification of the symmetry class by means of standard forms that we described in the previous section is only possible for ideal stiffnesses, i.e., stiffnesses without scattering. However, stiffnesses obtained by measurements or numerical experiments are in general non-ideal. A way to associate a symmetry class with non-ideal stiffnesses is to find the nearest symmetry class by calculating the distance to all symmetry classes. Since a symmetry class is a set of stiffnesses with certain properties, it is not trivial to define the distance between a stiffness and a symmetry class. After [15], we are defining the distance as the distance between a stiffness and the nearest stiffness in a given symmetry class. A formula for the latter is provided by [9], where |G| is the cardinality of G, i.e., the number of elements in the group. It is introduced as an average of the orbit of a stiffness in a symmetry group. Earlier GAZIS et al. [10] introduced a similar calculation procedure using an integral average. Equation (4.1) is applicable for the triclinic, monoclinic, orthotropic, tetragonal, cubic, trigonal, and hexagonal symmetry classes. The isotropic and the transversally isotropic case need a different solution.
In order to calculate the distance of a stiffness to the isotropic class, the isotropic part of the stiffness is used, In the transversally isotropic case, the nearest hexagonal stiffness is used since it is always equal to the nearest transversally isotropic stiffness. This is due to the fact that the maximal number of constraints on a stiffness due to rotations around one symmetry axis is reached by five rotations with the angle 2π /5. This means that the resulting nearest stiffness does not change if rotations around a symmetry axis with the angle 2π /n, 5 ≤ n ∈ N are considered. In order to avoid the introduction of a pentagonal symmetry class, the hexagonal symmetry class is used. Finally, the distance between a stiffness and a symmetry class is introduced by [15] as with the FROBENIUS-norm of an nth-order tensor The definition of the distance is unit-dependent and generally yields different results if the compliance is used instead of the stiffness. In order to obtain the same value in both cases, a logarithmic distance measure can be used, see [11]. A similar conclusion was already drawn by [6]. However, in the present context the absolute value of the distance does not play an important role. Rather it is the order that is induced by all such distances in exactly the same way. Equation (4.3) involves the nearest stiffness in the symmetry class G. If G is not the isotropic class this means that the symmetry transformations for the employed representation of the stiffness needs to be known. However, the symmetry transformations are generally unknown. Therefore, one needs to use the algorithm by BAERHEIM first in order to ensure that the stiffness is in a symmetry base. According to the symmetry base, the symmetry transformations are known and the distance can be calculated.

Finding the Symmetry Class of an Arbitrary Stiffness
The symmetry classification of stiffnesses has been the topic of discussions in the mechanics community for some decades. The recent paper by WEBER et al. [15] used a numerical optimization scheme to find the canonic base of the stiffness. However, non-linear optimization comes with some computational afford as well as the risk of finding second-degree minima. We therefore present here a completely algebraic solution for the symmetry classification problem.
Such a scheme is not convenient to be performed manually. We therefore employ a computational implementation.

Computational Implementation
In order to find the symmetry class of any given stiffness, the symmetry transformations of the symmetry groups are needed. To calculate these transformations the generators in Table 1 are used. We treat the transversally isotropic case as hexagonal since the restrictions concerning the stiffness are the same in both cases. Furthermore, Eq. (4.2) is employed to calculate the nearest isotropic stiffness. In order to classify the symmetry class of a stiffness 4 K BAERHEIM's algorithm is employed. This yields one of twelve equivalent canonic stiffnesses. These equivalent stiffnesses are connected via the transformations of the cubic symmetry group. The cubic transformations rotate a base in a way that the CARTESIAN base vectors are mapped onto one another or their negatives. The cubic transformations are given by P 1 = 1, P 2 = R π/2 e 1 , P 3 = R 3π/2 e 1 , P 4 = R π e 1 , P 5 = R π/2 e 2 , P 6 = R 3π/2 e 2 , P 7 = R π e 2 , P 8 = R π/2 e 3 , where d = (e 1 +e 2 +e 3 ) / √ 3. Since we do not know which orientation we are getting for any given stiffness, all twelve stiffnesses need to be considered. Furthermore, all orientations of the considered stiffness are normalised using the FROBENIUS-norm, i.e., The distance between all normalised orientations and stiffness is calculated by using Eq. (4.3). For isotropic and transversally isotropic stiffnesses the aforementioned workarounds are applied. The minimum of the distances between the different orientations of the stiffness to one symmetry class is considered as the distance of the stiffness to the symmetry class, i.e., The computational implementation is further illustrated by Fig. 1.

Examples
In the following paragraph purely academic examples for the identification of the symmetry class of stiffnesses are given. The unit of the stiffnesses is suppressed for brevity.

Stiffness with Scattering in a Non-canonic Base
We are going to look at an example that could be encountered while measuring a stiffness either experimentally or numerically. Consider the randomly generated stiffness It is obvious that there is no observable structure in the non-canonic form of the stiffness. In order to find the symmetry class of the stiffness we firstly calculate the residual part of the It is worth mentioning that this result is non-unique. One is free to choose the ordering of the calculated eigenvectors. As a result, the columns of the component matrix of the canonic transformation can be ordered in several ways. Using As explained in Paragraph 5.1, this is one of twelve equivalent representations. It is obvious that a structure is beginning to emerge. However, there appears to be some scattering in the stiffness. To identify the symmetry class of the stiffness, the distance to all the symmetry classes needs to be considered. The distances to various classes of symmetry are shown in Fig. 2. Since triclinic is the maximal symmetry class the distance of any stiffness to the triclinic class is always zero. Furthermore, it is noteworthy that the distances to the monoclinic and the orthotropic class are small and nearly the same. The distance to the monoclinic class is smaller since it includes the orthotropic class but is able to take more of the scattering into account. However, it can be concluded that 4 K ex1 is almost orthotropic and needs nine independent constants to describe the materials behaviour.

Ideal Cubic Stiffness in a Non-canonic Base
Consider the experimental measurement of the stiffness of a material which is known to be cubic. Due to a misalignment in the measuring apparatus the following stiffness is measured  It is obvious that this is the VOIGT-representation of a cubic stiffness.

Conclusion
We presented a method to assign a symmetry class to ideal as well as measured stiffnesses. The identification is based on the distance of a stiffness to a symmetry class to allow for taking scattering in measured stiffnesses into account. As an intermediate step, the stiffness is rotated into a special base by either the BAERHEIM or the eigenvalue method (cubic case). This special representation is then used to find the closest element to the stiffness in a given symmetry class. Finally, the relative distance between the stiffness and the closest element in the given symmetry class is calculated as distance of the stiffness to the symmetry class. The procedure has been illustrated by examples of relevant cases.
Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.