Weight trades in the design of a composite wing box: effect of various design choices

A process to efficiently design composite wing boxes is presented. It uses analytical and semi-empirical equations for failure modes such as material strength, plate buckling, stiffener column buckling and stiffener flange or web crippling. Laminate layups for the different components are selected in accordance with basic engineering rules and guidelines and are updated as necessary to meet the local loads. The emphasis is in allowing buckling of skins at any fraction of the ultimate load and allowing local load redistribution from buckled to non-buckled panels to save weight. The design process is automated and the design can be automatically transferred over to a commercial finite-element code for detailed design and validation. The effects on weight of number of spars, ribs, and stiffeners as well as the fraction of ultimate load at which buckling is allowed are examined and insight is gained to which of these the weight is most sensitive to. In addition, the effect of minimum gage on weight was found to be a driver.


Introduction
The preliminary design of a wing box can become very complex depending on the level of fidelity of the analysis methods used, the multiplicity of failure modes, and the load sharing between components. Consider a design consisting of upper and lower skin (upper skin predominantly under compression and shear, lower skin under tension and shear), spars, stringers, and ribs, as shown in Fig. 1. Any of these components can fail by material failure or stability-driven failure. Stability failure corresponds to global or overall buckling, and local buckling which can be broken into skin buckling including the stiffeners, skin buckling between stiffeners, column buckling of stringers, or crippling of stringers. As the design of one or more of the components is altered, the load in that component changes and load is redistributed away from or attracted to the component. Accounting for this load redistribution is paramount in generating a reliable and robust design. The increased computational capacity and speed in recent decades has helped a lot in increasing the level of detail in the model and the accuracy of the results. However, finding an approach where the level of detail and required computational effort are optimized, such that trade studies and optimization on the entire wing can be performed remain difficult.
One of the earliest works on the design of composite wing boxes was by Watts [1]. In that work, the potential of up to 25% of weight savings and significant life cycle cost savings were identified when compared with the life cycle cost of civil transports of that time. At the same time, the need for improved understanding of the performance of composites in service was mentioned. It was soon recognized that the weight reduction of composite wing boxes brought some challenges such as the manufacturing of post-buckled composite wing boxes which can undergo severe loading without delaminations [2]. More extensive efforts led to large-scale evaluation of composite wing-box designs. Shuart et al. [3] describe the detailed test and analysis effort used to gain insight on the failure mode 1 3 and load of the NASA Advanced Composites Technology (ACT) box beam. Here, the role of structural details such as cutouts, splices, and stiffener terminations was emphasized. A detailed test program to understand failure of impacted wing boxes, still within the ACT program, was undertaken by Jegley and Bush [4]. Here, the importance of impact damage was identified. Failure started in the impacted upper skin and proceeded to the spars. Secondary failures at stiffener terminations were also identified. With the increased understanding of load distribution and failure modes, came attempts for overall design of wing boxes. Schuhmacher et al. [5] combined structural analysis, aeroelasticity, and manufacturing constraints in a large multi-disciplinary optimization model. Kapania and Chun [6] used a simplified formulation combined with empirical relations to obtain the weight of a wing box under aerodynamic loads keeping the maximum angle of twist below a prescribed value. The relative importance of buckling and strength constraints was studied by York [7] with emphasis on accuracy of buckling loads even for non-rectangular skin portions. Load redistribution effects were examined by Ferreira and Almeida [8], where the number of ribs and spars was changed to reach lower weight designs with fewer ribs and more spars. Cost issues were addressed by the same authors [8] and by Klomp-de Boer [9]. The latter describes how reducing the part count and using automation can yield low-cost wing-box components for use on small aircraft. Chintapalli et al. [10] focused on the optimization of stiffened panels in a wing box using buckling for compression and fracture mechanics damage tolerance-based analysis for tension dominated loading. Stamatelos and Labeas [11] considered the effect of strain or stress allowables on the wing-box design and showed that multi-spar wings may lead to lower weight designs. Soloshenko [12] studied how various design concepts affect the stress allowables for a wing design with emphasis on local effects such as the wing root splice joint.
Recent improvements in computational enabled the study of advanced concepts in stiffened panel design such as the use of curvilinear stiffening members [13][14][15]. Furthermore, the use of parallel processing [13,16] or response surfaces [17] has been used to cope with the optimization aspects of individual stiffened panels [13,14,17] or entire wings [15,16]. Topology optimization has been used to determine the location of wing sub-structure under aeroelastic loads [18,19] or minimizing compliance of a wing modelled as truss structure [20]. All this work focuses on designs with buckling, strength, natural frequency, and aeroelasticity constraints, but does not account for crippling of flanges or webs of stiffeners ribs and spars. This is an important failure mode which drives a significant part of the design of compression loaded members [21] and is very rarely included in any of the optimization work in open literature [22]. In addition to frequently being a design driver, this failure mode is very localized as opposed to column or plate buckling and, in a good design, can be tuned to precede more global failure and thus come up with a failure sequence which gives a more robust and damage tolerant structure. Another important issue is the lack of post-buckling in designs of entire wing boxes. While there has been a very large amount of work on individual post-buckled panels, there is very little design optimization work for entire wings using postbuckling as a design driver. It is computationally intensive to do a post-buckling analysis, and for an entire wing, it is difficult to account for load redistribution from buckled panels to adjacent non-buckled ones. Yet, post-buckling is known to save significant amounts of weight (for example, [22][23][24][25]). Recent innovations in use of curvilinear stiffening elements and advanced optimization methods show a lot of promise, but do not account for crippling and post-buckling.

