A Digital Twin for Automated Layup of Prepreg Composite Sheets

Composites manufacturing methods form an integral part of the aerospace sector. Recently, the aerospace industry has gravitated towards an increase in the use of composites. The use of composites is projected to increase with a compound annual growth rate of about 10% year on year [1]. This growth in the utilization of composites has stirred a demand for automating functional tasks in the composite manufacturing processes. Out of the myriad of manufacturing methods, pre-impregnated with resin (prepreg) composite sheet layup method has been of key interest for the aerospace industry. Prepreg sheets are advantageous due to their controlled volume fraction, simple and inexpensive tooling, and ease of handling. These materials can be customized in terms of the types of fibers, weave type, and number of plies, allowing for diversification of component production and a broad range of applications [2,3]. Although the multi-component material systems delivered by prepreg sheets carry improved functionality, their manufacturing method is still susceptible to defects and quality issues. These defects can be classified as air gaps, wrinkles, bridging, etc. Currently, these layups are predominantly executed by skilled operators who place each ply (prepreg sheet) manually on the tooling and apply localized pressure using their hands or custom hand tools.

The high factor of human involvement in the layup processes introduces a potential variability and error in each subsequent layup. With an increase in demand for composite parts, it becomes crucial to increase production rates and improve quality to maintain process capability measures. This can be achieved through automation of the composite manufacturing methods. Presently, automated techniques such as automated fiber placement and automated tape layup are employed primarily for tooling which does not involve complex features. Automating a prepreg layup process introduces challenges that can only be addressed if the prepreg sheet behavior can be predicted accurately beforehand. For automated sheet layup, rectifying a defect in-process can be arduous as prepregs once draped are difficult to separate. Consequently, the separation can result in additional defects such as bridging, fiber misalignment, and impairment of previously laid up base sheets. As a result, executing a layup perfectly at once becomes a pivotal criterion for an automated prepreg layup process. This necessitates a need for a precise and low-latency simulation system that can accurately predict prepreg sheet behavior during the automated layup process. Such a simulation system acts as a digital twin of the prepreg and can be utilized to plan automation tasks.

In previous work accomplished at the Center for Advanced Manufacturing (CAM) at the University of Southern California (USC), Viterbi School of Engineering, the use of robots for performing automated prepreg sheet layup has been demonstrated [4–8]. Such a highly autonomous robotic cell needs estimation of the sheet parameters in situ for optimal and efficient planning. Along with accuracy, a high data transfer rate plays a critical role. In this research, the aim is to study the identification of different material parameters of a prepreg sheet, namely, bending stiffness, tensile stiffness, shear stiffness, damping stiffness, and surface density. These material parameters can help develop a FEM simulation model that can estimate and accurately simulate prepreg sheet behavior.

For developing the thin-sheet FEA model, we used the VegaFEM library [9]. In this paper, a modeling technique to accurately simulate prepreg sheets under fixed constraints is presented.¹ An experimental methodology of recording sheet characteristics that can aid in the determination of key material parameters is presented as well. These recorded data are used to train a model that can evaluate the parameters with an acceptable accuracy using the Vega FEM library. The evaluated parameters are subsequently used to construct a digital twin represented by a force, damping, and mass matrix with the ability to emulate prepreg sheet behavior under external fixed constraints. The study then focuses on model evaluation and testing for different conditions. A detailed comparison of the parameter model predictions and experimental data is presented as well. We test our methodology on other materials such as cotton cloth, felt, and canvas. Additionally, we introduce a real-time sheet tracking system with depth image sensing to track the sheet for in-process validation.

The proposed in situ simulation system can aid in the development of an autonomous human–robot collaborative cell for prepreg sheet layup. In order to test the modeling concept, we developed a robotic cell setup and the sheet simulation as depicted in Fig. 1. Such a simulation system can particularly assist in the development of a grasp planning system that can aid in the handling of the prepreg sheets by robots. The potential applications of this simulation system are also discussed in detail.

Mechanical simulation of composite prepreg sheets is extensively used in the context of predicting the draping of prepreg over a mold. The sheet is represented by a mesh and fiber alignment, shear, and bending of the sheet can be predicted. A simple and computationally fast approach is to use kinematic simulations which only consider the mold geometry [11,12]. More accurate models use elasticity theory to compute the strains within the fabric [13]. Advanced finite element analysis (FEA)-based models are used to capture the deformation mechanics of the cloth [14]. In Ref. [15], a robot was used to place the flexible material on a doubly curved guiding tool. They used FEA-based models to predict the conformity. While FEA models are more accurate, they are often computationally expensive. In our work, we utilize a variant of FEA developed in the computer graphics community, providing a good trade-off between simulation speed, stability, and accuracy [16,17]. As a result, our method is not just reasonably accurate but also runs at interactive rates which permits us to perform more design iterations.

