Cost-oriented process optimisation through variation propagation management for aircraft wing spar assembly

Abstract Overconstrained assemblies such as aircraft sub-assemblies present a challenge to production planners, as variations in parts and processes can make it difficult to achieve all assembly Key Characteristics (KCs) simultaneously. Despite assigning tight tolerances to sub-component manufacture, part variation propagation necessitates expensive and time-consuming variation management processes such as shimming in order to ensure the final assembly is within specification. This paper presents for the first time a variation propagation model for overconstrained assemblies, and develops a novel modelling method to connect variations with production costs. This facilitates a novel process optimisation method based on variation propagation, with the ability to analyse the trade-offs between the cost and achievable variation limits of the entire manufacturing chain in order to minimise the overall manufacturing cost. An overconstrained wing spar assembly is used as a case study to validate the methodology.


Introduction
In the highly competitive aerospace industry, cost reduction whilst ensuring high quality standards, has become an important method of gaining an edge in the market [1]. One area with cost reduction potential is the production of large overconstrained assemblies such as wing spar assemblies. It is commonplace for assemblies to be deliberately designed to be overconstrained by their parts in order to use locked-in stress to achieve the desired gross shape, for instance providing shape definition for flexible parts [2]. In other applications, overconstraint is beneficial in assemblies in order to withstand external forces when in operation, and is ensured by having more locators than needed to provide location [2]. In both cases, overconstrained assemblies present a difficult task to design engineers and process planners, as variations in sub-components make it problematic to assemble a final product that meets specification. This is because the achievement of each assembly Key Characteristic [3] (KC) is not independent in overconstrained assemblies. Therefore, improving one KC can degrade another KC, forcing rework or scrappage at the assembly stage [2]. To manage this variation, designers often assign high tolerance requirement in part fabrication, and process planners place stringent precision criteria in every manufacturing and assembly operation. Nevertheless, additional variation management processes such as shimming are still required at the assembly stage, which are expensive and time consuming [4,5].
Production systems are typically not optimised for variation management processes at assembly despite overconstrained assemblies often necessitating their inclusion. It is typical for manufacturers to outsource part manufacture to their supply chain [6], with only final assembly being completed in-house. This can result in a disjoint between the two processes, whereby part manufacturing tolerances are not appropriate for the assembly tolerances. However, cost reductions can only be confidently attained when the production chain is considered holistically. The cost for part manufacture and assembly may be estimated and minimised separately, but without a feedback loop between the various stages of production [6], the inter-relations between part manufacture and assembly are neglected, and the overall time and cost of assembly may therefore not be as low as is possible. For example, high precision requirements generally leads to high cost [7], but it is arguable that not all the components and processes should require high precision simultaneously. Relaxing the precision requirements for certain components and processes may not undermine the final precision requirements, but will result in a lower total cost of production. Therefore, it is important to provide a solution to the allocation of precision requirements in order to optimise the overall manufacturing process. This drives the need for a variation model that can map variation propagation from manufacturing through to assembly. field of robotics accuracy in the 1980s [8]. Since then, the area of variation modelling has grown to become a highly active research field, particularly with the development of the stream of variation model [9][10][11][12], which models variation transmission between stages in multistage production systems. Assembly variation modelling has been the subject of a large amount of research [13]. It has established methods for modelling connective assemblies [14,15], and has facilitated the analysis of assemblies, where parts are assembled through part-to-part mating surfaces [16,10,11]. There are a number of key challenges associated with variation, such as the large number of potential variation sources, disparate variation magnitudes, the stochastic nature of variation occurrence, and a lack of knowledge of variation propagation mechanisms [17,18]. The most recent research in the area of variation propagation has developed physics-driven variation models [19] and compliance-induced variation models [20][21][22]. To date, most assembly variation models have focused on the connective assembly relationships between parts. Overconstrained assemblies, however, present a unique challenge, as part relationships can be spatial instead of connective [23]. Cai et al. [24] noted this as a barrier to assembly variation modelling, and developed a variation propagation model for spatial assemblies of two parts with one part-to-part mating interface [24]. Large overconstrained assemblies typically feature many parts, with multiple part-to-part mating interfaces, which often cannot be achieved simultaneously with variations present. Recent work has derived how the variation level at KC features in these assemblies can be quantified [18]. This paper will develop this research further by relating the variation level at KC features to the cost of assembly in order to determine the most appropriate variation management process that will reduce the overall cost of production.
Each variation management process has a different associated cost, process capability and downstream process requirement. With so many possible processes available to production planners, determining the most appropriate process for a given production system is challenging. Process options range from manual shimming [4], to state-of-the-art processes such as in-process adjustment [25], adaptive assembly [26,27], compensating for deviations by machining [4,28,29] and predictive shimming operations [6]. Making informed decisions, therefore, requires an understanding of the relationship between variation and cost implications of each process in the production chain. Relationships between cost and variation are potentially nonlinear [30] and complex [31], and optimal production chain solutions are likely to be a compromise between manufacturing and assembly process cost.
Research has related cost to tolerance for one process or machine [32,33], but less focus has been given to the summation of costs from many processes in a production chain [30,34]. The cost impact of variation management processes has been only scarcely investigated [7]. This paper addresses these gaps in the literature, and proposes a methodology that relates cost and variation for an entire production chain, from part fabrication to final assembly. This paper is organised as follows. Section 2 introduces the overconstrained aerospace wing spar assembly which is used to demonstrate the proposed methods. Section 3 presents the variation representation and modelling method used throughout this paper. Section 4 deals with part manufacturing variation sources and assembly process variations. Section 5 presents a novel variation propagation modelling method for overconstrained assemblies. Section 6 describes how a cost model can be related to the variation model. Section 7 presents details on how the developed methodology can be used for manufacturing system analysis and optimisation. Section 8 concludes this paper with the key findings.

