Best theory diagrams for multilayered structures via shell finite elements

Composite structures are convenient structural solutions for many engineering fields, but their design is challenging and may lead to oversizing due to the significant amount of uncertainties concerning the current modeling capabilities. From a structural analysis standpoint, the finite element method is the most used approach and shell elements are of primary importance in the case of thin structures. Current research efforts aim at improving the accuracy of such elements with limited computational overheads to improve the predictive capabilities and widen the applicability to complex structures and nonlinear cases. The present paper presents shell elements with the minimum number of nodal degrees of freedom and maximum accuracy. Such elements compose the best theory diagram stemming from the combined use of the Carrera Unified Formulation and the Axiomatic/Asymptotic Method. Moreover, this paper provides guidelines on the choice of the proper higher-order terms via the introduction of relevance factor diagrams. The numerical cases consider various sets of design parameters such as the thickness, curvature, stacking sequence, and boundary conditions. The results show that the most relevant set of higher-order terms are third-order and that the thickness plays the primary role in their choice. Moreover, certain terms have very high influence, and their neglect may affect the accuracy of the model significantly.


Introduction
The finite element method (FEM) is one of the most common tools for the design of structures and makes use of three-dimensional (3D), 2D and 1D elements to solve a broad variety of linear and nonlinear structural problems. 2D and 1D elements, although less accurate than 3D, can lead to reduced computational costs. 2D models are referred to as shell and plate finite elements (FE) and can model metallic and composite thin-walled structures. 2D models available in commercial codes rely on the classical theories of structures [1][2][3]. In a 2D model, the primary unknown variables depend on two coordinates, x, and y. On the other hand, assumed fields define the unknown distributions along the thickness direction, z. A structural theory has a given expansion of the unknowns along z. Such expansions characterize the accuracy of a theory and its computational costs. For instance, in FEM, the expansion terms, referred to as generalized unknown variables, define the nodal degrees of freedom (DOF) of the model [4].
Richer expansions lead to higher accuracies and computational costs [5] but wider application scenarios. For a given accuracy level, the choice of a structural theory is problem dependent. In the case of composite structures, the following characteristics may require structural models with richer expansions than classical ones [6]: 1. Moderately thick or thick structures, i.e., a h < 50, where a is the characteristic length of the structure and h is the thickness. 2. Materials with high transverse deformability, e.g., common orthotropic materials, in which E L E T , E L E z > 5, and G E L < 1 10 , where E and G are the Young and shear moduli and L is the fiber direction of the fiber and T, z are perpendicular to L. 3. Transverse anisotropy due, for instance, to the presence of contiguous layers with different properties.
As well-known, such factors require the proper modeling of shear and normal transverse stresses, and variations of the displacement field at the interface between two layers with different mechanical properties, i.e., the Zig-Zag effect. The development of structural theories, i.e., the selection of the expansion terms, can follow two main approaches, namely, the axiomatic and asymptotic ones. The axiomatic method introduces expansions related to hypotheses on the mechanical behavior to reduce the mathematical complexity of the 3D differential equations of elasticity as in the case of classical theories [1][2][3]. The asymptotic method introduces a mathematically rigorous expansion having known accuracy if compared to the 3D exact solution [7,8]. Axiomatic models are easier than asymptotic ones to implement but may miss fundamental expansion terms. Asymptotic models are more rigorous but the simultaneous consideration of multiple problem parameters, e.g., thickness and orthotropic ratio, may be cumbersome.
Via the axiomatic/asymptotic method (AAM), this paper presents best theory diagrams (BTD) [37] providing the shell finite elements with the minimum computational cost and maximum accuracy for a given problem. In [37], the results stemmed from strongform solutions restricting the analysis concerning boundary conditions and stacking sequences. This paper is the first contribution based on shell finite elements allowing the generation of BTD for various boundary conditions and stacking sequences. Moreover, this paper presents a novel metric referred to as Relevance Factor (RF) to evaluate the influence of terms and outline guidelines for the proper choice of the expansion terms. This paper is organized as follows: the governing equations and the methodology are in "Finite element formulation" and "Best theory diagram" sections, then, the "Results" and "Conclusions" sections follow.

