Virtual shimming simulation for smart assembly of aircraft skin panels based on a physics-driven digital twin

A leading challenge in the assembly process of aircraft skin panels is the precise control of part-to-part gaps to avoid excessive pre-tensions of the fastening element which, if exceeded, impair the durability and the response under dynamics loads of the whole skin assembly. The current practice is to measure the gap in specific points of the assembly with parts already at their final location, and then be-spoke shims are machined and inserted between the mating components to fill the gap. This process involves several manual measurement-fit-adjust quality loops, such as loading parts on the assembly frame, measuring gaps, off-loading parts, adding be-spoke shims and re-positioning parts ready for the fastening operation—as a matter of fact, the aircraft is re-assembled at least twice and therefore the current practice has been proved highly cost and time ineffective. Additionally, the gap measurement relies on manual gauges which are inaccurate and unable to follow the actual 3D profile of the gap. Taking advantage of emerging tools such as in-line measurement systems and large-scale physics-based simulations, this paper proposes a novel methodology to predict the part-to-part gap and therefore minimise the need for multiple quality loops. The methodology leverages a physics-driven digital twin model of the skin assembly process, which combines a physical domain (in-line measurements) and a digital domain (physics-based simulation). Central to the methodology is the variation model of the multi-stage assembly process via a physics-based simulation which allows to capture the inherent deformation of the panels and the propagation of variations between consecutive assembly stages. The results were demonstrated during the assembly process of a vertical stabiliser for commercial aircraft, and findings showed a significant time saving of 75% by reducing costly and time-consuming measurement-fit-adjust quality loops.