Aircraft wing spar assembly
The methodology to relate variation and cost proposed in this paper is illustrated by an overconstrained aerospace case study. Specifically, the case of a wing section is considered, consisting of a structural spar and four hinge brackets [4], as shown in Fig. 1. There are two mating feet per hinge bracket, resulting in the spar being mounted at eight surfaces along its span, thus making it overconstrained. There are three KCs which must be achieved in this assembly, as highlighted in Fig. 1. KC features are those that are critical for performance, safety, or regulations and are at risk of not being achieved due to variations [35]. The first KC is the alignment of the four hinge brackets to each other along the hinge bracket datum hole x-axis. This is to ensure the correct operation of the aircraft control surfaces, which are mounted to the hinge brackets. The second KC that must be satisfied is the alignment of each of the four hinge brackets to the spar surface. This is to make certain that the assembly fits within the wing's aerofoil profile. Finally, the third KC that must be achieved is the mating of the hinge bracket feet to the spar surface, without gaps or clashes occurring. This is to ensure the structural integrity of the overall assembly.
Variations induced during part manufacture of the spar and hinge brackets, and variation induced during the assembly processes lead to gaps and clashes during assembly, resulting in the failure to achieve all KCs without additional variation management processes. These variation sources and effects will be presented in the following sections.

Variation representation and modelling
A feature-based variation modelling method is developed in this work. Parts consist of numerous features, while assemblies consist of a number of parts and assembly features. The part features significant for variation propagation modelling are classified into one of two categories: mating features or datum features. Mating features are contact features to be mated together between two parts, whereas datum features are the contact features used as a reference for assembly processes between a part and a fixture. Products may be assembled with or without the help of fixtures. Fixtures also consist of features. A standard naming convention, as proposed by Cai et al. [24], is used in this paper to define parts, features and fixtures of the spar and hinge bracket assembly, as shown below. The postures of features are represented using coordinate frames. Each mating or datum feature on an individual part is given a coordinate frame, as is each feature on the fixture. A global coordinate frame, G, is also defined. Coordinate frames are related to each other in 3D space using homogeneous transformation matrices (HTM) [36]. HTMs define the relative posture and position of each coordinate frame to each other. Each HTM consists of a 3 × 3 matrix, R, describing the orientation of the frame, and a 3 × 1 vector, P, which defines the position of the coordinate frame. Eq.
, & are small. The orientation and position transforms can be obtained from a part's CAD data, from physical measurement, or assumed if the part is at an early planning stage. The X, Y and Z values define the feature's position relative to another feature, and the θx, θy and θz values, which are XYZ Euler angles, define the feature's orientation.
Variations in the assembly are represented by deviations in the coordinate systems of part/assembly/fixture features. For example, variations of P 1 relative to the global frame, (3)