A survey work done in Ref. [18] presents an overview of all the techniques and models used to simulate fabrics. An overview of recent advancements in automated composite draping was presented in Ref. [19]. Progressive drape models which are a hybrid between kinematic and FEA simulations have also been proposed [20,21]. Such methods improve the accuracy compared to the purely kinematics model and also reduce the computational expense compared to FEA. Other commonly used models are particle systems [22–25] which nonetheless suffer from an accuracy loss due to their inherent particle discretization.

The parameters of the model need to be tuned to match real-world thin-sheet observations. For cloth, this can be done using Kawabata plots [26]. Prepreg composite sheets are sturdier than cloth, however, necessitating precise displacement and force measuring equipment, rendering Kawabata plots less applicable. Instead, we develop an approach which uses FEA simulation and optical motion tracking to tune the material parameters for a given sheet, without explicitly measuring any internal elastic forces.

We use thin-shell FEA to simulate the behavior of viscoelastic prepreg materials. Our thin-shell mathematical simulation model was described in computer graphics references [16,17], and we summarize it in this chapter for completeness. We note that these models were designed as a tradeoff between computational efficiency and accuracy, and have not been previously applied to real-world structures such as composite prepreg. We stress that accuracy alone is not sufficient for a successful material optimization system, as the speed by which design iterations can be explored also greatly affects one’s ability to obtain satisfactory results. Our system runs interactively at 10 fps, which enables rapid exploration of the material parameter optimization space. We incorporate layup domain knowledge, realistic constraint-handling, and real-time tracking of sheet deformation and optimize the material properties to match real-world prepreg observations.

The prepreg is represented as a triangulated mesh with x, y, z displacements of the vertices treated as simulation degrees-of-freedom (DOFs) of the model. The internal forces as a result of bending, shear, and stretching of the mesh govern the behavior of the material. Therefore, these internal forces need to be computed based on laws of elasticity for mechanical simulations. In our approach, we also compute the bending and tensile-shear force Jacobians to improve the computational efficiency of simulations.

Our algorithm takes the triangle elements of the mesh as an input. We parameterize the 3D surface of the material using 2D parameter space represented by parameters u and v. Corresponding 3D weft and warp direction vectors U and V of the mesh can be computed in terms of these parameters. These vectors need not be orthonormal to each other after deformation. Consider the undeformed and deformed states of a triangle element in the mesh to illustrate the concept shown in Fig. 2. The 3D vertices of the triangle P_a, P_b, P_c are computed as a function of u_a, v_a, u_b, v_b, and u_c, v_c.

The stiffness component governs how stress–strain relationship affects the forces acting on the vertices of the triangle. The new position of the vertices is then found by taking into account internal and external forces during the simulation. The force acting due to bending stress is also superimposed with the tensile and shear forces to improve simulation accuracy. We will now discuss the method to compute the bending force.

In this section, we discuss the methodology of estimating material parameters for composite prepregs, cotton cloth, felt, and canvas. The experimentation was accomplished in three stages: data acquisition, data and parameter training/testing, and parameter optimization. Data acquisition was conducted within a physical environment through a guided manipulation of the various material sheets and collecting point cloud data of the position of each sheet as it moved. The training data optimization was conducted in a purely computational space, repetitively running the acquired data through an optimizing simulation sequence to characterize optimal parameter outputs. The outputs were then run through the simulation testing sequence to produce the final results. This can be seen in Fig. 4. This method was completed for all four material types separately. The parameter estimation approach will be discussed comprehensively in Sec. 4.3.

The data collection was conducted using a Hexagon RS5 laser scanner and contact probe attached to the Romer absolute arm (87-Series). The class 2M laser scanner is hand operated and generates point cloud data. All scanning equipments are shown in Fig. 5. The sampling filter can be manually set to optimize the percentage of points recorded and exposure time based on the light in the sampling environment and color of the component being scanned to reduce noise and capture the target locations. The settings can be seen in Table 1.

The contact probe was calibrated and operated using a TESA TKJ 3 mm Ruby Ball Probe. The probe and laser were subject to accuracy specifications, designated by the manufacturer in Table 2.