Introduction and motivation
The assembly process of aircraft skin panels, such as wings and fuselage, requires joining of multiple skin panels through fastening and riveting operations. This requires on the order of thousands of custom shims to fill part-to-part gaps between the deformable skins and internal skeleton (i.e., ribs, spar). Shims can considerably reduce the assembly stress on mechanical fasteners and improve the joint stiffness and load capacity, and this effect is more remarkable with the increase in gap values. So, a strict control method must be implemented to ensure the part-to-part gaps are within the targeted control limits (typically below 1 mm) [1]. Part-to-part gaps are the consequence of accumulated manufacturing tolerances and variations incurred during the sequential stages of the assembly process [2,3]. The shimming operation undergoes multiple timeconsuming quality loops of measurement-fit-adjust that have been proved highly cost and time ineffective. These operations involve pre-positioning of parts, measurements of gaps, off-line fabrication of shims, off-loading of parts, placement of fabricated shims, re-positioning of parts and riveting [4]. It is worth noting that large aircraft sub-assemblies such as fuselage, wing and stabiliser, are built in different manufacturing plants often placed at far geographic locations. This implies that the quality loops may take weeks or even months before all the quality requirements are fulfilled. Quality loops, as indicated in [5], require between 10 to 14 days to complete and correct defects in the barrel assembly of commercial aircraft. This has detrimental effects in terms of increased cost, up to 25%, and elongated time-to-market [6].
A critical element is the measurement of the part-to-part gap which currently relies on manual gauges such as feeler gauges, as shown in Fig. 1a, or capacitive sensors. After measuring the gap in specific points the shape of the shim is developed by interpolating those points and then machining the mating surface in the CNC machining centre, or 3D printing [7]. The shim fabricated following this procedure is then positioned and glued (to avoid any unwanted drop or movement of the shim during the sub-sequent part re-positioning) before parts are fastened together- Fig. 1b illustrates a typical shim generated with this approach and retro-fitted to the nominal design of the skin panels. This approach has three challenges. Challenge C1, the measurement, being manual, is inaccurate and unable to follow the 3D profile of the gap. Challenge C2, only one surface of the shim is machined: the bottom surface of the shim is assumed flat, and this does not correspond to the actual gap distribution between parts. Challenge C3, gaps are measured in a "closed-configuration", where parts are mounted on the same assembly jig to form a closed structure, and access to each gap area may be cumbersome and therefore the measurement accuracy would be impaired.
Alternative solutions to shimming were also proposed and are discussed as follows: (1) in-situ printing of shimming-in [8] shims have been manufactured directly on the part being joined by a 3D printer attached to an anthropomorphic robot moving over the structure. (2) Liquid shimming-since cracking and delamination can appear around the mechanical fasteners of composite structures due to unpredictable re-distribution of loads, liquid shims (epoxy-based adhesive) have been introduced as alternative to mechanical fasteners [9]. They are commonly used in bonded applications, such as the composite rib-to-skin assembly. Although they offer high compressive strength properties, they suffer under shear and traction loads. (3) Local machining of panels-in [10] a method has been described to eliminate the shimming operation by developing automated processes to remove excessing material from the appropriate interfaces using rapid metrology techniques and robotic systems. This technique was claimed to be more advantageous than working with shims. On the same line, in [11], a possible solution to the high variability of the manufacturing process of an aircraft wing assembly was described by using a reconfigurable assembly process which was shown to be highly repeatable.
Research efforts have been devoted to model and predict both the variability of parts and the final shape of assembly. The variation propagation has been modelled via various techniques such as kinematics-based models, assuming rigid parts [12], and Finite Element Methods (FEM) approaches [13,14] to model parts deformations. However, these methodologies are limited by the fact that they assume parts are assembled in a single-stage operation. Other studies have dealt with the fixture and locating scheme optimisation to reduce part deformation occurring during the assembly phase. For example, determinant assembly holes have been used to locate parts to each other and hence to guarantee the correct hole-to-hole alignment; however, clamping forces cause local elastic deformation of parts and subassemblies. Using a statistical approach with Monte Carlo simulation, a method has been proposed in [15] to predict assembly variations, accounting for deformations occurring due to the part-to-part interaction during the assembly process. Nonetheless, the dimensional variation simulation alone appears not to be able to accurately estimate part-topart gaps. A possible route to enrich variation simulation approaches is via experimental scanning data. The work described in [16] combined scanning data and FEM simulations in order to predict the 3D shape of the shims for wing assembly process. The method was based on meta-heuristic optimisation and single-stage modelling of the assembly process. Manhoar et al. [17] proposed a methodology to predict part-to-part gaps in aircraft assembly with machine learning and sparse sensing.
In order to overcome the three aforementioned challenges (C1-C3), this paper proposes to measure in-line the mating surfaces of parts being assembled while they are in an "open configuration" (Fig. 2). This opposes the state-of-theart approach which manually measures gaps only when parts are positioned and aligned together on the same assembly jig. The proposed approach allows to predict and fabricate shims using in-line measurement data with no need to iterate the measurement-fit-adjust quality loops. This will undoubtedly reduce inspection and re-work costs, enabling operators to virtually test assembly operations before deployment and installation in the field. The proposed approach shifts the problem from a metrology-driven challenge to a prediction challenge. For instance, we propose a physical-driven digital twin in the sense that integration between in-line measurement data (physical domain) and physics-driven variation propagation (digital domain) is deemed. In our definition, the physics-driven digital twin makes use of high-fidelity simulation based on first-principle physics (i.e., deformation under static and linear regime) to predict the response of the system. The high-fidelity simulation is complemented by measurement data [18]. The technological challenge now is the need to predict the geometrical configuration of parts once they are moved from the open configuration to the closed one. The challenge is driven by the fact that (1) skin panels may deform under the effect of gravity; (2) and/or, be pushed away by the pressure exerted by the ribs subject to manufacturing tolerances. This is illustrated in Fig. 3, where GRS and LRS are the Global and Local Reference Systems, respectively. Once parts are scanned in-line (Fig. 3c), the gathered cloud of points (defined into LRS) is aligned to GRS. State-the-art methods, based on rigid rotations/translations (for example, iterative closest points or best fitting alignments) are uncapable to model physical effects, such as gravity and part deformations (i.e., the elastic deformation induced by partto-part contact). Finding the correct alignment is therefore the technological challenge. The fact that parts are moved from the open to the closed configuration brings the need to model the variation propagation in a multi-stage scenario. In order to clarify this point, we have used the case of the vertical stabiliser where the internal skeleton (front and rear spars, ribs located by rib-post and shear-tie angular flanges) is located on a large jig, while left and right skins are located on separate lateral structures. The assembly of internal skeleton comprises of 3 stages: (1) location of front and rear spars; (2) assembly of rib-post and shear-tie angular flanges; (3) ribs located and assembled by using the datum offered by the rib-post and shear-tie angular flanges. Once the internal skeleton is completed, it is married to the lateral structures (stage (4)-skins to internal skeleton). The need to precisely Proposed methodology for the prediction of part-to-part gaps and the generation of 3D shims model multi-stage assembly processes has been also highlighted in Franciosa et al. [19]. Findings indicated that in most cases the quality of the final assembly is strongly influenced by complex interactions which might generate defects within individual stages and their propagation from stage-tostage.
The remainder of this paper is organised as follows: Sect. 2 outlines the methodology. Section 3 describes the implementation of the methodology while the case study is reported in Sect. 4. Section 5 draws final conclusions and remarks.