Hinge bracket fabrication variation sources
The hinge bracket is manufactured by machining a billet of aluminium to the net hinge bracket shape. Variations are inevitable in the fabrication process of the hinge bracket. The key features of the part are the two feet, which mate with the spar, and two datum features, as shown in Fig. 2. The two datum features, DFeia and DFeib, locate the part in the fixture. DFeia defines the X, Y and Z location and rotation about the Y-axis and Z-axis of the hinge bracket in the assembly. DFeib defines the rotation about the X-axis, θx, of the hinge bracket. MFeia and MFeib represent the hinge bracket feet features which interface with the spar. Each of these features are assigned coordinate frames, allowing an ideal model of the ith hinge bracket, PHi, to be produced by linking the feature coordinate frames using HTMs.
To fabricate the part, first an aluminium billet is loaded into a machine centre, and is located using a 3-2-1 locating scheme. There may be discrepancies in the condition of work piece, such as defects in raw billet material and material property variations [15,37], which may cause variation in the loading of the billet into the machining centre. Cutting tools, the fixture used, the operator [15,37], the order of processes and the environment can also all lead to variations in the final part. In this case study, a roughing cut is performed on the billet to define the shape. Geometric feature variations (variations in shape of component features, such as orientation and location of features) [15] can occur during this manufacturing process due to machine tool elastic and thermal deformation [38,39], cutting tool wear [15], and machine vibrations [39]. A finishing cut can then be completed, which has greater control over the geometric accuracy of the cut, thus typically producing a final feature with less variation and a better surface finish. Finally, a datum hole is drilled. Variation can be induced in this process due to drill wear [38,40], incorrect calibration of the drill, and vibrations in the drill [39]. The overall dimensions for the hinge bracket part are shown in the engineering drawing in Fig. 3. Table 1 collates the variations which are inputted into the hinge bracket ideal model. The tolerances of typical machining processes [41] and engineering drawings of existing hinge brackets were used to define the likely variation input at each feature, taking into account the potential variation sources described above.
The variations described in Table 1 are combined into δ Phi (i = 1,2,3,4), a vector representing all variations associated with the hinge bracket.

Spar fabrication variation sources
The accuracy of the finished spar depends on the manufacturing process used. The coordinate frames of the key mating and datum features of the spar for assembly are presented in Fig. 4. There are two datum hole features, i.e. the king pin datum hole, DFe1a, and the queen pin datum hole, DFe1b. The two datums position the spar on the assembly fixture. DFe1a, as the primary datum, is located first on the fixture, and therefore defines the X, Y and Z location as well as the rotation about the X-axis and Z-axis of the hinge bracket in the assembly. DFe1b is located second, and defines the rotation about the Yaxis, θy, of the spar. There are four additional datum features, DFe1c, DFe1d, DFe1e and DFe1f, which align the four hinge brackets respectively to the spar during assembly. Finally, there are eight mating features, MFe1a, MFe1b, MFe1c, MFe1d, MFe1e, MFe1f, MFe1g and MFe1h, where the spar joins with the four hinge brackets.
In this paper, two fabrication methods are considered, i.e. manual layup of pre-impregnated Carbon Fibre Reinforced Plastic (prepreg CFRP), and Resin Transfer Moulding (RTM). All the CFRP features, including the hinge bracket alignment datums, will have some degree of variation due to the CFRP fabrication process [42]. Variation in CFRP manufacture comes from a variety of sources, including ply alignment tolerance, heterogeneous material properties, and undulation effects as additional plies are added [43]. These variation sources have been classified by Perner et al. [42]. Furthermore, due to the size of the spar component, and the fact that it is supported at each end by the assembly fixture, features in the centre of the component may suffer variation as a result of compliance [31]. This manifests itself in the form of sagging in the middle of the spar. Finally, the spar alignment datum features may also exhibit variation due to spring-out when it is removed from the manufacturing mould. Although both RTM and prepreg methods utilise a mould to lay up CFRP, each induces characteristic variation unique to that process. The accuracy of the final part is highly influenced by the accuracy and location of the mould. Typically, the CFRP laid on the mould, known as the tooled surface, will have a greater accuracy. The RTM process is generally considered to be more accurate, as all surfaces of the part are tooled, whereas with the prepreg method only one side is tooled. Two holes are also drilled in the spar in order to locate it onto the assembly fixture. Variation can be induced in this drilling process due to CFRP delamination, tool wear [38,40], and vibrations [39]. Operator setup, sequence and tooling can influence drill angle, and can induce angular variation into the feature. The dimensions of the wing spar are shown in the engineering drawing in Fig. 5. Table 2 collates the variations, which are inputted into the hinge bracket ideal model. The tolerances of typical CFRP processes [44] were used to define the likely variation input at each feature, taking into account the potential variation sources described above.
The variations described in Table 2 are combined into δ PS , a vector representing all variations associated with the spar.

Assembly process
The assembled spar and hinge brackets must be assembled spatially   Fig. 6. This is achieved through the use of an assembly fixture. The assembly can be classified as a type-II assembly, as part mating features are not defined until the assembly process [25]. A major difference between overconstrained assemblies and constrained assemblies is that overconstrained assemblies have an ambiguous KC delivery chain, i.e., there is more than one way to assemble the final product. Different KC delivery chains will result in different magnitudes of assembly variation, occurring at different locations in the assembly. Therefore, in order to model variation propagation in overconstrained assemblies, the order in which KCs will be achieved during assembly must be determined prior to modelling. In the spar and hinge bracket case considered in this paper, the current assembly technique used in industry is modelled. The spar and hinge bracket case is described in Fig. 7.
The KC delivery chain is a consequence of the assembly technique. In assembly step 1, KC1 is achieved. In the next assembly step, KC2 is then achieved, without diminishing KC1. However, although KC1 and KC2 can be achieved simultaneously, it has been found that KC3 cannot be achieved, as gaps or clashes between the spar and the hinge bracket feet occur, as shown in Fig. 8.
As a result, variations not only come from the spar or hinge bracket fabrication but also from the assembly processes or fixtures. Assembly process variation sources include operator reach, preference, stature and positioning [45], the sequence of operations, misalignment of parts [46], fit-up problems [15], as well as variations in fastening forces and clamping forces [28]. Fixture variation sources include profile variations in locating surfaces [16,47], clamping-force induced variations [48], variations in locating parts on fixtures [40], and the misassembly of fixtures [46]. Values for these variations were assumed to be an order of magnitude smaller than part variation, and are displayed in Table 3.