3
It is possible that a post-buckled design, with straight members, will give comparable weight savings with, for example, designs with curvilinear stiffeners and spars. In such a case, it would be preferred, because it has significantly lower cost in manufacturing. In addition, post-buckled designs with straight stiffening members have potential for lower maintenance cost than non-buckled designs with curvilinear stiffening members. Panels with curvilinear sub-structure are harder to repair, because repair patches usually attach to the next available and undamaged "hard" members in the substructure (stiffeners, ribs, and spars). If the available substructure is curved, the repair patch must match this curved shape increasing significantly its fabrication cost and time for installation.
The present communication is a contribution in the direction of finding a reliable, accurate, and efficient approach for preliminary design of composite wing boxes. The objective is to include important failure modes such as crippling and investigate design approaches such as post-buckling which are not often included in the preliminary design of entire wing boxes. This approach leads to very useful conclusions about the effect of different design choices such as the number and type of stringers, number of spars and ribs, choice of layup, buckling load level, etc.
The problem at hand is the determination of the number of spars, ribs, and stringers along with their respective stacking sequence and geometry (for example, stiffener flange widths and web heights), such that the weight is minimized under given loading. Emphasis is placed on efficiency and automation, so configurations can be changed easily making it possible to perform trade studies and optimization. In addition, it is required that a direct automated link with finite elements be possible for immediate validation of results and further design refinement.

Approach
Within the context of a preliminary design, the approach consists of minimizing the weight, while at the same time making sure that there is no failure under the applied loads and that various design guidelines and robust design practices, as selected by the designer, are met. The design focuses on skins, spars, ribs, and stringers. The number and location of spars, ribs, and stringers can change and are user inputs or can be determined as part of the optimization effort at a significant increase in cost. An effort is made to keep the approach simple and efficient. While this will not allow the exhaustive evaluation of more powerful and time-consuming approaches [26], it will quickly allow the creation of an accurate preliminary design.
To check for the various failure modes, the externally applied loads are translated to internal local loads and stresses which are used for strength and stability checks. This requires knowledge of a starting design in terms of layups and internal geometry. The layups for the different components are created and modified with a layup generator that was created for this purpose. Manufacturing and other constraints such as knockdowns for material scatter, environmental effects, and damage are imposed. If a component fails, the layup and/or geometry are adjusted accordingly and local loads are re-calculated to account for load redistribution. Schematically, the approach is shown in Fig. 2. The various steps are described below in some  detail. In Fig. 2, RF is the reserve factor and RF min and RF max are the minimum and maximum value for RF, respectively. It should be emphasized that the focus at this stage is for the creation of a starting design which does not include the effect of high point loads, coming, for example, from a landing gear, or aeroelastic effects. In addition, load cases such as bird strike or impact by runway debris are not explicitly included. As such, the load case used is the static (ultimate) load case used for the certification of the composite DC-10 vertical stabilizer [29,30].
The validation of the proposed preliminary design procedure requires the choice of a realistic aircraft structure. A well-documented aircraft structure is the DC-10 vertical stabilizer [29,30]. The geometry and loads are available in the open literature for the metallic [29] and composite [30] designs. This interesting practical case illustrates the potential gain of using composites as replacement of aluminium with the available technology during the 1980s. Another important point is that horizontal stabilizers do not have landing gear loads and the design is not driven by aeroelasticity. The simplicity of this structure facilitates the interpretation of the trends studied herein. However, it must be emphasized that the fact that these loads are not considered in the present work does not diminish the generality of the proposed preliminary design approach. Landing gear loads and aeroelastic effects may be modelled using simple and yet robust preliminary design techniques compatible with the objectives of the proposed approach.