Proposed methodology
The proposed approach for generation and prediction of shims is illustrated in Fig. 4. The paper combines the need to model and simulate the variation propagation of dimensional and geometrical errors in multi-stage assembly operations [20,21], along with a predictive tool to model the shimming condition between parts being assembled.
Two areas are central to our approach: (1) variation model of the individual skin panels generated via a morphing algorithm and calibrated through surface scanned data; (2) and, variation model of the multi-stage assembly process via a physics-based simulation which allows to capture the inherent deformation of thin and large skin panels. The morphing algorithm is used to complement measurement data and generate the "variational non-ideal" shape of the skin panels. This approach is particularly useful when measurement data are unavailable or incomplete.
The proposed methodology is presented and discussed using the assembly process of the vertical stabiliser for a commercial aircraft. The gap distribution is computed for several assembly scenarios considering the deformation of the skin panel and the variation arising from misalignments of both structural and fastening components. Then, the 3D shape of the shims is generated.

Step (1): Generate form errors before assembly
This step involves the development of the variation model of the individual skin panels (before parts are assembled). It is assumed that the skin panels, being flexible, are subject to form errors, whereas the internal skeleton (ribs and spars) is reasonably rigid. The variation of the internal skeleton is only associated to rigid rotations/translations of the individual components, caused by manufacturing tolerances of the mounting points (fasteners, locating holes, etc.). Therefore, form errors are neglected for the internal skeleton. State-of-the-art methods are used to model rigid rotations/translations. Conversely, we propose to complement measurement data with a data-driven model to generate form errors. If the measurement systems (i.e., camera-based or laser-based scanners) are available, surface scanning data is acquired and, in this case, deviation patterns are generated by mapping surface data to the nominal CAD geometry, G N . Otherwise, and when data are partially missing, deviation patterns are simulated using GD&T design tolerances. The morphing mesh method is used for the purpose of emulating the form errors of the skin panels. It allows generating a deviation pattern in accordance with geometrical tolerance specifications. This technique allows to parametrise any 3D CAD geometry and to embed deviation patterns.
The fundamental idea behind the method is to define a set of control points (see Fig. 5a) and corresponding influence hulls. Each control point determines a local deformation of the surface acting along the direction of deformation. Control points are taken from the quality inspection plan-production data (either scanning data or historical data) are used to extract the temporal-spatial correlation, which determines the dimensions of the influence hulls. The influence hull defines the 3D region within which any point is influenced by the related control point [22]. The deviations (and the related statistical distribution) assigned to each control points are experimentally calibrated using historical inspection data (control plans) of the skin panels. Hence, the contribution of each control point accounts for the natural variations incurred during the fabrication of the panels.
Step (2): Simulate multi-stage assembly process This step involves the physics-based simulation and variation propagation during multi-stage assembly process [23] considering part deformation and the deviation patterns computed in STEP (1).
The model follows the fundamental first-principle laws of deformation under static and linear regime. This assumption is deemed reasonable under small deviations and pre-defined load and boundary conditions as exerted by locating and fastening points. The model is solved via FEM. The mesh representation of the assembly is generated from the nominal CAD. Data about the spatial location and orientation of both the locating and fastening points is input to the model. The generalised multi-Place-Clamp-Fasten-Release (m-PCFR) approach is used to model each assembly stage (Fig. 6).
The proposed model for variation propagation in multistage assembly stages is conceptually represented in Eq. (2), being U A the cumulative deviation pattern generated by the assembly process. U A,k is the deviation (Eq. (2a)) accumulated throughout the assembly stages, from S[1] to the S[k] stage; Ω k is the ensemble of all parts belonging to the k-th assembly stage.
The variation propagation model has been implemented in the VRM simulation suite [24] which is a MATLAB-based finite element modelling software toolkit with a dedicated FEM kernel [25].
Essential part of the simulation of the multi-stage assembly process is the need to avoid part-to-part penetration. The approach integrates the FEM kernel with an iterative contact solver which makes use of the node-to-surface contact search approach. The contact solver converges if and only if the gaps of all the active contact pairs are greater (or equal to) then a pre-set gap tolerance, and all pairs are in compression (negative load). A penalty method has been implemented to enforce these conditions [26]. At each iteration of the iterative solver, part-to-part gaps, g, are computed with the source-todestination projection (Fig. 7), as expressed in Eq. (3), for any source point P G,s projected on the destination part, P G,d , along the local normal vector to the surface, N m .