Variation propagation modelling for overconstrained assemblies
To model the variation propagation which occurs during the assembly process, the spar, hinge bracket and fixture features are mapped, and represented by coordinate frames. Fixture features are also referenced using the standard naming convention. A complete network of features in an assembly is generated by creating chains of HTMs. Chains are significant for defining how parts are mated together. Depending on the final product architecture and assembly technique, there are open chains or closed chains. A chain of transforms, is generated by multiplying HTMs. When all features are connected in this way, the model will contain the spatial relationships between features within the same part, and between different parts in the same assembly.  KC1, KC2 and KC3, shown in Fig. 1, can be translated into a liaison diagram as shown in Fig. 9. The liaison diagram displays all the relationships between part features and fixture features within the assembly. An example of the numerous variation propagation routes (routes 1, 2 and 3) for one of the mating interface, KC3a, is also presented.
The assembly variation witnessed at the KC3 mating interface is due in part to manufacturing variation induced during spar or hinge bracket fabrication. Some of the assembly variation will be due to fixture variation of the spar and hinge bracket locators. Furthermore, variation is induced when the hinge bracket is rotated to contact the index in order to satisfy KC2. Whereas manufacturing variation can be readily measured during inspection, and controlled through tolerance allocation, the variation induced during assembly is not as intuitive to understand, and requires a mathematical model to quantify at the planning stage. A variation propagation model was therefore developed as shown in Fig. 10. Fig. 10 can be converted into chains of HTMs for each feature. By relating both the spar and hinge bracket mating features to the global coordinate system, G, the maximum gaps or clashes between the mating features at KC3 can be calculated. For example, the spar mating feature, MFe1a, can be related to a global coordinate system, G, by using Eq. (4).
When variations are considered, the chain of HTMs to represent the spar mating feature, MFe1a, is as shown in Eq. (5). This corresponds to variation propagation route 1 in Fig. 9.
In Eq. (5) corresponds to assembly process variation, such as misalignment of parts.
The hinge bracket mating feature can also be related to the same global coordinate system G. The hinge bracket is located by two datums, DFeia and DFeib. Therefore, two equations are used to define the location and orientation of the hinge bracket. For example, the first hinge bracket, PH1, is first located by the datum hole feature, DFe1a, by using Eq. (6).
When variation is considered, the variation model of the hinge bracket shown in Eq. (7). This is equivalent to variation propagation route 2 in Fig. 9.
The hinge bracket is then aligned to the spar using the fixture index blade, Fi1c. This defines the orientation of the hinge bracket. Furthermore, it allows the location of the mating feature, MFe1a, to be determined, as shown in the ideal model in Eq. (8).
Variation must also be considered in the chain of HTMs, and can be introduced as shown in Eq. (9). This corresponds to variation propagation route 3 in Fig. 7.
By relating both features to the same global coordinate system, the spatial relationship between the two mating surfaces can be determined by using Eq. (10).
Eq. (10) measures the gap or clash between the two mating surfaces and corresponds to KC3a from Fig. 9. δ KC3a includes all variation sources, propagated from part and process variation to fixture and assembly variation.

Manufacturing cost modelling on parts and assemblies associated with variations
Cost and variation are intrinsically linked to process selection. The following section describes a method of associating variation with cost for the overconstrained assembly, with process capability being the key link between the two models. Activity Based Costing (ABC) is utilised, because the cost of variation is primarily linked to activity drivers, and not resource drivers or costs drivers [7]. In practice, this means that the effect of variation on cost is most influenced by the activities required to complete production. For example, additional activities added to the process chain at the assembly stage will have a significant effect on variation and cost. In contrast, the consumption of resources during a given activity has been considered to be constant in literature [7]. Other costing methods such as parametric costing are not explicitly associated with each activity, and are therefore less suitable for this application. The ABC method has been previously applied to both tolerance allocation [34] and quality and cost optimisation [49].
Manufacturers have utilised the selective assembly method to group parts in order to minimise the effects of assembly and improve the quality of their products in a cost effective manner [50]. But this method is only applicable to medium to large volume production, and may not be suitable for small volume production such as aerospace manufacturing. Furthermore, the cost of WIP would in particular be prohibitive for such large and high value components.
To establish the cost model, the production chain was first divided into all its constituent operations. Each operation was then assigned a standard process time. The time was then multiplied by an operation wrap rate in order to estimate the operation cost. Wrap rate is also known as the fully burdened labour rate, and typically includes direct labour wage rate, overhead costs rate and other costs firms incur, such as general and administrative costs. Material cost was calculated based on the raw material volume and its unit cost. The summation of the operation cost and material cost for all operations provides the total cost for each process.

