Composite Structures 73 (2006) 229–236 www.elsevier.com/locate/compstruct
Influence of modelling and solution methods on the FE-simulation of the post-buckling behaviour of stiffened aircraft fuselage panels P. Linde
a,*
, A. Schulz b, W. Rust
b,c
a Airbus, Kreetslag 10, D-21129 Hamburg, Germany CAD-FEM GmbH, Schmiedestraße 31, D-31303 Burgdorf, Germany University of Applied Sciences, Ricklinger Stadtweg 120, D-30459 Hannover, Germany b
c
Abstract Stiffened fuselage panels with laminated constructions play an increasing role in aircraft design. The static behaviour through the buckling- and post-buckling regime until failure has to be established. Apart from analytical calculations, experimental tests for different load combinations are indespensible, both of which are expensive and time consuming. The virtual testing described here is based on a development project aiming at reducing the amount of experimental tests, and narrowing the numerical predictions to experimental results. A tool developed for parametric modelling and simulation of test shells is discussed. The numerical model is based on layered shell elements in ANSYS and for special purposes in LS-DYNA. It is outlined how far the behaviour of laminates (interaction of different and anisotropic materials, delamination and splices) can be simulated in this context. Results are given for welded panels and for fibre metal laminate panels. Comparison with experimental data is made. Recommendations for future research is provided. 2005 Elsevier Ltd. All rights reserved. Keywords: Stiffened shell structure; Metal fibre laminate; FE limit load analysis; Quasi-static dynamic solution; Virtual testing
1. Introduction
1.2. Objective
1.1. General
FML is finding use as fuselage skin material in civil aircraft. Static testing of such panels until failure is carried out in a specially developed test rig. These tests are complex, expensive and very time consuming. For the purpose of reducing the amount of these tests, i.e., not replacing them altogether, a virtual test rig has been developed [1], see also [2,3]. The numerical simulation of the tests in this test rig of laminated panels is the subject of this paper.
In modern aircraft design laminated materials obtain increasing importance. One such family of materials is made up of the fibre metal laminates (FML). This is a hybrid composite consisting of thin aluminium layers alternating with thin layers of glass fibre reinforced epoxy. This material displays improved fatigue behaviour versus monolithic materials, as well as improved fire resistance. The static damage behaviour is also improved compared with other laminated materials, it is however not easy to simulate realistically. *
Corresponding author. Tel.: +49 40 743 82353; fax: +49 40 743 84277. E-mail address:
[email protected] (P. Linde).
0263-8223/$ - see front matter 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.compstruct.2005.11.048
1.3. Scope Following the introduction this paper in the second section treats the stiffened panels, the geometry, and the material modelling. Physical aspects such as splices, delamination and boundary conditions are discussed. Section 3
230
P. Linde et al. / Composite Structures 73 (2006) 229–236
treats the computations. Solution methods discussed include implicit static, implicit dynamic (quasi-static), and explicit dynamic. Results are presented in Section 4 in which comparisons between numerical and experimental results are given. Section 5 finally provide a summary, conclusions and recommendations. 2. Model The panel mainly consists of the skin which is stiffened by stringers. The circumferential stiffness and radial support is obtained from frames which are connected to skin and stringers by a formed sheet metal called clip. Main connections are rivets and welds but bonding is also under consideration. In the test rig the panel is partially stiffened by extra metal layers being bonded to the skin (doublers) or by riveted metal profiles (below referred to as stiffeners). The boundaries are supported by clamps being guided by a very stiff beam which on the other hand is used to apply shear loads. The skin may consist of aluminium, metal laminates or fibre metal laminates (Glare). Other composites are under consideration. In case of homogenous metal regions of smaller thickness are created by removing some material (pocketing). In case of laminates the same effect can be obtained by laminating layers with holes. 2.1. Geometry model Supported by a graphical user interface designed for this special purpose the fully parametric modelling of the geometry of the whole panel is done by the pre-processor of the general purpose FE-tool ANSYS. Afterwards the FE-mesh is created with the density level chosen by the user. This model is mainly used to create a FE-mesh for the ANSYS solver. However, in case of convergence problems an input file for the explicit dynamic FE solver LS-DYNA can also be generated.
Fig. 1. Varying thicknesses.
Furthermore, the element can account for layered cross sections. This is necessary if the behaviour of the laminates should be modelled in detail. The shear correction factors are calculated internally. Each layer can have its own material law. Welds are modelled via common nodes, but with small elements to get an idea of weld stresses. Timoshenko beam elements serve as rivets. Contact is checked wherever necessary. The edge beams are considered rigid, they form rigid targets for contact with shells representing parts of the clamps. The motion of such rigid bodies is constrained at so-called pilot nodes which also are used for single-force loading. Rigid targets are also used for the other parts of the clamps. They are fixed to the edge of the panel by bonded contact, but use normal contact for the support of the skin in the adjacent regions (Fig. 2). 2.3. Model for explicit analysis ANSYS offers a pre-processing tools for LS-DYNA. However, some features are not supported or different concepts are followed. In these cases some programming in the ANSYS Parametric Design Language (APDL) has been done to complete the transfer of the model. Particular differences are
2.2. Model for implicit analysis The static limit-load analysis including all kinds of nonlinearities is done in ANSYS. With only very few exceptions the FE-model consists of the 4-node Reissner–Mindin shell element SHELL181. The bending part is formulated according to Bathe/Dvorˇkin (assumed strain), the membrane part offers the choice between reduced integration (at least used for the skin) and full integration with incompatible modes (mainly used for open-shaped profiles like the stringers). This element allows for a free location of the reference plane (not only mid-plane). This feature is essential because the outer skin is the reference in the real airplane but the thickness can vary for several reasons such as doubling or pocketing, see Fig. 1.
Fig. 2. Rigid targets for load application and constraints.
P. Linde et al. / Composite Structures 73 (2006) 229–236
• Layered cross sections and reference-plane offsets are modelled by user-defined integration rules: the location and the weighting factors of integration points together with material pointers describe the lay-up of the composites; the offset from the mid-plane requires a layer of ‘nothing’, i.e., the sum of the weighting factors is not 1. • Rigid targets are transferred to rigid bodies, then contact can be defined. The pilot nodes can become a user-defined center of gravity (for constraints) or an extra node placed anywhere (for force application). • Contact definition is different: in LS-DYNA the standard is automatic overall contact or part related input, but segment based definition like in ANSYS is also possible. The shell elements used are 4-noded, usually reduced integrated for solution time reason. Like in ANSYS model full integration is chosen for open-shaped profiles. The special LS-DYNA features for spotwelds being useful for the representation of rivets are not used because ANSYS determines the reference idealisation. 2.4. FE-model for LS-DYNA As long as the panel and its components consist of single-layered aluminium sheets anisotropies are not distinct and a usual isotropic hardening plasticity law with a piecewise linear stress–strain curve can be used. The fibre metal laminate consists of alternating layers of aluminium and prepreg. For the prepreg only in situ properties are of importance, i.e., derived from measurements of a complete lay-up. Nevertheless, the prepreg is extremely anisotropic, especially if the 45-direction is compared to 0 and 90, where are still slighter differences. It could be shown by experiments that the properties for different number of layers can be determined from the volume fraction of each material type. Two levels of detailing have been considered: • a smeared model considering the material as homogenous over the cross section and, • a layered model.
231
The Mij terms are derived from the condition that the plastic deformation is incompressible. For shear 3 2 3 0 0 2 7 6 C 12 7 6 7 6 3 shear 2 T6 0 0 7 ð4Þ ðreqv Þ ¼ s 6 2 7s C 23 7 6 4 3 5 0 0 C 213 holds where pffiffiffi 3syi C ij ¼ ry0
ð5Þ
Cij again is a scale factor for different directions but accounts for the different influence of shear and normal stress on yielding. Thus, the scale factor for normal stress in the 45-direction can directly be used instead of shear. Hill’s criterion assumes that the stress–strain curves are scalable which is not exactly the case. Therefore, the reference curve has been chosen after the main loading direction, either shear or axial compression. In a layered model of the laminate it is possible to get good accordance with measurement in both the axial and
Fig. 3. Stress–strain curve for fibre metal laminate, measured and calculated in the FE-model, axial and circumferential direction.
In the two cases the anisotropy is modelled by using Hill’s equivalent stress qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 shear 2 rHill ¼ ðrnorm ð1Þ eqv eqv Þ þ ðreqv Þ 2 1 3 M 12 M 13 C2 6 1 7 1 6 2 T M M 23 7 12 ð2Þ ðrnorm 7r C 22 eqv Þ ¼ r 6 4 5 1 M 13 M 23 C2 3
where Ci is a scale factor for the reference stress–strain curve to determine the curve for the specific direction: ryi Ci ¼ ð3Þ ry0
Fig. 4. Stress–strain curve for fibre metal laminate, measured and calculated in the FE-model, 45 direction.
232
P. Linde et al. / Composite Structures 73 (2006) 229–236
Fig. 5. FE-shell and volume model of a splice: metal (dark blue), prepreg (green), glue (magenta). (For interpretation of references in color, the reader is referred to the web version of this article.)
the shear direction (see Figs. 3 and 4 for a laminate of five aluminium and four prepreg layers), although the aluminium is nearly isotropic whereas in the case of the prepreg the 45 direction is significantly weaker than the 0 and 90 directions. 2.5. Splices
aminations influence load capacity and if they grow. In ANSYS these regions were modelled by two layers of elements with coincident but separate nodes and section definitions to model offsets from the common reference plane. Thus their deformations are independent, Fig. 8. At the boundary of the delamination zones the nodes of one row are connected to the corresponding nodes of the regular region by stiff resp. rigid beams.
In conjunction with Glare splices are used to join different sections of the fuselage. Within spliced regions two metal layers are glued together whereas one prepreg layer steps to the next aluminium sheet (see Fig. 5 lower part) and so forth. Since the skin must be modelled with shells for time and accuracy reasons a layered shell model is used where the diagonal prepreg is replaced by a horizontal one at an average position (see Fig. 5 upper part). 2.6. Delaminations The skin has a laminated construction if consisting of Glare. In some tests the different layers are intentionally not connected to each other in small regions, so-called initial delaminations, see Fig. 6. The question is if these del-
Fig. 7. FE-meshing of delamination zones.
Fig. 6. Intentional delaminations in the skin.
Fig. 8. FE-model of delaminations.
P. Linde et al. / Composite Structures 73 (2006) 229–236
Furthermore, two rows of nodes of the separate regions are also connected by these beams, Fig. 8. The inner row is used to determine interlaminar forces and thus estimate stresses at the crack tip. The coordinates are oriented in such a way that the three components correspond to the three crack modes. The beam elements of the inner row can be removed (kill option). From the following opening (gap) and the former forces the energy release rate can be estimated. These values should lead to a statement about the probability of delamination growth. The progress itself cannot be simulated by this model. However, in physical test the delamination growth appeared just before reaching the limit load. Since these regions represent beginning delamination their diameter is small and in the order of magnitude of the regular mesh size. Inside the circle the mesh size must be small enough to allow for out-of-plane bending. Therefore, a mesh size transition zone in the regular region surrounds the delamination, see Fig. 7. In the numerical simulations the effect on the ultimate load is low because very small relative motion (interference is prevented by contact elements) and thus small influence on the stiffness is considered. This may change when the regions are positioned at inflection points of the buckles. The main information, however, is expected about the probability of growth. For that purpose the criteria are still under development.
Goals of this concept are moment free shear loads, reducing boundary effects and avoiding to influence compressive behaviour by shear load attachments. Same shear loads on all application points are technically provided. The clamps are fixed to a certain skin area (modelled by a rigid target of size and shape of the clamp surfaces with bonded contact), support the edge of the panel in an adjacent region (modelled by standard contact) and get forces at a defined point (modelled as pilot node). For computational time reasons a slightly simplified model is also used as shown in Fig. 10. Since contact between two rigid bodies is not supported the shells perpendicular to the skin are made nearly rigid by a large thickness (Fig. 11). This model allows for the same motions of the skin edge as the more detailed one and is proved to be sufficient for good results. For aircraft-like conditions, the behaviour of the neighbouring panels is to be estimated and here symmetric, antisymmetric or cyclic conditions may be used.
2.7. Boundary conditions For simulation of tests carried out in a physical test rig, it is proved to be of particular importance to realistically model the boundary conditions. Already small deviations in this modelling resulted in large deviations in the results. The test shell is fixed to the test rig firmly along the upper edge, as well as to a moving table at the lower edge. Along the sides, shear forces may be introduced into the shell by means of steel clamps, mounted on the sides of the shell’s edge, supplied with doubler plates, see detail in Fig. 9.
Fig. 9. Detail of boundary conditions in numerical model, reproducing a test rig.
233
Fig. 10. Slightly simplified model of the boundary conditions.
Fig. 11. Detail of simplified boundary conditions.
234
P. Linde et al. / Composite Structures 73 (2006) 229–236
3. Computation 3.1. General The virtual test rig is designed to carry out a numerical computation of a stiffened fuselage panel through the buckling and post-buckling regime until failure. Non-linear material and geometric behaviour is considered at all times during the computation. Arbitrary shear–compression ratios can be chosen. Both static and quasi-static dynamic solutions are feasible, depending on numerical demands of the computed example. In the latter case the choice between implicit dynamic with ANSYS and explicit with LS-DYNA is offered. 3.2. Static solution As default the implicit solver of ANSYS is used for a static analysis. It uses a static load step with the default incremental-iterative solution. Full Newton–Raphson iteration is used. Residual tolerances are kept at the default settings of ANSYS. The load is applied force-controlled. Displacement control is not possible because the ratio of shear and compression forces is prescribed according to the physical testing. In most cases the limit load is approached to a sufficient level with these setting. Sometimes, however, non-convergence appears before signs of global failure are to be seen. Main reasons are sudden changes in the buckling patterns as shown in Figs. 12 and 13 which lead to local instabilities. Tests with arc-length method were not satisfying because of the following reason: ANSYS uses Crisfield’s method which always determines two solutions. Because of the plastic material behaviour loading and unloading paths are not the same. Thus the algorithm did not detect that a decrease in applied load was a simple change in load direction rather than post-failure (Fig. 14).
Fig. 13. Buckling pattern before sudden change.
Fig. 14. Load–displacement-curve in an arc-length method.
3.3. Implicit dynamic solution If there are doubts that the load-level at non-convergence is the ultimate load, the analysis is restarted at a certainly lower level and continued by a transient dynamic solution using the implicit Newmark scheme. Care must be taken that inertia effects remain negligible. Some betadamping is used to increase the time-step size in the autotime-stepping procedure. Small steps must be allowed close to the limit point. In a number of cases this method helped to increase the point of non-convergence until results indicate global failure (Fig. 15). 3.4. Explicit dynamic solution
Fig. 12. Buckling pattern before sudden change.
For implicit solutions encountering divergence prior to physical test failure because of numerical convergence problems, a fast transfer to the LS-DYNA software is possible in the virtual test rig. Because of the explicit time integration similar to the central difference scheme no convergence check is necessary, i.e., no convergence problem can appear.
P. Linde et al. / Composite Structures 73 (2006) 229–236
235
Fig. 17. Load–displacement-curve in an explicit dynamic analysis. Fig. 15. Continuation of non-convergent static solution by transient dynamics.
Another purpose could be to analyse the post-buckling behaviour of the panel (Fig. 16). Firstly this definitely shows that the limit-point is exceeded, secondly one can estimate how sensitively the system will react on imperfections. The latter should influence safety consideration. The computation of the ultimate load is a quasi-static problem, whereas LS-DYNA is a transient dynamic solver. That is why it needs to be carefully chosen in which time
Fig. 16. Post-buckling pattern obtained by LS-DYNA.
the load is applied. Experiences with this method are good, the static equilibrium (same force at the same time at different sections through the panel) is fulfilled until reaching the limit load (see Fig. 17). Further experiences and hints are given in [4]. 4. Results Following computation results may be viewed in the dedicated post-processor of the virtual test rig, which constitutes one of the main items in the selection tree. Particularly it addresses aircraft-specific results. Several of these may be displayed in conjunction with test data, which may be read from file. Results can be divided into global and local results described in the following. Figs. 18 and 19 show the skin deformation (buckling field) of the FE-analysis in comparison to the test results. The scaling of both contour plots is equal; it shows that both the amplitudes and the principal direction of the buckles agree very well. As a comparison a photograph of the buckling pattern from the test is shown in Fig. 20.
Fig. 18. Measured buckling pattern in the physical test.
236
P. Linde et al. / Composite Structures 73 (2006) 229–236
development: expensive and time consuming experimental tests can be partially replaced. The finite element based numerical models of the panels were presented. The implicit (ANSYS) and explicit (LS-DYNA) solvers were discussed. Aircraft specific results obtained in the virtual test rig were presented divided into global and local results, of which several may be presented jointly with test data. 5.2. Conclusions
Fig. 19. Buckling pattern in FE-analysis.
The simulation results of the developed numerical model agree very well with the test results. However, during solutions it was noticed that the structure reacts very sensitively to the modelling of the boundary conditions. The behaviour of the layered shell elements seem to be adequate in the here calculated examples. Concerning delamination growth only certain estimates of the probability can be achieved. The simulation of this effect requires much more effort, but for the specific systems considered here the expected benefit is low because the points of delamination growth and ultimate load are close to each other. 5.3. Recommendations For future research it is recommended to extend the current capabilities of the virtual test rig with the following: • inclusion of further test rigs and panels, • inclusion of further tests, such as frame bending tests, • interactive link with test data.
Acknowledgements
Fig. 20. Buckling pattern in the physical test (Photo).
The authors wish to thank Ju¨rgen Pleitner of Airbus, Hamburg for commenting on the manuscript. References
Further results are available such as various load-deflection curves, local stresses in individual elements, etc. Of particular interest is the formation of interaction diagrams in shear and compression, used in the design of aircraft. 5. Summary, conclusions 5.1. Summary It was shown how a virtual test rig for simulation of stiffened aircraft fuselage panels, consisting of laminates in major parts, was developed. The motivation for the
[1] Linde P. Feasibility study, Virtual structural test analysis system, internal report, Airbus Deutschland GmbH, 2002. [2] Linde P, Pleitner J, Frank C, Gotthold S, Rust W, Schulz A. Development of a virtual test rig for stiffenened fuselage panels. In: Proceedings of the 21st CAD-FEM Users’ Meeting, International Congress on FEM Technology, Potsdam, 12–14 November 2003. Grafing: CAD-FEM GmbH. [3] Linde P, Pleitner J, Rust W. Virtual testing of aircraft fuselage stiffened panels, ICAS 2004. In: Proc 24th International Congress of the Aeronautical Sciences, Yokohama, 29 August–4 September 2004. [4] Rust W, Schweizerhof K. Finite element limit-load analysis of thinwalled structures by ANSYS (implicit) and LS-DYNA (explicit) and in combination. Thin-Walled Struct 2003;41:227–44.