Step (3): Generate 3D shims
This step entails the calculation of part-to-part gaps and generation of 3D shims. The idea is that, once the simulations in STEP (1-2) have been completed for all assembly stages S[1] to S[n], a dense set of points (namely cloud of points) is sampled from the surface of the mating surfaces. The approach is illustrated in Fig. 8 and is discussed as follows: (a) the planar envelope of the shim is imported via reference CAD surface. This is a design input and is usually included in the product specifications; (b) a local reference system is built using the standard principal component decomposition applied to the nodes of the reference CAD surface; (c) a regression algorithm seeks to find the analytical surface which best fits the sampled surface data (dots in Fig. 8c). The regression is computed in the local reference system. Multiple regression techniques can be used for generating the analytical surface (i.e., polynomial fitting, spline, nurbs, etc.) pending on the distribution of the sampled data, accuracy and geometry topology. It is worth noting that the regression algorithm tends to smooth out local disturbances/noise in the cloud of points and generate a continuous surface. Furthermore, being analytical, designers can extrapolate the thickness of the shim for any point belonging to the envelope of the shim itself; (d) the analytical surface is exported in polygonal mesh format (i.e., STL format) and ready to be fabricated (i.e., 3D printed via additive manufacturing). The mesh resolution is a user parameter which is selected pending upon the approximation level of the geometrical topology.

Implementation
The predictive simulation tool for part-to-part gap has been implemented as a new module in the VRM suite, and named Virtual Shimming Simulator (VSS).
The tool allows to interactively define the assembly sequence, control points, locating and fastening points. VSS can parametrise any input variables and therefore large sets of "what-if" assembly scenarios can be "virtually" investigated. Figure 9 shows the Graphical User Interface (GUI).

Case study
The proposed methodology has been tested and validated on the vertical stabiliser of aircraft assembly system, shown in Fig. 10. The aim of the study is to control the part-to-part gap between the internal skeleton (spars and ribs) and the two skins (left-and right-handed). The following assumptions have been made: (1) spars, ribs and rib posts are assumed rigid compared to the skins; (2) skin panels are assumed deformable; (3) the material of the skin panels is CRFP; since the fibre layout is randomly distributed, the material is assumed isotropic with an equivalent Young's Modulus experimentally estimated and equal to 90 GPa.
Two cases have been considered: case [1]-sensitivity study of variation propagation to inform design decisions; case [2]-validation using physical demonstrator.