Cost modelling of part fabrication considering variations
The cost of part fabrication is dealt with first in this paper. The geometric variations of each part will lead to cost variations because of different manufacturing processes. Two methods of hinge bracket fabrication and two methods of spar fabrication are considered. In each case, one high accuracy process is considered, alongside one lower accuracy process. The fabrication process used for hinge bracket fabrication affects the variations induced into the hinge bracket, and thus affects the values of Table 1. A Monte Carlo simulation was developed in order to determine the magnitude of the assembly gap or clash. Similarly, the method of spar fabrication affects the range of variations from Table 2 which will be inputted into the Monte Carlo simulation. Again, this will have an influence on the magnitude of assembly variation, and thus will affect the cost of assembly.

Hinge bracket fabrication cost considering variations
Two methods of hinge bracket fabrication were considered. The first method considered was a high accuracy method. The hinge bracket was machined from a solid billet, using both a roughing cut and a finishing cut. The cost of the high accuracy process was estimated by using Eq. (11).    Assembly positioning variation, Fi1a-j ± 0.0100 ± 0.0080 Fixture geometric variation, Fi1a-j ± 0.0100 ± 0.0080 administrative costs. The second method considered was a lower accuracy method that omits the finishing cut step of the first method. This reduces the cost of the hinge bracket fabrication, but also reduces the quality of the finished part, i.e., it increases the variations associated with the part features.

Spar fabrication cost considering variations
Two distinctly different processes were considered for spar fabrication in this paper: prepreg and RTM. Both have different total spar fabrication cost, C PS , with differing material, labour, machine, and overhead costs. The manual prepreg method cost is estimated by using Eq. (12).   . The second spar fabrication process considered, RTM, has a lower material cost, as dry CFRP is used instead of prepreg. However, in this paper it is estimated to be more expensive overall, as it requires more tooling and a higher capital expenditure on equipment. The RTM method cost is estimated using HYPERLINK Eq. (13).

Part fabrication cost summary
As a result of considering two methods of hinge bracket fabrication, and two methods of spar fabrication, there are four possible combinations of processes. These are numbered Scenario 1-4, as shown in Table 4. The costs of each of the scenarios were calculated by using Eqs. (11)- (13). The values of these variables are listed in Appendix A.

Cost modelling of assembly processes considering variations
The cost of assembly C AS , is related to the assembly variation level δ KC3 . Whether a gap or a clash exists at assembly KC features depends on the variation propagation from all the possible variation sources. In this paper, two methods of variation management process are considered: shimming-based assembly and milling-based assembly.

Cost modelling of shimming-based assembly
The first assembly variation management process considered in this paper is a manual shimming process. This is the removal of gaps at assembly KC mating pairs, by using either polymer or aluminium shims. The cost of this process is dominated by the labour cost. However, the material used for shimming is also a major cost driver, and it is directly related to variation in the assembly. As regards for small gaps (≤0.2 mm), the assembly is within tolerance and no shimming is required. For gaps between 0.2 mm and 1.27 mm, peelable polymer shims are utilised. However, for large gaps (≥1.27 mm), aluminium shims are required, which have a significantly higher material cost. Eq. (14) presents a method for estimating the cost of shimming. . The material cost, C _ e material , has a large influence on the cost of the process, as shown in Fig. 11, where the cost of the shimming-based assembly, C AS , is plotted against the assembly variation level, δ KC3 .
It can be seen that there is a large step increase in cost, as aluminium shims are used instead of polymer shims. The cost of shimming gradually increases as the gap size, δ A , increases. This is because the volume of shim, and thus the material cost, increases as the gap size increases.