To prepare each sheet for scanning, quarter-inch markers were placed at even intervals along the sheet in a 17 × 17 grid, for a total of 289 markers. The markers used were white 3M double-sided foam tape squares, cut to the correct size. The markers were raised from the surface, allowing for easier detection by the scanner. One of the sample materials used was white in color and the markers had to be colored in black to be picked up by the scanner, otherwise the color difference was considered optimal for scanning purposes with the settings previously mentioned. The four locations of the clamps obscured the markers and were recorded separately via the contact probe for use within the simulation as a fixed point (Fig. 6).

These data were collected during two trials for each material tested. Each edge of a sheet was grasped using four vertically fixed clamps, allowing the sheet to rest suspended between them, as shown in Fig. 7. The clamping locations were chosen at differing and unique distances from adjacent corners. The collection period for each trial featured a five-stage setup, with each stage introducing a new modification to sheet positioning by moving the clamp with respect to the room, while retaining the clamping position on the sheet. The translational movement of the clamp was chosen arbitrarily, but all movements had a change in distance of no more than 350 mm.

The first stage required placing the sheet under the initial boundary conditions of the first stage, resulting in the sheet suspended in a relaxed, horizontal position. The next three stages involved isolated movement of only one grasping location to a new position, with a general increase in the z-axis, coinciding with x- and y-axis movement toward the center of the sheet. The last stage for every trial was a release of the clamp, allowing the material to settle into a hanging position. Within each stage, each clamp changed location only once and was allowed to settle into position, until no visible movement could be detected. At the beginning of each stage, after clamp relocation had occurred, the new clamp locations were collected via the contact probe and the sheet was scanned via the laser scanner. This procedure can be observed in Fig. 8.

The point cloud data collected from each stage were processed using pc-dmis and exported as an XYZ file for post-processing. From each stage, each sheet contained thousands of data points, with each marker averaging 300 points. The simulator requires only one point from each of the markers, so each XYZ file was processed through blender software using an original python script to isolate the center point of each marker and export them as 3D points, relative to the world coordinate system of the scanner base. A comparative image of the same sheet is shown in Fig. 9, demonstrating the reduction from the point cloud marker clusters to single points.

On successful construction of the initial mesh, the boundary conditions of the sheet are changed by releasing one of the clamps. Note that the remaining clamps should not be moved to maintain consistency in the entire process. The sheet state after changing the boundary conditions is defined as released state. Figure 8 shows the difference between the initial state and the released state.

Since no mesh is required to generate from the released state, the point cloud data of such a state are clustered to represent each marker. The data set acquired in this process is further defined as the observed data. Figure 9 shows the point cloud data clustering process.

The crucial element of the parameter estimation process is the composite sheet simulation system. Figure 10 gives the block diagram for the proposed simulation system. The simulator system utilizes the initial mesh as geometric information input and applies the model parameters to construct the composite sheet model.

The model parameters are categorized into two types: (1) material parameters and (2) integrator parameters. Material parameters consists of the surface density and the internal force parameters, which are tensile stiffness, shear stiffness, and bending stiffness. On the other hand, integrator parameters include damping stiffness and damping mass.

Once the sheet model is constructed [27], the numerical integrator applies the boundary conditions and external forces, such as gravity, to the sheet model and solves the deformation equation [28].

After the predicted mesh is generated, the prediction error is obtained by comparing the predicted mesh and the observed data. The prediction error, E, is a function of model parameters, P, initial mesh, M, and observed data, O. Algorithm 1 is used to calculate the prediction error function.

Let d be a data point in the observed data, O. After getting the predicted mesh, M_pre, by applying P and M to the simulator, we query each d to find the closest vertex, v, on M_pre. The prediction error is then defined as the average displacement between v and d plus half of the maximum displacement across the entire mesh. The prediction error is used for parameter optimization, which will be discussed in Sec. 4.3.3.

The model parameter identification uses a nonlinear optimizer to compute the optimal parameters for simulating the composite prepreg. However, the optimizer is not required to compute all model parameters. Some of these parameters can be measured directly.

As mentioned in Sec. 4.3.1, model parameters comprise of material parameters and integrator parameters. Since integrator parameters are not related to model construction, we set damping stiffness and damping mass to 1.0 and 0.0, respectively. The surface density can be measured directly by scaling the sheet and dividing it by the surface area. Therefore, the remaining parameters, tensile stiffness, shear stiffness, and bending stiffness, are the parameters that require optimization.