Finite element formulation
The Carrera Unified Formulation (CUF) defines the displacement field of a 2D model as where the Einstein notation acts on τ . u is the displacement vector, (u x u y u z ) T . F τ are the thickness expansion functions. u τ is the vector of the generalized unknown displacements. M is the number of expansion terms. In the case of polynomial, Taylor-like expansions, a third-order model, hereinafter referred to as N = 3, has the following displacement field: The third-order model has twelve nodal unknowns. The order and type of expansion is a free parameter. In other words, the theory of structure is an input of the analysis. This paper makes use of the Equivalent Single Layer (ESL) formulation and N = 4 as reference model to build the BTD. The metric coefficients H k α , H k β and H k z of the k-th layer of the multilayered shell are R k α and R k β are the principal radii of the middle surface of the k-th layer, A k and B k the coefficients of the first fundamental form of k , see Fig. 1. This paper focused only on shells with constant radii of curvature with A k = B k = 1. The geometrical relations are where The stress-strain relations are The governing equations make use of the Principle of Virtual Displacements (PVD) and the finite element formulation exploits the MITC technique via nine-node shell elements. The displacement vector and its virtual variation are u τ i and δu s j are the nodal displacement vector and its virtual variation, respectively. Considering the constitutive and geometrical equations, and the PVD, the following governing equation holds The 3 × 3 matrix k k τ sij is the fundamental mechanical nucleus whose expression is independent of the order of the expansion. p k sj is the load vector. More details regarding the finite element formulation are in [4].

Best theory diagram
One of the CUF extensions is the AAM as a tool to analyze the influence of expansion terms starting from a full axiomatic theory [38,39], in this paper, the N = 4. Via the AAM, asymptotic-like results related to the relevance of each variable are obtainable by varying problem parameters, e.g., thickness, orthotropic ratio, stacking sequence, boundary conditions. The AAM can follow various approaches, the one used in this paper has the following steps: 1. Definition of parameters such as geometry, boundary conditions, materials, and layer layouts. 2. Axiomatic choice of a starting theory and definition of the starting nodal unknowns.
Usually, the starting theory provides 3D-like solutions. 3. The CUF generates the governing equations for the theories considered. In particular, the CUF generates reduced models having combinations of the starting terms as generalized unknowns. 4. For each reduced model, the accuracy evaluation makes use of one or more control parameters, in this paper, the maximum transverse displacement.
The number of active terms and the error identifies each theory on a Cartesian plane in which the abscissa reports the error and the ordinate reports the number of active terms. The best theory diagram (BTD) is the curve composed of all those models providing the minimum error with the least number of variables, see Fig. 2. Given the accuracy, models with fewer variables than those on the BTD do not exist. Given the number of variables, models with better accuracy than those on the BTD do not exist. The graphic notation makes use of black and white triangles to indicate active and inactive terms, respectively. In this paper, the control parameter for the error evaluation is the maximum u z , that is,

Results
The numerical results focus on cases retrieved from [40]. The shell has a = b and R α = R β = R. The load is bi-sinusoidal and applied on the top surface, p z =p z sin(πα/a) sin(πβ/b). The material properties are E 1 /E 2 = 25, G 12 /E 2 = G 13 /E 2 = 0.5, G 13 /E 2 = 0.2, ν = 0.25. The finite element model of a quarter of shell has a 4 × 4 mesh as this discretization provides sufficiently accurate results [40]. In all cases, the BTD vertical axis ranges from 5   to 15 since, more often than not, models with 4 or less DOF provide very high errors and are not of practical interest.