Cost modelling of milling-based assembly
Milling-based assembly is an alternative method for assembly variation management. Milling-based assembly removes any clashes between the spar and hinge bracket by measuring the degree of clash using a laser scanner, and machining the excess material away from the hinge brackets using a PKM [4,51,52]. In contrast to shimming, the milling method is highly automated, and the main cost driver is the high capital cost of scanning as well as milling equipment. The process is typically much quicker than the shimming process. The milling-based variation management is similar to that described by Morgan et al. [4] and further setup and process details may be found in this reference. Eq. (15) presents a method for estimating the cost of milling-based assembly.
3 , Load fettled parts into assembly fixture Load fettled parts into assembly fixture 8 Apply sealant Apply sealant 9 , and = C T · f Load spar into assembly fixture 10 W Load spar into assembly fixture .
For an assembly clash, δ KC3 , between 0.2 mm and 0.8 mm, only one milling pass is required to remove the excess variation. There is a step increase in cost at a clash of 0.8 mm, as at this point a second milling pass is required to remove the excess variation from the part. For clashes smaller than 0.2 mm, the assembly is deemed to be within tolerance, and sealant is applied between the parts only. The cost of the milling-based assembly, C AM , can be related to the assembly variation level, δ A , as shown in Fig. 12. The main variation-cost driver when two milling passes are required is the time taken to mill the excess material from the hinge bracket, as larger gaps require more time, and therefore more cost to remove.

Monte Carlo simulation for entire wing spar assembly production chain
As manufacturing variations are stochastic by nature, a different combination of variations are induced each time a new part is manufactured. The Monte Carlo method has been used to simulate the effect of random variation in literature [53], and was therefore utilised. The simulation algorithm was repeated 100,000 times to robustly determine the effects of process variation in an assembly. To carry out a Monte Carlo simulation, the variation ranges from Tables 1-3 must be converted into probability distributions. In this work, the process capability, Ĉ p , is utilised to determine the probability distribution of variations of individual manufacturing processes. Ĉ p evaluates the Table 4 Hinge bracket and spar fabrication methods.

Scenario
Hinge bracket fabrication method C PHi ($) Spar fabrication method C S ($) Total fabrication cost statistical variation of a given process, such as the drilling and milling processes used to fabricate the hinge bracket. Ĉ p indicates the statistical range of variations that will be induced by each manufacturing process, which can then be analysed simultaneously using the Monte Carlo method. Ĉ p can be calculated using Eq. (16) [54].
where USL is upper specification limit, LSL is lower specification limit, and^is the estimated variability of the process expressed as a standard deviation.
Variations from each process are assumed to be approximately normally distributed, with the mean centred on the ideal feature dimension. This is for ease of representation of the developed methodology. Normal distribution of process variations have been assumed by a number of processes in literature, such as variation in carbon fibre [55], fixture errors [11], shape errors [19] and machining process errors [56]. In reality, each and every process will have a unique distribution of variations. Non-normal distributions can be used instead of normal distributions without modifying the overall methodology. The standard deviation, σ, can then be calculated using Eq. (16), where the upper and lower specification limits are obtained from the likely variation range for each process.
In this analysis, each variation range considered in Tables 1-3 is assumed to have a Ĉ p of two. This corresponds to a six sigma quality process. This means there are six standard deviations between the process mean and the nearest specification limit. A decrease in the value of Ĉ p represents a less accurate process. According to Eq. (16), a decrease in Ĉ p will reduce the standard deviation of the process probability distribution. This will decrease the percentage of runs that fall within the upper and lower specification limits. When this new probability distribution is inputted into the Monte Carlo simulation, it will result in a larger assembly variation being outputted, therefore leading to a higher cost during the assembly stage. If Ĉ p decreases and all other inputs remain constant, the assembly cost and overall production will subsequently increase.
As variations occur randomly, a Monte Carlo simulation was completed to model the effect of part, process, assembly and fixture variation occurring simultaneously. A flowchart of the Monte Carlo simulation is shown in Fig. 13. First, hinge bracket variations, δ Phi , were picked at random from the probability distribution of hinge bracket variations, which is related to the expected variation ranges shown in Table 1. Then, the spar variations, δ PS , were chosen at random from the probability distribution of spar variations, which are related to Table 2. The spar and hinge brackets were then considered together in the assembly stage, which adds fixture variation and assembly-positioning variation into the analysis. The fixture variations and assembly-positioning variations were also randomly selected from their respective probability distributions, and are related to the variation ranges shown in Table 3. Then, the assembly variation propagation model described in Fig. 10 was used to capture the effects of the assembly process itself. The output of the simulation was the gap or clash at each of the mating interfaces between the spar and hinge brackets, which relate to KC3. When the Monte Carlo simulation was repeated a 100,000 times, a distribution of assembly variation at KC3, δ KC3 , was outputted for each of the mating pairs.
The output of the Monte Carlo simulation, δ KC3 , determines the variation management process required (whether a gap or clash occurs), and also the magnitude of the variation, which affects the cost of managing it. Thus, the distribution of assembly gap or clash magnitudes graph is utilised to connect the variation and cost models together, as the variation magnitude determines the variation management process that will be required at assembly. This is illustrated in Fig. 14, which shows the logic behind the decision to shim or mill excess variations, depending on their magnitudes.