To ensure effective performance of the optimizer, appropriate upper bound, lower bound, and initial values for the parameters are required. Figure 11 shows the input–output diagram for parameter boundary selection. Tensile stiffness and shear stiffness are sampled into three categories: (1) 500, (2) 5000, and (3) 50,000. For bending stiffness, the parameter is sampled into (1) 0.01, (2) 0.001, and (3) 0.0001.

After getting three samplers for each parameter, we shuffle them and get 27 candidates to test for parameter feasibility. All candidates that cause system failure are eliminated, and those who have the smallest prediction error, or are comparable to the smallest one, are highlighted. Candidate {5000, 5000, 0.001} has the best overall performance among 27 candidates, and therefore, it is selected to be the initial values set for the parameter optimization. Then, the initial values are used as the median for the parameter boundary. Thus, the upper boundary for the {tensile stiffness, shear stiffness, bending stiffness} is {9000, 9000, 0.0011}. The lower boundary for the parameter set is {1000, 1000, 0.0009}.

The optimization library used in this work is NLOPT [29], an open-source library for nonlinear optimization. The selected algorithm is improved stochastic ranking evolution, a global-gradient-free optimization algorithm [30]. Figure 12 shows the block diagram for the optimization process.

The optimization problem can be expressed as Eq. (11). Let M be initial meshes in the training set, which contains mesh data, m_i. O is the observed data set, which contains observed data o_i. Recall that the prediction error, E_i, is defined in Algorithm 1. Since five samples are used for model training, k is set to 5. The optimizer tries to find the parameters that reduce both the average displacement difference and the maximum difference between the predicted mesh and the observed data for each training set. This error value was set to (average error + 2* Max error). Once the error converges, the identified parameters, P_opt, are used to simulate the composite sheet in real-time.

Four materials are considered in this work: prepreg composite sheets, common cotton cloth, felt, and canvas. The composite sheet was supplied by Boeing Inc. and came as 3 ft × 4 ft sheets, while all other sheets were purchased locally and were 2 ft × 2 ft in dimension, shown in Fig. 13. As discussed previously, the prepreg composite sheet contains viscoelastic properties due to the resin contents, while all other sheets contain no resin or viscous content. The elastic materials all have a standard/uniform weave. Details of the elastic materials are shown in Table 3.

The composite prepreg material provided by Boeing came with manufacturer density and thickness specifications (refer Table 4). The elastic materials came with no specifications and thus required density calculations and thickness measurements, shown in Table 3.

The physical experiment process described previously was repeated twice for each material type for a total of eight trials. In the case of the composite sheet, two different sheets were used, as it was determined that the combination of experimental movements and exposure to air may degrade the sheet, possibly providing poor results. The elastic materials all used the same sheet twice. Of the two trials for each material, the data collected from one trial were used for training purposes, while the other was used for testing purposes.

Estimation of sheet parameters was accomplished through the simulator, as previously discussed. The internal parameters were considered to be in one state throughout the simulation procedure, but research has shown that linear elastic woven fabrics show a variance in many of their internal parameters due to the anisotropic behaviors of the fabric [31]. As such, manual manipulation of the parameters from the initial internal parameters was needed to find an ideal starting range for the optimization. The initial parameters were found experimentally and were then run through the optimizer, providing lowest error results.

Tables 5 and 6 highlight the training and testing performance of the simulator in terms of the average error and maximum error for the training and testing data sets. Table 5 gives the results for the training of fabric materials and Table 6 gives the results for the testing of fabric material.

Figures 14, 15, and 16 show the configuration of stage 5 for the felt, cloth, and canvas sheets, respectively. In each figure, the upper image is the observed position photograph, the middle image is the mesh generated from the scanned point cloud data, and the lower image is an image of the simulation mesh at the end of the stage.

Figures 17 and 18 depict the initial and final positions for the composite sheet trials respectively. Within the initial stage shown in Fig. 17, consistency can be seen between the photo, observed data, and simulation data.

Within the final stage, there was consistent error across all material types, however this error was higher than all other stages indicating that this position was the most difficult for the simulator to emulate. Figure 18 depicts this position and displays results from the observed data and simulated data.

The error from the training data is shown in Table 7. The optimal parameters were a bending stiffness of 4.77682 × 10⁷ N/m² and a shear stiffness of 2942.93 N/m. The optimization process began with initial values based on the data sheet supplied by Boeing in Table 8.