Simply-supported, 0/90/0
A simply-supported shell with symmetric lamination is the first numerical case. The analysis aims to study the influence of the thickness and curvature on the BTD. R/a and a/h vary to consider deep, shallow, thick and thin shells.   Table 1 presents the transverse displacement with comparisons against a 3D solution and an analytical model based on a fourth-order layer-wise model. As well-known, the accuracy of the present N = 4 model decreases for thicker shells. However, given that the present work aims to investigate the role of higher-order terms and build BTD, the present N = 4 model accuracy is satisfactory. Figure 3 presents the accuracy of all models stemming by the 2 15 combinations of the N = 4 model. In other words, each dot provides     sake of brevity, R/a = 2 is not reported since does not present any significant changes if compared to R/a = 5. Each row shows the model providing the minimum error for a given number of DOF. For instance, for R/a = 5, a/h = 100, the 7 DOF BTD has the following displacement model: The last row of each table shows the relevance factor (RF) of given order terms in the BTD. The RF is the ratio between the number of active instances and the total number of cases. For instance, RF 0 = 1 indicates that the zeroth-order terms are always present in the BTD. RF 4 = 0.27 because fourth-order terms are in the BTD nine times out of 33 cases. The RF provides a metric to measure the influence of a set of variables, higher the RF higher the relevance. Table 5 reports the error from all the 14 DOF models. Each row refers to a model indicated by the inactive term. The results suggest that • In all cases, no more than six DOF are necessary to provide errors lower than 1%.
• The analysis of all combinations shows that for thin shells there is a significant gap between models providing acceptable accuracies and those with errors larger than 70%. On the other hand, as the thickness increases, the distribution has fewer accuracy gaps. As shown in Table 5, the zeroth and first-order terms affect the gap width to a great extent. In thin shells, their role is predominant, whereas, in thick shells, higherorder terms gain relevance. A more regular accuracy distribution is an indication of more relevance of higher-order terms. • According to the distributions of accuracy from all combinations, the introduction of new terms in an expansion is ineffective if a very relevant term is not present. For instance, u x4 gains significance as the thickness increases. • For thin shells, the FSDT provides higher accuracy with less DOF than the BTD due to the correction of the Poisson locking. For moderately thick shells, a/h = 10, the FSDT matches the BTD but with moderate accuracy. The use of 6 DOF improves the accuracy to a great extent. As a/h decreases further, the FSDT is no longer on the BTD. • The N = 3 is always on the BTD, whereas the N = 2 is a BTD only for thin shells.
• The thickness ratio influences the BTD more than curvature. • The zeroth-order terms are active in each BTD independently of the thickness, i.e., RF 0 = 1. • The relevance of first-and second-order terms decreases as the thickness increases.
• The influence of third-order terms increases as the thickness increases. • The fourth-order terms are the least influential, although, at a/h = 5, the RF increases considerably to the level of second-order terms. • Most of the zeroth-, first-and third-order terms present a regular pattern along a BTD table, i.e., as one of these terms becomes inactive, it does not appear in the BTD anymore. On the other hand, second-and fourth-order terms have a more irregular pattern indicating that their influence depends on the activation or deactivation of other terms.

Simply-supported, 0/90/0/90
The second numerical case deals with a different stacking sequence to investigate the effect of an asymmetric lamination on the BTD. All other parameters remain as those of the previous case. Moreover, this section considers two additional R/a values, 100 and 50, for a more comprehensive analysis on the effect of the curvature. Table 6 presents the transverse displacement values with comparisons with other models from literature, when available. The all combination accuracy plot is in Fig. 6, whereas, Figs. 7 and 8 present the BTD for given R/a values and varying a/h, and Fig. 9 shows the BTD for a given a/h and varying R/a. The BTD models are in Tables 7, 8, 9, 10, and 11. The results suggest that • The present case has more uniform accuracy distributions than the previous one indicating higher relevances of the higher-order terms. For the thick case, there are no relevant gaps up to 60%, and the proper choice of terms can provide any accuracy level. For a/h = 10, there is an accuracy gap between 20 and 35% meaning that there are not structural models that can provide such level of accuracy. • As in the previous case, the FSDT validity is confirmed for the thin case, whereas its accuracy is not sufficient from a/h = 10 and below. • Unlike the previous case, from a/h = 10 and below, some ten DOF are necessary to have errors lower than 1%. • As the thickness increases, the RF distributions are similar to the previous case with a slightly higher influence of the higher-order terms and lower for zeroth-and firstorder ones. • For a/h = 5 the influence of higher-order terms is of particular relevance. For instance, the 5 DOF BTD differs significantly from the FSDT and requires third-order terms. • The variation of the curvature leads to less significant modifications of the BTD than the thickness.