Cost modelling for the entire wing spar assembly production chain
The total cost of the entire wing spar assembly production chain, ∑C, can be obtained by summating the sub-component fabrication costs and the assembly costs together, as described in Eq. (17).
The cost of the processes C PHi and C PS are dependent on which method of hinge bracket or spar fabrication is used (respectively). Each fabrication process introduces random variations of δ PHi and δ PS into the assembly at different magnitudes. As a result, the assembly variation, δ KC3 , has the form of a probability distribution of gaps and clashes at assembly KCs. There are eight KC3 mating pairs between the spar and hinge bracket, corresponding to the eight hinge bracket feet, which locate the spar. As the variations are random and come from many sources, some KC mating pairs in the same assembly may clash whilst other mating pairs may have a gap between them. The size of the gap or clash will also be different at each of the eight assembly mates. Therefore, the cost of removing variation will be different at each assembly KC due to both the size of KC variation, and the process that is used. The result is a different total assembly cost for each run of the Monte Carlo simulation, thus creating a distribution of total production costs. A Monte Carlo simulation of 100,000 runs was completed. Fig. 15 shows the results of the Monte Carlo simulation in the form of a cumulative distribution function (CDF). This shows the distribution of costs for each of the four manufacturing scenarios described in Table 4. Fig. 15 shows the total cost of production for each of the four scenarios plotted against the number of simulations completed. Each of the four scenarios exhibits the same graph shape. For Scenario 1, a finishing cut is used to fabricate the hinge bracket, and the spar is manufactured using a prepreg process. The lowest 43% of costs corresponds to the case where all eight mating pairs are within tolerance ( ± 0.2 mm), therefore no shimming or milling processes are required. The step increase in cost at simulation number 43,000 corresponds to one of the eight mating pairs being out of tolerance, therefore one of the four hinge brackets will require variation management. At simulation 60,000, the cost increases again. This is because two of the four hinge brackets now require shimming or milling: therefore two mating pairs in two different hinge brackets are out of assembly tolerance. There are four step increases in total between simulation 43,000 and simulation 68,000, corresponding to the four hinge brackets in the assembly. At simulation 68,000, there is a step increase of cost of approximately    $1,000. This increase is due to both shimming and milling being required in the same assembly. This dramatically increases the cost of assembly, as the costs of both methods are incurred. The same shaped curve can be deduced from Fig. 15 for each of the remaining three scenarios. It can be seen that the lowest mean cost of production was Scenario 1, when the hinge bracket is manufactured with a finishing cut, and the spar is manufactured using the prepreg method. However, like each of the other three scenarios, it has a very large standard deviation of cost due to the potential of requiring both shimming and milling. Scenario 2 has the most accurate hinge bracket and spar fabrication processes, and therefore has the most cases where all eight mating pairs are within tolerance limits (83% of cases were within tolerance). However, given the considerably higher cost of the RTM process compared to the prepreg method in this case, it is estimated to have a higher mean value than Scenario 1. Scenario 3 has the lowest accuracy (no finishing cut + prepreg process), and thus has the most cases where variation management is required.

Optimisation method of manufacturing processes
In order to optimise the overall production chain, and to minimise ∑C, the optimal variations distribution for δ Phi and δ PS must be found. In the following section, shimming-based assembly and milling-based assembly will be discussed in turn before a comparison is drawn between the two assembly methods.

Optimisation of manufacturing processes using shimming-based assembly method
In the previous cost analysis displayed in Fig. 15, each of the four scenarios resulted in a large standard deviation of costs. This is due to gaps and clashes occurring in the same assembly, meaning the setup and capital costs of shimming and milling are incurred simultaneously. The analysis was therefore repeated but with only shimming allowed to occur. This can be achieved by manufacturing the hinge bracket to be undersized, thus ensuring a gap between it and the spar. For each of the four scenarios in Table 4, the largest magnitude clash between the spar and hinge brackets was found through the Monte Carlo simulation. The hinge brackets were then made to be undersized by this magnitude, thus ensuring that no clashes will occur, and only shimming processes will be required. Fig. 16 presents the CDF of costs for each of the four scenarios when the hinge bracket design has been optimised for shimming-based assembly, Scenarios 5-8.
In this case, gaps and clashes no longer occur simultaneously in the same assembly. As a result, there is no large step increase in costs for each of the four scenarios. Subsequently, when only shimming processes occurs, the standard deviation of all four scenarios is greatly reduced. Once again, Scenario 5 has the lowest mean cost, and Scenario 8 has the largest mean cost. The higher cost of Scenario 6 compared to Scenario 5, despite the RTM process having less variation impact than prepreg, is due to the large difference in cost (∼$800) between the two methods of spar fabrication. It can therefore be concluded that Scenario 6 would be the lowest cost method of spar and hinge bracket manufacture and assembly if the cost of the RTM process was reduced ∼$400 relative to the prepreg process.
The CDF for Scenario 5 and Scenario 6 in Fig. 16 shows the same characteristic steps as Fig. 15, which are associated with each subsequent mating pair being out of tolerance. However, Scenario 7 and 8 exhibit a smoothly increasing line. This is due to all eight mating pairs being out of tolerance for all instances: therefore they all require a shimming process to be carried out.