Determination of internal stresses
The external loads can be forces in any direction applied at any location on the outer surface of the wing box. The wing box is divided into any number of segments (up to a limit set in the input file) in the spanwise direction and the forces and moments are determined at the end of each segment. The equivalent system consists of an axial force, two shear forces, two bending moments and a torque (Fig. 3).
These are then used to determine normal and shear stresses in the cross section. The normal stress is given by with A c the area of the cross section, I xx , I yy , and I xy the moments of inertia of the cross section, and the remaining quantities defined in Fig. 3.
Since the layups differ from one member to the next (skin, stiffener, spar, and rib), Eq. (1) is rewritten to compute strain ε z by multiplying with the corresponding stiffness. Thus, membrane stiffness EA and bending stiffness EI become the quantities of interest: Note that the Young's modulus E appearing in each EI product in Eq. (2) is not the same and depends on the neutral axis location, the axis about which the calculations are done and the (usually different) layup or stacking sequence of segments involved in the calculation.
The shear stress anywhere in the cross section due to the shear forces S x and S y is calculated from the relationship: where t is the thickness at location x, y of the cross section and q so is the constant shear flow corresponding to a pure torque generated by the shear forces S x and S y . This constant shear flow is determined by requiring that the total torque created by the shear flows is the same as the applied torque at the cross section of interest: where q b, the basic shear flow, is the quantity in brackets in Eq. (3) excluding the q so term, p is the normal from the point of intersection of shear forces S x and S y and the tangent to the cross section of the wing box at the point of interest, the integration is taken around the entire cross section, and A e is the enclosed area of the cross section.
The shear stress due to torque applied on the cross section is obtained from For multi-spar designs, the rate of twist is the same for each cell defined by two adjacent spars and the corresponding upper and lower skin and is given by with q the running shear flow, A ec the enclosed area of the cell in question, and A 66 the corresponding entry of the A matrix for the composite laminate at the point in question.

Design guidelines and manufacturing constraints
Several design guidelines and manufacturing constraints are implemented, most of which are at the user's discretion to enforce or relax. From a manufacturing constraint perspective, only 0°, 45°, − 45°, and 90° plies are used for ease of manufacturing. In addition, a minimum allowable thickness is imposed to facilitate handling and avoid damaging of thin skins. This minimum thickness also protects against moisture ingression during service and improves residual strength under impact. From the perspective of a robust design, the laminates are always symmetric and balanced with the + θ and − θ plies close to each other to minimize bending-twisting coupling (small D 16 and D 26 terms in the D matrix of the laminate). In addition, any user-defined knockdown factors for environmental effects, material scatter, and impact damage are included in the allowable stress and strain levels of the materials used in the laminate. In addition, if deemed necessary, a correction on the stiffness of the material can be made in the materials' input file.

Layup generator
For each component, skin, stringer, spar, or rib, the thinnest possible layup is created which satisfies the requirements of Sect. 2.2. This starts with the minimum number of 45, − 45, 0, and 90 plies keeping the layup symmetric and balanced with the 45 and − 45 next to each other to minimize bending-twisting coupling. In addition, if selected by the user, a minimum number of plies, e.g., 10%, in any principal direction is implemented. As a result, the thinnest laminate used is: [45/− 45/0/90]s. Maintaining the same minimum layup everywhere and increasing the number of plies on the basis of local critical direction (see below) allow continuation of a minimum number of plies from one panel to the next (blending). At any location, the current laminate is checked for strength and stability. The lowest reserve factor (RF) is where the multiplier RF mult should be greater than 1. If RF mult is close to 1 (between 1 and 3), the optimizer takes RF threshold = RF min RF mult , 1 3 longer to converge. If it is too large (greater than 10), a lot of un-needed weight is left in the design.
Considering the large number of skin, spar, rib, and stiffener elements, a rigorous optimization to determine the layup would be too time consuming. The approach selected here for the determination of the layups of different members is a compromise which is very efficient and facilitates layup blending between adjacent panels.

