## Laminated composite shells: buckling of a cylindrical panel with a circular hole with Abaqus

This example illustrates modeling a thin, laminated composite shell in the presence of buckling.

*Animation-1.*

**Geometry and model**

The structure analyzed is shown in Figure 1 and was originally studied experimentally by Knight and Starnes (1984). The test specimen is a cylindrical panel with a 355.6 mm (14 in) square platform and a 381 mm (15 in) radius of curvature, so that the panel covers a 55.6° arc of the cylinder. The panel contains a centrally located hole of 50.8 mm (2 in) diameter. The shell consists of 16 layers of unidirectional graphite fibers in an epoxy resin. Each layer is 0.142 mm (.0056 in) thick. The layers are arranged in the symmetric stacking sequence {45/90/0/0/90/45} degrees repeated twice. The nominal orthotropic elastic material properties as defined by Stanley (1985) are

where the 1-direction is along the fibers, the 2-direction is transverse to the fibers in the surface of the lamina, and the 3-direction is normal to the lamina.

The panel is fully clamped on the bottom edge, clamped except for axial motion on the top edge and simply supported along its vertical edges. Three analyses are considered. The first is a linear (prebuckling) analysis in which the panel is subjected to a uniform end shortening of 0.8 mm (.0316 in). The total axial force and the distribution of axial force along the midsection are used to compare the results with those obtained by Stanley (1985). The second analysis consists of an eigenvalue extraction of the first five buckling modes. The buckling loads and mode shapes are also compared with those presented by Stanley (1985). Finally, a nonlinear load-deflection analysis is done to predict the postbuckling behavior, using the modified Riks algorithm. For this analysis an initial imperfection is introduced. The imperfection is based on the fourth buckling mode extracted during the second analysis. These results are compared with those of Stanley (1985) and with the experimental measurements of Knight and Starnes (1984).

The mesh used in Abaqus is shown in Figure 2. The anisotropic material behavior precludes any symmetry assumptions, hence the entire panel is modeled. The same mesh is used with the 4-node shell element (type S4R5) and also with the 9-node shell element (type S9R5); the 9-node element mesh, thus, has about four times the number of degrees of freedom as the 4-node element mesh. The 6-node triangular shell element STRI65 is also used; it employs two triangles for each quadrilateral element of the second-order mesh. Mesh generation is facilitated by specifying node fill and node mapping, as shown in the input data. In this model specification of the relative angle of orientation to define the material orientation within each layer, along with orthotropic elasticity in plane stress, makes the definition of the laminae properties straightforward.

The shell elements used in this example use an approximation to thin shell theory, based on a numerical penalty applied to the transverse shear strain along the element edges. These elements are not universally applicable to the analysis of composites since transverse shear effects can be significant in such cases and these elements are not designed to model them accurately. Here, however, the geometry of the panel is that of a thin shell; and the symmetrical lay-up, along with the relatively large number of laminae, tends to diminish the importance of transverse shear deformation on the response.

**Relation between stress resultants and generalized strains**

The shell section is most easily defined by giving the layer thickness, material, and orientation, in which case Abaqus preintegrates to obtain the section stiffness properties. However, the user can choose to input the section stiffness properties directly instead, as follows.

In Abaqus a lamina is considered as an orthotropic sheet in plane stress. The principal material axes of the lamina (see Figure 3) are longitudinal, denoted by L; transverse to the fiber direction in the surface of the lamina, denoted by T; and normal to the lamina surface, denoted by *. The constitutive relations for a general orthotropic material in the principal directions (**) are*

In terms of the data required to define orthotropic elasticity by specifying terms in the elastic stiffness matrix in Abaqus these are

This matrix is symmetric and has nine independent constants. If we assume a state of plane stress, then * is taken to be zero. This yields*

where

The correspondence between these terms and the usual engineering constants that might be given for a simple orthotropic layer in a laminate is

The parameters used on the right-hand side of the above equation are those that must be provided as part of the definition of orthotropic elasticity in plane stress.

If the (*) system denotes the standard shell basis directions that Abaqus chooses by default, the local stiffness components must be rotated to this system to construct the lamina’s contribution to the general shell section stiffness. Since ** represent fourth-order tensors, in the case of a lamina they are oriented at an angle ** to the standard shell basis directions used in Abaqus. Hence, the transformation is*

where *Qij* are the stiffness coefficients in the standard shell basis directions used by Abaqus.

Abaqus assumes that a laminate is a stack of laminae arranged with the principal directions of each layer in different orientations. The various layers are assumed to be rigidly bonded together. The section force and moment resultants per unit length in the normal basis directions in a given layer can be defined on this basis as

where *h* is the thickness of the layer.

This leads to the relations