Case [1]-Sensitivity study
Six assembly scenarios have been simulated in order to investigate the impact of manufacturing tolerances on the part-to-part gap between skin panels and the structure of the skeleton. Each scenario is constituted by three assembly stages: S[1]-fabrication of skins; S[2]-assembly of internal skeleton (ribs and rib-post with shear-tie angular flanges); S[2]-assembly of skins to internal skeleton. Form errors are modelled for the skins in S [1]. Rigid rotations are modelled for S [2] to account misalignment of both ribs and rib-post. Variations in the orientations are used to emulate manufacturing tolerances of the locating holes. Table 1 summarises the simulated assembly scenarios.
The orientation error of the rib is generated with respect to the reference system (Fig. 11a) centred in the middle of the rib. This assumption has been made to take into account the fact that the rib goes in contact with the two spars and therefore any error tends to be symmetric around the middle axis of the rib. For each locating hole of the rib post a Nominal skin panels Fig. 11 Definition of the coordinate reference systems for the simulated assembly scenarios local reference system is defined as shows in Fig. 11b, and attached to the locating hole of the rib post itself. The form errors of the skins have been generated with the morphing mesh approach- Fig. 12 shows the adopted control points and related influence hulls. The deviations assigned to each control points have been experimentally calibrated using historical data. The first assembly scenario, S1, highlights the effect of the rib misalignment. Thus, the only simulated variation is Scenarios S2 to S6 focus on the effect of the misalignment of the rib posts. The gaps are computed between the rib-posts and the skin in the middle of the rib post. 100 different orientations around the Z-axis where randomly generated assuming a normally distributed error. Furthermore, S4 and S5 are carried out with a local deformation around the interface with the 7-th rib post. Conversely, S6 assumes that the skin panels are nominal (no form error). Figure 13a,b show the representative results of the three assembly stages. It is worth noting that the part-to-part interference that happens in S [2] (the misalignment has been magnified 100 times to appreciate the theoretical interference with the skins) is then resolved in the clamping stage, S [3], after solving the physics-based simulation of the assembly. Part-to-part gaps have been measured on the 4 rivets shown in Fig. 11a-Rivet [1] & [2] and Rivet [3] & [4] corresponding to the left and right skin, respectively. The maximum acceptable gap has been set to 1.0 mm. Figure 14 shows the results of the sensitivity study for the assembly scenario S1. The charts have been generated by sampling the rotation α and γ in the range [−1, 1] degree. For each sampled point, the multi-stage simulation has been run (STEP (2) of the proposed methodology) and deviation pattern computed. A piece-wise polynomial regression has then been computed to Fig. 13 a First assembly scenario, S1: deviation patterns (in mm) for the vertical stabiliser. Deformation has been magnified 100 times. The colour code represents the deviation normal to the nominal surface. b Second assembly scenario, S2: deviation patterns (in mm) for the vertical stabiliser. Deformation has been magnified 100 times. The colour code represents the deviation normal to the nominal surface Fig. 14 Results of scenario S1-impact of the misalignment of the rib on the part-to-part gaps between skins and rib. The colour code represents the part-to-part gap in mm. Rivet [1][2][3][4] are in Fig. 11  Fig. 15 Distribution of part-to-part gaps for scenarios S2 to S6 fit the sampled data (total 100 data points) with R 2 95% and RMSE 0.05 mm. This result is significant since informs the designers about possible un-feasible design decisions. For example, assuming α 0.5 and γ 0.2, a gap of approx. 0.8 mm is predicted for Rivet [1]. Results suggest that part-topart gaps cannot be compensated only by re-orienting the rib. For example, though Rivet [1] shows a feasible region (gap < 1 mm) for α [0, 1] degrees and γ [−1, −0.5] degrees, the same does not apply for Rivet [3], which is placed on the opposite side of the rib. Same conclusion applies to Rivet [2]. As consequence, corrections made on the left skin may tend to further increase the gap on the right skin. Designers may consider to change the assembly sequence to avoid the detrimental impact of the variation propagation from S[1]-S [3].
The misalignment of the rib posts has been modelled with a normal distributed error with a standard deviation of 0.1 degrees (and zero mean). The chart in Fig. 15 reports the gap distribution for the 9 rivets placed on the rib posts, and for all the simulated assembly scenarios, S2 to S6. Comparing S2, S3 and S6, the results suggest that the form error of the skin panels gives a minor contribution towards the gap values, which, instead are significantly affected by the manufacturing tolerances of the rib posts. This result is backed up by the fact that the 1-way ANOVA gives a p-value higher than 10%, thus suggesting that the variation among the analysed scenarios is neglectable. However, for scenarios S4 and S5, the gap values are the highest on average. This confirms the fact that the additional local deviation in the skin panels (error in the control point ξ 1 ) contributes significantly to the variation in the part-to-part gap. With those results in mind, designers may consider imposing stricter control of the form errors of the skins to guarantee better control of the part-to-part gaps.