Strength checks
The locally applied internal stresses from Sect. 2.1 were used to determine ply-by-ply stresses which were used in a Tsai-Hill failure criterion: where σ x , σ y , and τ xy are ply stresses (x parallel to the fibers) and X, Y, and S are the ply strengths parallel and perpendicular to the fibers and in shear. If the left-hand side of Eq. (8) exceeds 1, there is failure. Then, the RF is calculated as one over the quantity on the left-hand side of this equation.

Stability checks
A number of stability checks were done. These included: Note that if the portion between adjacent stringers is considered, D ij refer to the skin only. If it includes the stringers (panel buckling between spars and ribs, for example), D 11 also includes a term EI/d s, where EI is the bending stiffness of the stringer, and d s is the stringer spacing.
In addition, the expressions shown in Table 1 were used. For stringer column buckling the standard Euler equation is used for the buckling force F crit assuming simply supported ends: with variable a being the length of the stringer and EI the bending stiffness accounting for different stacking sequence for the different members making up the stringer cross section. As the length of the stringer, the portion between ribs is taken with the assumption that there are attachments (e.g., shear ties) between each rib and stringer.
For local buckling or crippling of webs and flanges, the force per unit width N xi applied on each flange or web is calculated using strain compatibility: Buckling of web under bending moment  where N xtot is the total force per unit width acting on the cross section of interest, b i is the long dimension of the flange or web in question, and (a 11 ) i is the 11 entry of the A matrix (membrane stiffness) of the ith flange or web. This applied force is compared to the critical crippling force per unit width: for one-edge-free flanges and for no-edge-free flanges or webs.
"One-edge-" and "no-edge-free" refer to whether one of the edges of the flange is free (not connected to another stabilizing member). In Eqs. (11) and (12), a is the length of the flange, and b is the long dimension of the flange cross section.

Updating the design
Strength requirements are enforced first. If any RF for strength is greater than 1, the local laminate is updated following the procedure, as described in Sect. 2.3. This is done observing all the constraints and guidelines discussed in Sect. 2.2. Buckling and crippling checks are done next with the layup updated accordingly (see Sect. 2.3) for cases, where the RF is greater than 1. An optimized blending approach, maximizing the number of common plies, while minimizing the overall weight would require a detailed extra optimization step [27,28] updating of local layups, an effort is made to enforce some blending between adjacent areas. A minimum number of plies are kept common to all areas and the locally added plies are only 45°, − 45°, 0°, or 90° resulting in a significant number of common plies across different areas.
If a skin panel buckles, its stiffness is reduced by a factor based on the buckling reserve factor. A minimum stiffness multiplier can be set by the user. By reducing the stiffness, the load taken by the buckled panel is reduced and redistributed to other elements which have not yet buckled.