The results obtained from the material parameter estimation are then utilized to develop a simulation system that can predict the behavior of the material under varied constraints. This section focuses on one of the potential applications of such a high-fidelity simulator. In composite prepreg sheet layups, as discussed earlier, it is important for the sheet to be held at appropriate locations such that it does not result in any potential defects during the layup. A sheet simulation system that can predict the appropriate locations for grasping the sheet can be instrumental in process planning. The key elements of such a simulation system would be (1) a sheet simulation model generated using the estimated material parameters and (2) a real-time sheet tracking system that can generate a mesh of the sheet in real-time at high frame rates, and that can then be used for comparison with simulated result. To demonstrate the feasibility of such a system, an experiment was designed using the carbon fiber reinforced epoxy sheet provided by Boeing Inc.

The methodology for material parameter estimation proposed earlier was used to generate an appropriate model for the Boeing composite sheet. The simulation system was built using the VegaFEM library [9]. Figure 19 shows the block diagram of the proposed real-time sheet tracking and simulation process.

The dimensions of the composite sheet used in this experiment were 4 ft × 3 ft, which required a multi-camera system to track the entire sheet. The real-time sheet tracking system proposed in this study consists of three Realsense D415 sensors. The entire sheet is captured by fusing the red-green-blue-D feed from all the three cameras. Color-based filtering techniques are employed to achieve filtering of the prepreg sheet from rest of the scene. Re-sampling is performed on resulting points to obtain a uniform distribution of the filtered points. The face normals are then recomputed using these points. Surface reconstruction is performed in scale space by implementing advancing front surface reconstruction [32]. A scale space describes the point set at a dynamic scale, and this additional dimension allows us to control the degree of smoothness required for the reconstruction [33]. After post-processing, a surface mesh with around 6000 triangles was obtained.

To test the system performance, we constrained the composite sheet in four different states. Initially, the sheet is constrained at four grasping locations. A human and two 7-DOF robotic manipulators KUKA iiwa R7 are used to apply the fixed constraints. The real-time sheet tracking system captures the mesh and grasping locations for the corresponding state. In the next step, one of the grasping locations is released. The sheet is allowed to settle to its minimum energy state and the real-time sheet tracking system captures the sheet in the respective state. This exercise is repeated for three other states by changing the grasping locations and capturing the sheet behavior. These data are then compared to the predictions generated by the simulation system. Figure 20 contrasts the four states in which the composite sheet was suspended and the corresponding simulator predictions.

The proposed simulation and real-time sheet tracking system can play a pivotal role in predicting optimal grasping locations for the draping of composite sheet. These grasping locations can then be used to deploy a human robot collaborative cell where the robots can aid the human in the draping process by holding the sheet appropriately [34]. Furthermore, the sheet tracking data can act as a rectifying input for the grasping location in case of sub-optimal simulator predictions. Such a conjunctive system can aid in streamlining the prepreg layup process and achieve an overall higher degree of automation.

In this work, we presented a methodology to predict material parameters to generate an energy-based model of various materials. We showcased how these models aid in simulating the material behavior. In addition, we discussed the implementation and development of a high-speed thin-shell simulator based on the VegaFEM library. The procedure followed for data collection for model training and testing is highlighted as well. The inputs to the parameter model are elaborated along with the boundary conditions pertaining to fixed constraints and external forces. Validation of the developed model for various test cases is presented and discussed in depth.

This work has provided detailed insights into the simulation of multi-material composite material-based systems, particularly prepreg sheets along with other non-viscoelastic materials. Emphasis has been laid on predicting behavior of prepreg composite sheets under different types of constraints and on simulating how the sheet will behave under these conditions. The prediction capabilities of our system open new avenues in the domain of automated composite layup processes, particularly in planning for the development of completely autonomous and human–robot collaborative cells. This system can be employed to plan and execute a defect-free layup using robotic assistants, using the simulator to generate optimal locations to grasp the prepreg sheets during a layup process. In the future, we aim to explore the capability of employing this simulator in trajectory and motion planning of collaborative robots mainly used for manipulation of prepreg sheets, thus enabling the automation of composite layup tasks.

This work is supported in part by National Science Foundation Grant No. 1925084. Opinions expressed are those of the authors and do not necessarily reflect opinions of the sponsors.

There are no conflicts of interest.

The authors attest that all data for this study are included in the paper. Data provided by a third party are listed in Acknowledgment.