Optimisation of manufacturing processes using milling-based assembly method
The production chain can also be designed specifically to suit a milling-based assembly method, by manufacturing the hinge brackets to be oversized using the new emerging parallel kinematic technology [4,51,52]. For each of the four scenarios in Table 4, the largest magnitude gap between the spar and hinge brackets was found by means of the Monte Carlo simulation. The hinge brackets were then made to be oversized by this magnitude, thus ensuring that interference will occur, and only milling processes will be utilised. Fig. 17 presents the CDF of costs for each of the four scenarios when the hinge bracket design has been optimised for milling-based assembly, Scenarios 9-12.
In this analysis, only clashes occur. Therefore, there is no step increase in costs greater than $400 for each of the four scenarios. The standard deviation of all four scenarios is subsequently greatly reduced compared to Fig. 15. The CDF for Scenario 9 and Scenario 10 in Fig. 17 shows the same characteristic steps as Fig. 15, which are associated with each subsequent mating pair being out of tolerance. However, Scenarios 11 and 12 exhibit a smoothly increasing line. This is due to all eight mating pairs being out of tolerance for all instances. Therefore, they all require a milling process to be carried out. This highlights that the inclusion of a final finishing cut has a large influence on the final variation level of the assembly despite the small relative cost of the process.
It is noted that the mean cost of Scenario 9 is slightly higher for a purely milling process compared to Fig. 15, where both shimming and milling occurs. This is because the variation distribution of the hinge bracket is no longer centred on the ideal value but is instead centred at the oversized value to ensure there are no gaps at assembly. Hence, there are less cases where all eight mating pairs are within tolerance. However, in the production chain optimised for milling, the standard deviation is significantly smaller than in the non-optimised production chain for each of the four scenarios. Table 5 presents a comparison of costs between the shimming-based and milling-based assembly methods, and the original mixed assembly method.

Comparison between shimming-based and milling-based assembly methods
Scenarios that utilise a prepreg spar manufacture and hinge bracket fabrication with a finishing cut (Scenarios 1, 5 and 9) were found to give the lowest mean cost. This indicates that the cost of spar fabrication is one of the main drivers of overall cost. One difference between milling-based and shimming-based assembly however, is that the second lowest mean cost for shimming-based assembly was found to be when a finishing cut and RTM is used. In contrast, no finishing cut and the prepreg method is utilised in the case of milling-based assembly. This can be attributed to the smaller step increase in cost when variation removal requires two passes instead of one in milling-based assembly (which can be seen in Fig. 12), compared to the larger step increase in cost when polymer shims are replaced by aluminium shims for shimming-based assembly (shown in Fig. 11).
It is noted that in Table 5 the mean cost of Scenario 5 is slightly higher for a purely shimming process, compared to Scenario 1 in Fig. 15 where both shimming and milling occurs. This is because the variation distribution of the hinge bracket is no longer centred on the ideal value but is instead centred at the undersized value to ensure there are no clashes at assembly. Therefore, there are less cases where all eight mating pairs are within tolerance. However, in the production chain optimised for shimming, the standard deviation for Scenario 5 is more than twenty times smaller than in Scenario 1. This is a significant difference, and must also be considered by production planners.
It was found that the shimming-based assembly method was cheaper than milling-based assembly for each of the four scenarios. In this paper, the assembly milling processes were assigned a high wrap rate to represent the large capital expenditure of milling equipment, compared to the manual shimming processes. However, different companies will have different costs and therefore wrap rates. Thus, conclusions made in this paper are based only on this exemplar case.

Conclusion
This paper presents a methodology for holistic process optimisation based on a variation propagation model. The methodology is demonstrated by a real life aerospace case study, i.e. wing spar assembly, in order to minimise production cost. The key deliverables of this paper are: • A variation propagation model for overconstrained assemblies has been developed for the first time, using a KC achievement algorithm to overcome the challenge of ambiguous KC delivery chains.
• Production cost was linked to the variation propagation model using the ABC method, providing a robust method to evaluate cost and variation simultaneously for an overconstrained assembly. As a result, the cost variations resulting from the geometric or process variations can be easily identified.
• A methodology for process optimisation for an entire production chain has been developed, based on variation propagation and Monte Carlo simulation. This provides the ability to analyse the trade-offs between the cost and achievable variation limits of the entire manufacturing chain in order to minimise the overall manufacturing cost.
• The developed methodology will be useful for product designers to allocate suitable dimensions and tolerances, and production planners to arrange the fabrication and assembly processes in an optimised manner. Table A1, Table A2.