Implementation
The geometry is defined by specifying the root and tip airfoil shapes using the coordinates of 101 points for each with the first and last point coinciding. A sweep angle is defined along with number of spars, number of stringers between adjacent spars and number of ribs. Spars are modelled as shear webs only with the skins taking up the bending loads. In addition, an end cap can be added which closes the wing box and acts as a rib at the tip location. Spar and stringer locations are specified as percentages of chord length. Rib locations are specified as percentages of span length. The ribs can have any orientation specified by the angle of each rib with the root chord. This allows ribs to be perpendicular to a particular spar or parallel to the windstream.
Stringer shapes are defined by a type and a maximum of four parameters which define the length (or angle) of each element in the stringer. The available stringer types are: "Z", "Hat", "L", "C", and "J". These are shown in Fig. 5 and the required parameters for their definition are indicated.
Using linear interpolation, the wing box is divided into a number of spanwise cross sections, the amount of which can be defined by the user. The applied loads for each load case of interest are transformed into a normal force, shear forces, bending moments, and a torsional moment on the corresponding cross section. This system enables the user to define forces and moments acting on any point in 3D space.
Any number of composite (or metallic) materials can be used. A material library has been set up and the layup at each location can combine any number of materials and orientations as long as the design guidelines and constraints of Sect. 2.2 are observed.
The approach is implemented in a Matlab ® program which starts by calculating the properties of each spanwise cross section defined by the user. These properties include centroid and shear centre locations, inertia terms and stiffness terms. Using the cross-sectional properties and applied loads the stresses and strains are calculated for each point on the cross section, as described in Sect. 2.1, while at the same time, a stability check is performed using the methods, as shown in Sect. 2.5. Using the material and layup allowables, reserve factors for strength and stability are calculated and compared to the requirements set by the user. If the laminate at a point on the wing box does not meet the requirements for strength or stability, the layup is adjusted using the layup generator. Once the laminates for the entire wing box are checked and adjusted, the program recalculates the cross-sectional properties and repeats the process until the requirements are met at each point on the wing box. Within the context of the assumptions made, this process leads to an optimized design for the given configuration and design limitations.
In addition to the optimized configuration, which is both in numerical data (dimensions, layups of every geometric point of the wing box), and in plot form, strains, stresses, RF values, local applied loads, local stiffness, and laminate thickness are provided by the software presenting a complete overview of the design.
An important feature of the software is the ability to directly and automatically transfer the geometry created to a finite-element software. In this case, NASTRAN ® was selected. Each geometric point or node in the Matlab ® software becomes a grid point in the finite-element software. Every four non-colinear adjacent grid points are used to create quadrilateral (QUAD4) plate elements. Then, for each element, a corresponding composite layup card (PCOMP) is created using the layup from Matlab ® . Stringers are represented by bar elements with the associated area and moments of inertia calculated with the Matlab ® code. Applied point forces and moments in the Matlab ® code are translated to distributed loads in the finite-element model to avoid unrealistically high stresses and deformation at grids, where high point forces act. These distributed loads are applied on both upper and lower skins of the wing box. Investigation into three different approaches of constraining the root of the wing box (skin nodes fixed, spar nodes fixed, and skin and spar nodes fixed) showed that fixing the spar nodes only gives unrealistically high stresses at the root of the wing, while the other two approaches were more or less equivalent.
For ease of implementation, the first one (skin nodes fixed) was chosen. A typical (four-spar six-rib) NASTRAN ® model contained of 4000 degrees of freedom.
The advantage of this seamless link with NASTRAN ® is the ability to transition immediately from the preliminary design generated by the Matlab ® software to detail design of the down-selected configuration. This would allow tweaking the geometry to reflect more accurate results obtained by the finite-element program and designing additional structural details such as joints, cutouts, doublers, etc.

Validation of the design approach
The automatic transfer feature of the Matlab ® code was used to transfer one of the designs to NASTRAN ® and compare the results. This would give a feel for the accuracy of the preliminary design tool described here and indicate areas, where more attention should be paid during the detail design. As already mentioned, the DC-10 vertical stabilizer was used and the critical static load case that proved and certified the structure [29,30] was imposed. This is a convenient choice, because a lot of the data for the vertical stabilizer are available in open literature and because it was a composite design. While the applied loading, overall geometry, airfoil, and loading were kept the same with the flying hardware, for the purposes of the code validation, different internal configurations (number of spars, stringers, and ribs) were used. All results in NASTRAN ® were obtained with a geometric non-linear static run (SOL 106) which is more accurate than the linear run (SOL 101) when the deflections are large. The applied load was the same as the failure load during the test [29], which was 167% of the design bending moment.
Four different configurations were designed and compared to the NASTRAN results. These range from 2 to 4 spars; 3 or 4 ribs, 0, 3, or 6 stringers; and buckling at 30, 60, or 100% of ultimate load. They are shown in Table 2.
It is important to note that higher numbers of stringers, ribs, or spars are possible and, in some cases, show lower weights. However, the associated cost increase due to increased tooling complexity and increased assembly time was considered prohibitive and results with higher number of stiffening elements are not reported here.