where the components of this section stiffness matrix are given by

**Results and discussion**

The total axial force necessary to compress the panel 0.803 mm (0.0316 in) is 100.2 kN (22529 lb) for the mesh of S9R5 elements, 99.5 kN (22359 lb) for the mesh of S4R5 elements, and 100.3 kN (22547 lb) for the mesh of STRI65 elements. These values match closely with the result of 100 kN (22480 lb) reported by Stanley (1985). Figure 5 shows the displaced configuration and a profile of axial force along the midsection of the panel (at . It is interesting to note that the axial load is distributed almost evenly across the entire panel, with only a very localized area near the hole subjected to an amplified stress level. This suggests that adequate results for this linear analysis could also be obtained with a coarser mesh that has a bias toward the hole.

The second stage of the analysis is the eigenvalue buckling prediction. To obtain the buckling predictions with Abaqus, an eigenvalue buckling prediction step is run. In this step nominal values of load are applied. The magnitude that is used is not of any significance, since eigenvalue buckling is a linear perturbation procedure: the stiffness matrix and the stress stiffening matrix are evaluated at the beginning of the step without any of this load applied. The eigenvalue buckling prediction step calculates the eigenvalues that, multiplied with the applied load and added to any “base state” loading, are the predicted buckling loads. The eigenvectors associated with the eigenvalues are also obtained.

The buckling predictions are summarized in Table 1 and Figure 6. The buckling load predictions from Abaqus are higher than those reported by Stanley. The eigenmode predictions given by the mesh using element types S4R5, S9R5, and STRI65 are all the same and agree well with those reported by Stanley. Stanley makes several important observations that remain valid for the Abaqus results: (1) the eigenvalues are closely spaced; (2) nevertheless, the mode shapes vary significantly in character; (3) the first buckling mode bears the most similarity to the linear prebuckling solution; (4) there is no symmetry available that can be utilized for computational efficiency.

Following the eigenvalue buckling analyses, nonlinear postbuckling analysis is carried out by imposing an imperfection based on the fourth buckling mode. The maximum initial perturbation is 10% of the thickness of the shell. The load versus normalized displacement plots for the S9R5 mesh, the S4R5 mesh, and the STRI65 mesh are compared with the experimental results and those given by Stanley in Figure 7. The overall response prediction is quite similar for the Abaqus elements, although the general behavior predicted by Stanley is somewhat different. The Abaqus results show a peak load slightly above the buckling load predicted by the eigenvalue extraction, while Stanley’s results show a significantly lower peak load. In addition, the Abaqus results show rather less loss of strength after the initial peak, followed quite soon by positive stiffness again. Neither the Abaqus results nor Stanley’s results agree closely with the experimentally observed dramatic loss of strength after peak load. Stanley ascribes this to material failure (presumably delamination), which is not modeled in his analyses or in these.

Figure 8 shows the deformed configurations for the panel during its postbuckling response. The plots show the results for S4R5, but the pattern is similar for S9R5 and STRI65. The response is quite symmetric initially; but, as the critical load is approached, a nonsymmetric dimple develops and grows, presumably accounting for the panel’s loss of strength. Later in the postbuckling response another wrinkle can be seen to be developing.

**Animation-2.**

**Tables**

Mode 1 | Stanley | 107.0 kN (24054 lb) |

S9R5 | 113.4 kN (25501 lb) | |

S4R5 | 115.5 kN (25964 lb) | |

S4 | 114.3 kN (25696 lb) | |

STRI65 | 113.8 kN (25579 lb) | |

Mode 2 | Stanley | 109.6 kN (24638 lb) |

S9R5 | 117.6 kN (26429 lb) | |

S4R5 | 121.2 kN (27244 lb) | |

S4 | 116.5 kN (26196 lb) | |

STRI65 | 117.8 kN (26492 lb) | |

Mode 3 | Stanley | 116.2 kN (26122 lb) |

S9R5 | 120.3 kN (27049 lb) | |

S4R5 | 124.7 kN (28042 lb) | |

S4 | 124.1 kN (27889 lb) | |

STRI65 | 121.1 kN (27217 lb) | |

Mode 4 | Stanley | 140.1 kN (31494 lb) |

S9R5 | 147.5 kN (33161 lb) | |

S4R5 | 156.1 kN (35092 lb) | |

S4 | 152.3 kN (34247 lb) | |

STRI65 | 146.9 kN (33015 lb) | |

Mode 5 | Stanley | 151.3 kN (34012 lb) |

S9R5 | 171.3 kN (38512 lb) | |

S4R5 | 181.5 kN (40800 lb) | |

S4 | 184.2 kN (41413 lb) | |

STRI65 | 172.8 kN (38843 lb) |