Clamped-free, 0/90/0/90
The last numerical example proposes the 4-layer shell with two edges parallel to β clamped, and the other two free. The aim is to provide some insights into the effect of the geometrical boundary conditions on the BTD. All the other parameters are as in the previous case. Table 12 presents the transverse displacement from the N = 4 model. The BTD for the present case are in Tables 13, 14, 15, and Fig. 10. The results show that the most relevant   effect from the new set of boundary conditions is an increased relevance of higher-order terms at a/h = 10.

Analysis of the relevance of generalized displacement variables
This section aims at investigating the role of each generalized unknowns in the BTD and how their relevance changes with varying parameters. To this purpose, the RF restricts      to each variable as shown, for instance, in Fig. 11. In this case, RF = 1 means that a given variable is present in each BTD of the shell configuration considered. For instance, for the 0/90/0 case with R/a = 5, u x1 , u y1 and u z1 are in all BTD independently of the thickness ratio. Each set of figures presents the RF for the three terms of a given order, see Figs. 11, 12, 13, 14, and 15. The discussion for each order follows.
• Zeroth-order terms As expected, these terms have very high influence and are almost always present in BTD. Just u x1 presents RF lower than unity in three cases in which the 5 DOF BTD requires higher-order terms as discussed in previous sections. • First-order terms In-plane components have unitary RF in most cases. On the other hand, the out-of-plane component has lower relevance and is consistent with the appearance of the FSDT model as 5 DOF BTD for thin and moderately thick shells. • Second-order terms The influence of these terms varies consistently. u x3 has little relevance in the 0/90/0 case but higher in the asymmetric case, and such relevance tends to increase for higher thickness, and the curvature does not influence it. The u y3 relevance has smaller variations due to the thickness change. The thickness strongly influences the out-of-plane component and its influence decreases for thicker shells. • Third-order terms The in-plane components have significant influence which increases for thicker shells. The out-of-plane influence is relevant in the symmetric case and increases for thicker shells. • Fourth-order terms These terms are the less influential except for u x5 in the clampedfree case. The relevance of these terms should increase as soon as the BTD considers stress distributions.

Conclusions
This paper presented results on the accuracy of higher-order generalized displacement variables for composite shells. Investigations used the CUF and the AAM. The former provided the finite element matrices for any-order models, and the latter led to the analysis of the relevance of each generalized variable. The combined use of these tools generated the BTD and relevance factor diagrams. The BTD provides the minimum number of DOF for a given accuracy level. The RF diagrams measure the importance of a variable, or of a set of variables, as various parameters change, e.g., thickness, curvature, stacking sequence. All results considered the maximum transverse displacement as the control parameter. The analysis led to the following guidelines and recommendations: • For the cases considered in this paper, the thickness and stacking sequence are the most important factors for the choice of the primary variables. For thin shells, six DOF are sufficient to obtain errors lower than 1%. For thick shells, ten DOF are necessary. • In most cases, the accuracy level obtainable from combinations of a given set of variables is not continuous as the DOF decrease. In other words, there are no structural models that can satisfy certain accuracy of the solution. • Accuracy gaps indicate the presence of very effective terms that must be present in the expansion to ensure satisfactory accuracies. For instance, for thin shells, these terms coincide with the FSDT expansions. However, as the presence of non-classical effects due to asymmetries or high thickness increases, the relevance of higher-order terms increases and the accuracy gaps tend to disappear. • The FSDT and second-order model are BTD only for thin shells. The third-order model is close to the BTD in most cases. • As the thickness increases, the relevance of third-order variables increases significantly, and these terms can be the most relevant together with the zeroth-order ones. • The out-of-plane displacement variables tend to have less relevance than in-plane ones. Such a relevance should increase significantly as soon as the analysis considers stress distributions as control parameters. • The set of variables composing a BTD model depends on the boundary conditions; however, such a dependency is weaker than the thickness one.
Most immediate future developments should deal with the inclusion of all displacement and stress components as control parameters and the analysis of more complex configurations. In fact, for the boundary conditions adopted, the use of the transverse displacement is the minimum requirement for a BTD. The inclusion of other control parameters, e.g., transverse shear and axial stresses, may modify the BTD concerning accuracy and set of active terms with higher relevance of higher-order terms.