3
The results of the comparison are shown as plots of the difference between Matlab ® and NASTRAN ® in Fig. 6. Results are shown for the first three cases of Table 2 (the last one is very similar to case 3). Green shows perfect agreement, red shows maximum difference when the Matlab ® software is higher than NASTRAN ® , and black shows maximum difference when it is lower than NASTRAN ® . In these figures, the upper skin has been removed for clarity.
It can be seen from Fig. 6 that the stresses along the 0° direction (set of plots to the left of the figure) are in goodto-excellent agreement with NASTRAN ® with the exception of the root region and localized skin regions near the tip, where differences are larger. This can be put in perspective if the maximum difference in each case is compared to the maximum stress in the wing box. For the first case (Fig. 6a), the maximum difference of 50 MPa is 24% of the maximum value. For the second case (Fig. 6b), the maximum difference of 50 MPa is 18% of the maximum value. For the third case, (Fig. 6c), the maximum difference of 20 MPa is 49%. Similarly for the shear stresses, the corresponding percentages are 24, 20, and 33%. The differences are attributed to a number of reasons: (a) boundary condition differences between Matlab ® and NASTRAN ® which would explain why the differences are maximized near the root of the wing box; (b) simplifications in the analytical approach used here, Eqs. (1)-(6), while the NASTRAN ® model has a more accurate distribution of stresses; this includes neglecting warping effects at the root of the wing box; and (c) non-linear effects present in the NASTRAN ® model due to moderate deflections which are not accounted for in the present approach. Nevertheless, the differences are confined to a small area of the wing box and/ or are at regions, where the stresses are only a small fraction of the peak stresses (see shear stresses with maximum difference in Fig. 6) suggesting that the present approach is quite reliable and accurate over the most part of the wing box. Separate finite-element runs for different cases using NASTRAN ® reproduced exactly the same trends in terms of relative position of the different designs and showed less than 5% changes in the maximum stresses and the absolute weights compared to the results of the Matlab ® runs. The present method, therefore, provides a good starting point for further refined analysis. In any case, the present approach is not expected and is not meant to replace detailed finiteelement analysis. It is expected to provide a good preliminary design reasonably optimized over thousands of potential candidates and thus "jump-starts" the finite-element modelling effort. Considering it only takes 5-10 min on a laptop (Intel i7, 2.2 GHz core [20% loaded by Matlab ® ]; 8 GB RAM [500 MB used by Matlab ® ]) for the complete preliminary design, this is a promising approach allowing quick trades of number of spars, stringers, and ribs.

Case study: application to the DC-10 vertical stabilizer
The approach described above was applied to the DC-10 vertical stabilizer [29]. The load case corresponds to maximum bending which is the critical ultimate load requirement. This is a good example first because of openly available data and second because it went through a weight reduction study transitioning from a metal design at 456 kg to a composite design at 352.7 kg. If fittings, access doors, and miscellaneous items not included in the present study are excluded from the DC-10 composite wing box, the result is 275.8 kg [29]. Depending on the number of spars, ribs, and stiffeners and the percentage of ultimate load at which buckling of skins is permitted, the present approach yields a composite weight in the range 244.5-272.6 kg (see Sect. 5.4 and Table 3 below). It should also be noted that the composite design has gone through detail design which, typically, during the last stages, reduces weight below the one achieved during preliminary design. The present approach does not include the detail design step and may not capitalize on additional savings from the detail design phase. Such a comparison, however, does show whether the present approach is in the right ball park and it will show if it is promising when  applied to existing hardware. Furthermore, given that the baseline composite design [29] is a sandwich structure, the comparison with the results of this approach which is a skinstiffened design will give valuable insight into differences between sandwich and skin-stiffened designs.
The results in this section used a combination of knockdown factors (see Sect. 2.2) leading to strength values equal to 42% of the mean, undamaged room temperature strength values.
The applied load, maximum moment distribution was the one corresponding to the failure load during test which was 167% of limit load. Two types of stringers were considered, "Z" and "Hat" with the geometry shown in Fig. 7.
The number of spars and ribs was varied from 2 to 6. The number of stringers was varied from zero to a maximum, which is a function of the minimum allowable spacing, set at 150 mm at the tip of the wing box. This led to a maximum number of six stringers per skin when four spars were used. The material used in the design was AS4/3501-6. Results are presented for the effect of buckling ratio, number of ribs, and number of spars. In all cases, the effect of varying the number of stringers is also shown.