Case [2]-Validation using physical demonstrator
The proposed methodology has been tested and validated on the physical demonstrator shown in Fig. 16. The aim of the study is to compare the results obtained with the current best practice approach using manual feeler gauge, against the predictions obtained using the proposed methodology. Total 40 shims were generated following the methodological STEP (1) to (3) between skin and ribs/spars. Results are discussed as follows: • Accuracy of the prediction-results are illustrated in Fig. 17 and findings show a strong correlation (approx. 98%) between the average shim thickness measured by feeler gauges and predicted using the proposed methodology, with mean squared error of just 0.08 mm. • Time saving for measurement of part-to-part gaps-it was reported that 2 h were spent to gather the scanned data (measurement made in the "open-configuration" of the Fig. 17 Average part-to-part gap (shim thickness) measured by feeler gauges and proposed Virtual Shimming Simulator (VSS) assembly) of both ribs/spars and skin (the computational time to generate the shims using the proposed methodology was only few minutes, thus neglectable comported to the scanning time). Whereas, near 8 h were required to collect individual measurements using manual feeler gauge (measurement made in the "closed-configuration" of the assembly). This equates to a significant time saving of 75%.

Conclusions and final remarks
Motivated by the need to accurately predict the part-to-part gap between aircraft skin panels, and therefore minimise the need for manual quality loops, this paper has presented a novel methodology which combines two main elements: (1) variation model of the individual skin panels generated via a morphing algorithm and calibrated with historical data; (2) and, variation model of the multi-stage assembly process via a physics-based simulation which allows to capture the inherent deformation of the panels. The impact of the research is discussed as follows: • Academic: state-of-the-art approaches are only limited to single-stage models for assembly systems of skin panels. The proposed methodology presents a novel and unique approach to model variation propagation in multi-stage assembly systems with deformable parts. Non-ideal parts can be generated either using scanning data or only based on design tolerances. • Industrial: the proposed approach will help to reduce number of quality loops necessary to reach zero-defect manufacturing. The developed predictive tool offers the unique capability to eliminate costly and time-consuming measurement-fit-adjust quality loops. Designers and manufacturing practitioners will benefit of this tool to automatically predict the part-to-part gap and will be able to "virtually" generate the 3D geometry of the shims. The experimental validation corroborates the fact that the proposed approach is a viable solution for handling assembly of aircraft skin panels with a significant time saving of 75%.
This work informs on the broader spectrum the possibility to deploy a "digital twin" representation of the physical assembly process to support interactive design and achieve right-first-time, thus reducing manual and costly trial-anderror iterations. Future work will be devoted to automating the process of fabrication of 3D shims using the deviation pattern calculated with the Virtual Shimming Simulator. Optimisation of assembly sequence will be also taken into account.