Effect of buckling ratio
The effect of the percentage of the ultimate load at which skin buckling is allowed is shown in Fig. 8 for different number of Z stringers, from no stringers to four stringers between spars. In this case, two spars and two ribs (one at the root and one near the tip) are used.
As might be expected, reducing the buckling ratio, reduces the weight. The lower buckling ratio reduces the bending stiffness required for the skin which leads to lower skin thickness. This weight reduction is gradual up to buckling at 60% of the design ultimate load and then is steeper. There even appears to be a minimum for designs buckling at around 30% of ultimate load for cases with stringers, however the increase in weight for lower buckling ratios may be a combination of enforced guidelines (minimum thickness and 10% rule for layup). It is important to note that buckling ratios below 50% of ultimate are very challenging during detail design because of large deflections, additional failure modes such as skin/stiffener separation, and fatigue concerns. Therefore, buckling ratios lower than 50% of ultimate should be chosen with care and only if the extra weight savings potential justifies the risk and the extra effort. Similar trends of decreasing weight with decreasing buckling ratio were observed for other configurations with multiple spars and ribs.

Effect of number of ribs
The effect of adding ribs is shown in Fig. 9. Buckling is enforced to occur at 60% of the design ultimate load and two spars are used. Two curves are shown in Fig. 9, one for the case of no stringers and one for two stringers. Going from two to three ribs adds a small amount of weight for all cases. This was observed also for cases with one and three stringers. The reason is that the weight saved by increasing the buckling load of the skin by going to lower aspect ratio designs is less than the weight added by the rib because of the minimum thickness requirement (minimum set to 2 mm). The intermediate rib is lightly loaded and could be designed with lower thickness. Further increasing the number of ribs beyond three reduces the weight.

Effect of number of spars
The effect of increasing the number of spars is shown in Fig. 10 for two-and three-rib designs with no stringers. Again, buckling is enforced at 60% of ultimate load.
As can be seen from Fig. 10, increasing the number of spars decreases the weight monotonically. The drop is significant at the beginning, up to four spars and less drastic for five and six spars. This trend has been observed by other authors also [8,11]. The present method allows quantifying this trend when the skins are allowed to buckle below ultimate load at any preselected percentage. As was observed in Fig. 9, increasing the number of ribs from 2 to 3 increases the weight slightly. This is shown by the 3-rib curve being above the 2-rib curve in Fig. 10. In addition to the effect of the number of stringers shown in Fig. 8, this is an alternative design approach. A multi-spar multi-stringer design promises to yield the lowest weight. It should be kept in mind, however, that the cost of such a design should be considered carefully because of the large part count and the complexity of the associated tooling.

Best design and comparison to baseline
Based on the trends discussed in Sects. 5.1-5.3, for the lowest weight configuration, a design with low buckling ratio and a relatively high number of stringers would give the lowest weight. It is not clear whether a two-spar or a multispar design would be more efficient. Increasing the number of spars unloads the skins which makes achieving lower buckling ratios easier. At the same time, each added spar adds weight and it is not always guaranteed that the weight removed from the skins is greater than the weight added by the increased number of spars. For this reason, a number of results varying also the number of ribs and types of stringers are shown in Table 3. This is done for two different buckling ratios, 30% of the design ultimate load, which gives the lowest weight, and 60% of ultimate load, which is considered more realistic. The lowest weight of 244.5 kg corresponds to a twospar, two-rib design with four Z stringers between spars buckling at 30% of the design ultimate load (45% of design limit load). As mentioned earlier, such a buckling ratio is very aggressive considering the relatively large number of manoeuvers during which the skins would buckle leading potentially to skin-stiffener separation concerns and fatigue   issues. Therefore, buckling at 60% of ultimate load (90% of limit load) is considered more realistic and achievable. In this case, the lowest weight is 272.6 kg corresponding to a two-spar, three-rib, six-Z stringer design. There is, however, a fairly wide variety of slightly sub-optimal designs shown in the second half of Table 3. Within 8.7% of the 272.6 kg there are six other designs ranging from two to four spars, two to six ribs, and three to six-Z or hat stiffeners. This can be very useful as near-optimum alternatives can be selected based on design preferences and/or cost issues. For example, designs with fewer spars and stringers are expected to have lower cost. The tool developed here allows determination of all these near-optimum designs and gives a good weight estimate for each alternative.
The baseline composite design [29] had sandwich skins and weighed 275.8 kg. The best current design obtained here at 272.6 kg is lighter, by 1.2%. Or, to be consistent with the no-buckling requirement imposed in the baseline, the lowest weight obtained here is 277.1 kg which is half a percent heavier than the baseline. This is preliminary design and refinements based on more accurate simulations using detailed finite elements and further optimization, for example, not limiting the layups on combinations of 45, − 45, 0, and 90 plies, are expected to further improve the weight. Furthermore, allowing buckling earlier, as already mentioned, leads to weights (244.5 kg in Table 3) as much as 11.3% lower. The thickness distribution for the optimum design of 272.6 kg is shown in Fig. 11.
Additional insight into the optimum design is gained by studying Fig. 12, where the reserve factors for buckling are shown. Regions with high reserve factor, such as the two spars near the root and the rib closest to the root, suggest that the design there is not driven by buckling but by strength requirements.

Equations for weight trends
The efficiency of the method allows multiple runs to examine in detail the effect of changing the number of stringers, spars, and ribs. For the case of stringers, see, for example, Fig. 13, the weight savings as a function of number of stringers is found to be given by where ΔW 1 is the percent weight savings in going from no stringers to one stringer, N is the number of stringers and the right-hand side of Eq. (12) is summed as a geometric series. The different points in Fig. 13 at each number of stringers correspond to different buckling ratios. Equation (12) is found to be within less than 1% error in predicting the best fit curve, as shown in Fig. 13. Repeating for other spar and rib configurations shows that Eq. (12) is still less than 4% off from the best fit curve. It should be noted that for multi-spar designs, the highest accuracy is obtained when averaging the results obtained from Eq. (12) when N is the number of stringers per cell, and from Eq. (12), again, but now, with N, the total number of stringers across all cells. The accuracy of Eq. (12) has been verified over a large number of cases. However, it still is a simple "curve-fit", and as such, it is not expected to always apply. For one thing, different load cases including aeroelastic effects, not accounted for here, may lead to results violating this equation. Its applicability should be checked on a case-by-case basis. It is very useful for trending and performing trade studies. Unlike the trends on stiffeners, where the effect on the weight of each additional stiffener is non-linear, the effect of adding ribs is linear, but the accuracy of the linear trend is not as high as with stringers, especially for the cases of two and three ribs. For four ribs or more, for a two-spar design, every time a rib is added, 2.55% of the wing-box weight is saved, while for three-and four-spar designs, adding a rib increases (as opposed to decreasing) the weight by 0.67 and 0.39%, respectively. The reason for this non-monotonic behaviour is a combination of competing factors. Adding a rib increases the weight, but at the same time changes the aspect ratio (length along span/width) of the skin and spar panels. Depending on the number of spars used, this change in aspect ratio may be an overall reduction, leading to a decrease in buckling load and thus weight increase in the skin or spar, or an overall increase, which, for low enough aspect ratios, increases the buckling load meaning the skin and/or spar thickness can be reduced thus saving weight.
The effect of adding spars is very similar to that of adding stringers shown in Fig. 13 and Eq. (12). Equation (12) holds up to five-spar designs with two or three ribs.

Summary and conclusions
An automated and efficient approach to perform preliminary design of composite wing-box structures has been presented. It accounts for strength, buckling, and local buckling failures and incorporates a variety of design guidelines and practices. A layup generator is used to determine the layup of skins, spars, and stringers at each location adding or subtracting plies as necessary to avoid failure and minimize the safety margin. The approach is fully automated in Matlab ® and is directly linked to NASTRAN ® for further fine-tuning and detail design. The number of spars, ribs, and stringers can be varied at will and the speed of calculations allows extensive trade studies and optimization. The results of the method were compared to NASTRAN ® and found to be within 20% 1 3 for the most part with higher deviations near the root of the wing box. The biggest weight reductions are achieved by reducing the buckling ratio (skin buckling load/ultimate load) and increasing the number of stringers. Buckling ratios as low as 30% (skin buckling at 45% of limit load) were used. It was also found that multi-spar designs can be quite efficient if combined with the optimum number of stiffeners. Z stringers were as weight-efficient as hat stringers. Increasing the number of ribs beyond three reduces the weight but by a relatively small amount. Simple and accurate equations for the weight reduction as the number of stringers, spars, or ribs changes were also determined.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creat iveco mmons .org/licen ses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.