Towards a digital twin for mitigating void formation during debulking of autoclave composite parts

Journal Pre-proofs Towards a Digital Twin for Mitigating Void Formation in Autoclave Composite Parts Guillaume Seon, Yuri Nikishkov, Andrew Makeev, Lauren Ferguson PII: DOI: Reference:

S0013-7944(19)30655-1 https://doi.org/10.1016/j.engfracmech.2019.106792 EFM 106792

To appear in:

Engineering Fracture Mechanics

Received Date: Revised Date: Accepted Date:

20 May 2019 18 November 2019 22 November 2019

Please cite this article as: Seon, G., Nikishkov, Y., Makeev, A., Ferguson, L., Towards a Digital Twin for Mitigating Void Formation in Autoclave Composite Parts, Engineering Fracture Mechanics (2019), doi: https:// doi.org/10.1016/j.engfracmech.2019.106792

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

© 2019 Published by Elsevier Ltd.

Towards a Digital Twin for Mitigating Void Formation in Autoclave Composite Parts Guillaume Seon1, Yuri Nikishkov2, Andrew Makeev1 and Lauren Ferguson3 1

Department of Mechanical and Aerospace Engineering, University of Texas at Arlington, Arlington, Texas, U.S.A. [email protected] 2 3

Numerical Technology Company LLC, Dallas, Texas, U.S.A.

Air Force Research Lab AFRL/RXCC, Wright Patterson Air Force Base, Ohio, U.S.A

Abstract High-performance polymeric composites are increasingly used in the design of aircraft structural components; however, susceptibility to manufacturing irregularities, including porosity/voids, remains a primary challenge delaying the implementation of advanced composites in modern aircraft. Voids may be precursors to structural damage significantly affecting structural integrity and remaining useful life of the aircraft. High-performance aerospace composite parts are commonly manufactured from resin-saturated pre-impregnated plies of uncured material laid-up over a rigid tool and consolidated and cured in an autoclave. In many applications, especially in thick and curved composite sections common in aerospace designs, autoclave consolidation alone is not sufficient to remove the air, or bulk, that might have been entrapped within the laminate during the layup process. Vacuum consolidation, or “debulking”, has been a standard practice extensively used by manufacturers to reduce the amount of bulk prior to autoclave curing. The debulking process typically takes long time and requires intensive manual operations. The underlying physical principles governing the formation and evolution of voids during these early stages of the manufacturing process are not yet well understood, and controlling parameters remain largely determined empirically or based on prior experience. Driven by the need to improve such understanding, there has been a recent interest in using newly available high-fidelity nondestructive inspection (NDI) techniques, such as X-ray Computed Tomography, for in situ observations of composites internal structure during the early stages of manufacturing. In situ observations are allowing researchers to identify the driving mechanisms involved during defect formation, and develop improved predictive models. The objective of this work is to show that high-fidelity NDI data and new developments in numerical modeling might be combined into the creation of a Digital Twin for mitigation of void formation in composite parts. In particular, X-ray CT data is used to extract bulk content and distribution in uncured carbon-epoxy curved-beam specimens after manual layup, and such information is transferred into a finite element (FE) model for simulation of debulking. The FE model uses a fracture-based approach that relies on porepressure cohesive zone modeling recently proposed for the discrete representation of entrapped air pockets in uncured resin-saturated prepregs. Preliminary results support a promising ability of the digital twin concept for optimizing the debulking process of autoclave composites towards mitigating void formation.

Keywords Digital Twin; Autoclave Composites; Void Formation; Debulking; Vaccum Consolidation; Pore-pressure Cohesive Zone 1. Introduction High-performance polymeric matrix composites are increasingly used in the design of aircraft structural components and their application has progressed throughout the past five decades from secondary to primary structures. Composite aircraft start dominating the commercial aircraft market with the production of Boeing 787 and Airbus 350 and expected upgrades of legacy systems. On the United States DoD (Department of Defense) application side Lockheed Martin has been ordered production of more than 3,000 F-35 aircraft. Similarly, the US Army and Rotorcraft Industry are facing the Future Vertical Lift (FVL) aviation challenge to replace more than 6,300 military vertical lift aircraft [1]. However, susceptibility of composite structures to manufacturing irregularities, including porosity/voids, remains a primary challenge delaying the implementation of advanced composites in modern aircraft. Voids may be precursors to structural damage and significantly affect structural integrity and remaining useful life [2-7]. Recent advances in high-fidelity non-destructive inspection of composites, such as X-ray Computed Tomography (CT), have also shown that individual small voids at critical locations could significantly affect structural performance, despite low porosity volume content typically less than 1% in aircraft parts [5-7]. For example, subsurface CT scan inspections and finite element stress analysis in [5,6] revealed that reduction of strength and fatigue life in overall low porosity content curved beam specimens could be related to shape, size, and location of individual critical voids. Similar conclusions are presented in [7] based on statistical considerations. Continuous improvement of the capabilities of industrial CT scanners as well as ongoing development of new methods for 3D volume reconstruction, such as limited angle X-ray computed tomography [8], increase the confidence in the ability to break through the current limits of X-ray CT and enable high-fidelity NDI of larger composite structures. Once the condition of a fabricated part is known at the appropriate scale and fidelity level, and the structural response can be accurately captured through numerical models that integrate part condition, a fundamental shift from statistics-based to condition-based structural substantiation could be achieved. Such transition will be central to mitigating the effects of defects through the deployment of the Digital Twin paradigm as advocated in [9,10], with promising outcomes in terms of lowering maintenance costs and increasing reliability and safety of aircraft. High quality of polymer-matrix carbon fiber-reinforced composite aircraft structures is typically achieved through autoclave processing of resin-saturated thermoset pre-impregnated (prepreg) plies. Prepreg plies are commonly laid-up over a rigid tool and consolidated in the autoclave oven under controlled temperature and pressure, as defined by the curing cycle. The use of an autoclave allows application of high external pressure on the open surface of the laminate (6-7 bars), which is usually encapsulated in a vacuum bag. High external autoclave pressure, combined with the application of elevated temperature, results in an increase in resin hydrostatic pressure. Resin hydrostatic pressure enables moistures, or other volatiles, to remain in solution and be removed from the laminate as the excess resin bleeds out during compaction, leading to an overall low porosity content. Creating resin hydrostatic pressure is crucial to achieving good part quality; however, it is generally not sufficient to evacuate the large air pockets that may be 2

entrapped at ply interface during layup prior to autoclave consolidation [11-13]. Air is easily entrapped during the layup process due to tackiness of the pliable prepreg material, which might be regarded as a pressure-sensitive adhesive material, coupled with the fact that out-of-plane air permeability of typical resin-saturated raw thermoset autoclave prepregs is very low. Debulking, or vacuum consolidation, has been a common practice extensively used by aircraft manufacturers to reduce the amount of “bulk” or air entrapped during the lay-up of composite parts prior to curing. The debulking process can span several hours or days and it requires intensive manual interventions. Multiple debulking operations are typically required during layup of complex composite parts, such as composite helicopter components, which include thick and curved sections. Tight radii combined with thick sections inherently result in non-uniform consolidation pressure, making curved composites especially prone to formation of large voids, or pockets of entrapped air, at ply interface. As highlighted above, controlling porosity and void content in composite aircraft parts can be a complex and labor-extensive task that drives high costs and is time consuming. Prepreg manufacturers typically provide recommendations on processing and curing their materials, such as recommended temperature and pressure cycles for autoclave consolidation. While the recommended cycles might not be optimal for a given part geometry and tooling configuration, it is at least a starting point for limiting porosity development during curing. It is also worth mentioning that analysis models, such as approaches based on the Kardos model which considers stability and growth of small spherical voids in a viscous resin [14], have shown potential on providing guidance in the selection of optimized temperature and pressure profiles to reduce volatiles-related porosity [15-17]. Homogenized air extraction models based on Darcy’s Law have also been used to predict air evacuation of partially impregnated panels in out-of-autoclave manufacturing [18,19]. On the other hand, and to the best of the authors’ knowledge, debulking in autoclave resin-saturated prepregs has received much less attention and parameters are mostly determined empirically in current industry practices. The efficiency of debulking is bound to be related to many factors and heavily depend on operator experience. Initial distribution of entrapped air during manual layup, presence and geometry of air pathways, proximity to laminate edges and vacuum ports, conflicting effects from application of bag pressure, frequency and duration of vacuum application, leaks in the vacuum bag system, as well as many other conditioning and environment factors, can potentially affect the efficiency of the debulking process. Accordingly, numerous best practices, such as number and positioning of vacuum ports, choice of bagging material, vacuum pump dimensioning, use and placement of breather cloths, release films and perforated films, must be established when new part geometries or new materials are introduced. Advanced vacuum bagging practices, including vacuum bag pleating [20] and double-bagging techniques [21], are also worth mentioning. Driven by the need to improve current understanding of the manufacturing process of composites, and in particular the mechanisms driving part quality, there has been recently a growing interest in employing advanced NDI techniques such X-Ray CT technology for inspecting and monitoring the evolution of composite microstructure during the manufacturing cycle. In [22], Centea and Hubert interrupted the curing cycle of woven out-of-autoclave prepregs and used XRay CT NDI to inspect and visualize the evolution of resin impregnation at different stages during curing. In [23], Torres et al used similar interrupted tests to quantify and analysis void distribution using X-Ray CT measurements throughout the curing cycle of carbon/epoxy IM7/M56 out-ofautoclave prepreg laminates. Such studies directly foster the development of improved models for process simulation. For instance, in [13] the authors of this work performed X-Ray CT inspections of unidirectional curved-beam specimens made of raw carbon-epoxy IM7/8552 prepregs before 3

and after vacuum bag debulking at room temperature. Data suggested that the amount of initial bulk entrapped at ply interface and its evolution during debulking are determining factors in the formation of voids in autoclave contoured composites. These observations prompted the development of a new simulation concept based on the discrete modeling of the voids entrapped at ply interface during debulking [13]. Very recently, continuous reduction of X-Ray CT acquisition times has allowed researchers to performed in-situ X-Ray monitoring at 1μm voxel resolution with 90 seconds acquisition time using state-of-the art synchrotron X-Ray facilities for dynamic evaluation of void content during pressure-assisted resin impregnation of unidirectional dry fiber preforms [24]. Reference [24] showed that the data could enable improved prediction of material permeability. Such technological advances support extending the Digital Twin concept to the very early life stage of the aircraft, including the manufacturing process of its composite parts. In particular, the Digital Twin might be a key enabling technology to mitigate formation of voids in autoclave composites. As mentioned previously, air entrapment at ply interfaces during layup is most likely a stochastic process, as it depends on too many uncertain factors, including operator skill and experience. If the initial distribution of bulk after layup of a given part can be rapidly evaluated and such information transferred into an efficient physically-based predictive model, a Digital Twin, such model could be used to provide guidance on the best way to debulk the part. After debulking, the Digital Twin could be updated based on the new part condition and interrogated again for selection of the best curing parameters to further mitigate part irregularities. Contributing to such development is the main objective of this work. In particular, this work focuses on simulation of air removal during debulking and mitigation of the large voids that might have been entrapped during layup. A Digital Twin model is proposed that relies on X-Ray CT data of the composite laminate obtained after layup and the discrete modeling of air pockets using a new FEM-based approach first introduced by the authors in [13]. The discrete model relies on cohesive elements enriched with pore-pressure degree of freedom inserted at ply interface to represent large pockets of entrapped air. Such approach uses simulation concepts originally developed for the analysis of hydraulic fracturing in geomechanics problems. In the case of composite manufacturing, air pressure within the cohesive gap is coupled with the interfacial debonding process where cohesive forces are a representation of tackiness of the uncured prepreg material. Air pressure gradients might also lead to tangential air flow within the pore-pressure cohesive interface, which is representative of the slow in-plane air seepage process that occurs in fully-saturated prepregs from application of vacuum at part boundaries. The rest of the manuscript is organized as follows. Section 2 summarizes the main concepts and governing equations of the discrete pore-pressure cohesive approach. In Section 3, a procedure for automatic quantification of entrapped air pockets using image segmentation of X-Ray CT data is presented and illustrated for assessment of initial bulk content in uncured unidirectional curved beam specimens after manual layup. Creation of the Digital Twin with transfer of the entrapped air information into a FE model developed for Abaqus FEM commercial software is introduced in Section 4 and illustrated using the data from curved beam specimens. In this preliminary work, the methodology is demonstrated in two dimensions only and several key properties can only be inferred from reasonable engineering considerations in the absence of available test data. In Section 5, further modeling concepts are introduced and the Digital Twin is used to evaluate efficiency of void removal under different possible debulking scenarios. Concluding remarks are presented in Section 6.

4

2. Pore-Pressure Cohesive Zone Modeling 2.1 Concept and Governing Equations The main concept for discrete modeling of an air pocket entrapped at the interface between two impermeable layers of tacky material is illustrated in Figure 1. The reader might also refer to [13] for additional details. As shown in Figure 1, an air-filled pocket is represented by a fully broken cohesive zone where transverse fluid flow and air injection is permitted. At air pocket boundaries, the cohesive zone is partially broken and surface tractions are nonzero due to tackiness of the interface between the two layers, which resists the reaction forces due to distribution of air pressure.

Figure 1. Pore pressure cohesive zone modeling for an air pocket entrapped at the interface between two impermeable media. The relationships between the deformation of the surrounding material, the pressure pf inside the pocket of entrapped air and the cohesive law representing the bonding and debonding process at the interface is found by solving a system of coupled constitutive equations which are summarized below for reference. The reader might refer to Refs. [25-27] for more background and recent developments on the derivation and implementation of the constitutive equations with applications to fluid-driven fracture problems in rock mechanics. Assuming that the air entrapped between the two layers of material behaves as an incompressible fluid with Newtonian rheology, the tangential flow within the broken cohesive zone is governed by the equation of mass conservation: ∂𝑤 + ∇.𝒒 = 𝑄(𝑡)𝛿(𝒙,𝒙0) ∂𝑡

[1]

where q is the fluid flux of the tangential flow, w is the fracture opening, Q is the fluid injection rate and 𝛿(𝒙,𝒙0) is a spatial Dirac function equal to unity at the point of fluid injection. The flux of tangential flow q can be related to the gradient of fluid pressure ∇𝑝𝑓 along the cohesive zone according to the Poiseuille’s law: 𝒒𝑤 = ― 𝑘𝑡∇𝑝𝑓

[2]

5

where kt is the tangential permeability to air flow, which might be defined as a function of the gap opening w according to Reynold’s equation: 𝑤3 𝑘𝑡 = 12𝜇

[3]

with  the fluid viscosity. Continuity of displacement on the crack faces imposes: 𝑤 = 𝒖.𝒏

[4]

Where u is the displacement field in the surrounding material and n is the local unit normal vector to the crack faces. The fracture process at the crack tip is modeled using a traction-separation law that relates the effective traction tensor' to the separation  across the cohesive surfaces: 𝝉′ =

∂𝜙 ∂𝛿

[5]

 is a cohesive potential function that describes the bonding and debonding process at the interface between the two layers of surrounding material, including description of the tacky behavior. The effective cohesive tractions can be related to the stress tensor 𝝈 in the surrounding material from the system of equations resulting from force equilibrium at the crack faces: 𝝈.𝒏 + = ―𝝈.𝒏 ― = 𝝉′ ― 𝑝𝑓𝒏

[6]

Deformation in the surrounding media is described by their respective mechanical constitutive relations. For an air pocket entrapped between two layers of the same material, only one set of constitutive equations is needed. In the general framework of hyperelasticity, which may be used to represent deformation of the uncured prepreg material [28,29], the constitutive relation can be defined by relating the Cauchy stress tensor 𝝈 to a strain energy potential function 𝝍 and the hydrostratic pressure p: 2 ∂𝝍 𝝈 = 𝑭 𝑭𝑇 ― 𝑝𝑰 𝐽 ∂𝑪

[7]

where F is the deformation gradient, C is the right Cauchy-Green deformation tensor such as 𝑪 = 𝑭.𝑭𝑇, J = det(F) and I is the second-order identity tensor. Finally, the stress tensor 𝝈 is also related to the deformation of the material and the applied body forces b through the equation of conservation of momentum: ∇.𝝈 + 𝜌𝒃 ― 𝜌𝒖 = 0

[8]

Where  is the density of the material and 𝒖 is the acceleration vector.

6

The uncured prepreg is assumed to be fully resin-saturated and impermeable to air flow in the out-of-plane direction. Such assumptions are representative of the behavior of uncured autoclave carbon/epoxy prepregs during vacuum consolidation or debulking. The Abaqus/Standard coupled pore fluid flow and stress analysis solver for transient analysis is used to solve the set of coupled constitutive Equations (1-8) in a weak form in the time domain [30]. 2.2 FE Implementation The pore-pressure cohesive model is implemented using zero-thickness cohesive pore pressure elements available in the Abaqus/Standard element library [30]. In this work, six-node twodimensional (2D) pore pressure cohesive elements COH2D4P are used, as illustrated in Figure 2. As shown in Figure 2, the element includes four nodes with displacement and pore pressure degree of freedom and two additional pore pressure nodes defined in the gap interior.

Figure 2. Six-node COH2D4P elements available in Abaqus/Standard. The cohesive traction-separation law used in this work is implemented within an Abaqus/Standard User Material subroutine UMAT. The cohesive model is based on the mixedmode bilinear traction-separation laws introduced by Turon et al for simulation of interlaminar delamination failure in composites under mixed-mode loading [31]. It is worth noting that Turon et al bi-linear traction-separation cohesive model is available in Abaqus/Standard built-in material models [30]. However, the available built-in models for definition of pore-pressure cohesive properties do not allow for the modeling of tackiness. 2.3 Modeling Tackiness The tackiness, or pressure sensitive behavior, is included in the cohesive constitutive model by the introduction of a history state variable (SDV1) within the UMAT subroutine. The state variable is set to unity when a previously open cohesive element is closed under compressive loading and set to zero when a previously closed element is re-opened. The state variable is used within UMAT to fully restore the cohesive properties of broken cohesive elements that undergo compressive loading. The state variable is also passed to user subroutine USDFLD to assign solution dependent properties to the tangential air flow permeability kt and prevent tangential flow in cohesive elements under compressive loading until further tensile loading occurs. The flowchart diagram for the update of history state variable SDV1 and the interactions between the UMAT subroutine, USDFLD subroutine and Abaqus/Standard pore fluid flow and analysis solver is summarized in Figure 3.

7

Figure 3. Flowchart for implementation of the effects of tackiness.

3. Quantification of Initial Bulk Content in Curved Beam Specimens 3.1 Curved Beam Specimens An automatic procedure was developed for extraction of the 3D attributes of pockets of entrapped air at ply interface in unidirectional curved beam specimens using X-ray CT inspection. Curved beam specimens were fabricated by manual layup of 36 unidirectional plies of raw IM7/8552 prepreg material on a special 3D-printed tool, as shown in Figure 4. The tool was 3Dprinted from a low density thermoplastic to limit obstruction during X-Ray CT scanning of the specimens. Teflon release tape was applied on the surface of the tool to ensure limited friction at the part/tool interface during consolidation and eventually facilitate separation from the tool if the part is cured.

Figure 4. a) 3D-printed tool and b) uncured 36-plies unidirectional IM7/885 curved beam specimen after manual layup on the tool. After manual layup, the curved-beam specimen and tool assembly was enclosed in a vacuum bag and mounted on the rotary stage of the X-ray CT facility for in situ CT scanning of the specimen before and after debulking. A North Star Imaging X5000 industrial CT system with a 8

225 kV micro-focus X-ray tube and Varian 4030E series flat panel detector was used in this work. Figure 5 shows the setup for in situ X-ray CT scanning.

Figure 5. Bagged curved beam specimen and tool assembly installed in the X-ray CT facility for in situ inspection It is worth noting that only CT data measured prior to debulking is used in this work. However, the setup depicted in Figure 5 allows for in situ CT NDI during vacuum consolidation. At any time during debulking, which typically spans over several hours, a valve in the tubing leading to the vacuum pump can be closed and the specimen can be CT scanned under sealed conditions corresponding to the current debulked state. Comparison and quantification of bulk distribution before, during and after debulking is outside the scope of this work, but will be the object of a future publication by the authors. Figures 6 and 9 show examples of CT slice images obtained in curved beam specimens after layup. 3.2 Methodology for Void Quantification The methodology for quantification of the bulk content in curved beam specimens follows feature extraction principles common in modern computer vision algorithms. In this paradigm, a density field for three-dimensional cloud of points provided by the CT reconstruction is analyzed using a set of fixed assumptions about the specimen, such as overall structure, coordinate system, number of plies etc., and a set of varying spatial properties, such as specimen thickness, ply thicknesses and void location and sizes. Fixed assumptions are usually easy to define but difficult to deduce from the CT scan analysis alone, and have major global impact on extracted variable properties; while variable properties may have complex spatial distributions but only local effect on the extracted feature set. Due to the presence of noise in CT scans and potential foreign material 9

inside specimens, additional smoothing and outlier identification is required after the feature extraction is completed. Input data and user-specified fixed assumptions used in the algorithm are detailed in Table 1: Table 1. Input data and user-specified fixed assumptions used in the algorithm for void quantification Data or Parameter

Type

CT scan slice images approximately Input Data along the width of the specimen Angle beam radius Maximum beam thickness Number of plies

Fixed assumptions

Comments Slice images in a desired direction can be exported by most CT scan inspection software Parameters characterizing angle beam geometry

Fixed assumption Image coordinate system for the first and Corner section coordinate system in based on last image using a pre-defined ordering CT images manual user with an example in Figure 6 input Voxel size of CT reconstruction Figure 7 shows examples of air and tool Density threshold for border CT scan border thresholds defined from density (specimen versus tool) parameters profiles using the approach detailed in Density threshold for voids (air [12] versus material) Number of extraction points Start point offsets of density profiles Widths of moving average for the specimen shape and plies Number of outliers

Feature extraction algorithm settings

User-defined discretization and algorithm parameters that may be changed depending on CT data quality

Void merge parameters

An example of CT slice image with definition of the coordinate system used in the extraction algorithm is shown in Figure 6. Figure 7 shows an example of threshold selection for detection of the boundary between the tool and the part and the part material and air pockets, based on the approach introduced in [12].

10

Figure 6. An example of CT slice image in a unidirectional curved beam specimen after layup with overlaid coordinate system definition.

Figure 7. Desnisty profile of a CT slice. Air, specimen, tool, Teflon and voids are indicated on the profile, and threshold selections for air and border are shown. The algorithm used for extraction of the void information in a 2D CT slice image of a curved beam specimen is detailed below. All steps are repeated for each CT slice image along the width of the curved beam specimens. 1) Interpolate the arm inner surface points specified by user and find coordinate system for this corner section, including center of the corner and corner angle. 2) March over the length of the arms and the corner section and analyze density profiles over the thickness of the beam. Approximate surface location is determined by the interpolated coordinate system, angle beam geometric parameters and start point offsets. Figure 6 demonstrates an example profile, from which the algorithm defines thickness of the beam and void thickness and location, if any. 3) Remove specified number of outliers and apply smoothing to the raw data from step 2 over the length of the beam using the moving average method. 11

4) Assuming that voids are located at ply interfaces only, estimate interfaces that contain voids along the length of the beam by allocating equal ply thicknesses between the voids to match total thickness. If contrast in CT scan allows detection of different plies (uncommon for Carbon/Epoxy composites), this step can be improved to detect ply thickness variations. 5) Combine void thickness measurements from adjacent profiles to make 2D voids (thicknessbeam length). 6) After analysis of all slice images is completed, 2D voids adjacent in width are combined into three-dimensional voids. Figure 8a shows an example result of the algorithm with thickness profiles and total void content plotted over the curved beam length for a sample CT slice in an 8-ply curved beam specimen. Curved beam length is defined in the local curvilinear coordinates illustrated in Figure 6. Values measured by the algorithm are shown by points and solid lines show interpolated data after smoothing is applied. Figure 8b shows smoothed ply boundaries that indicate presence of voids at the interfaces. Note that in both plots thickness is not to scale with respect to beam length.

Figure 8. Total shape thickness and total thickness of voids over beam length (a) and ply boundaries with voids at interfaces (b). Note that thickness is not to scale with length. 3.3 Results of Initial Bulk Measurements For verification and easier interpretation, the void information was used to build twodimensional FE meshes from CT slice images of 36-plies unidirectional IM7/8552 curved-beam specimens after manual layup. The 2D mesh was generated using a mesh morphing technique similar to the approach detailed by the authors of this work in [32]. Such method allows automatic generation of segmented FE mesh for complex geometry. 2D mesh were generated automatically using Abaqus CAE Python scripting interface. Figure 9 shows an example of 2D visualization mesh and compares it with the CT data from which it is automatically generated. Large pockets of air entrapped at ply interface during manual layup are illustrated in Figure 9.

12

Figure 9. Comparison of a) CT data and b) corresponding 2D visualization mesh generated using mesh morphing As illustrated in Figure 9, bulk distribution in the specimen is captured with good fidelity by the void extraction algorithm. It is worth noting that in the mesh morphing approach used for data visualization, ply boundaries are simply morphed to match the smooth lines generated from CT data using the void extraction algorithm. In this case, the gap between the plies is a visual representation of the bulk content; however, it lacks physical meaning. For instance, if external compressive loading were to be applied on this mesh, the gap would collapse with no resistance. In reality, external loading should result in increased air pressure within the gap, which in return resists to gap closing within some limits related to leak-off properties and in-plane permeability. Discretization of air pressure distribution within the ply gaps using pore-pressure cohesive elements and implementation of the constitutive model describing the interaction between air pressure and deformation of the surrounding material are key ingredients to the discrete modeling approach introduced previously in Section 2. 4. Creation of a Digital Twin The Digital Twin proposed in this work is a high-fidelity FE representation of the curved-beam specimen in its as-manufactured state after manual layup, including the initial distribution of entrapped air pockets, or bulk. The digital twin is created using inputs from the void extraction algorithm after processing of CT slice images and relies on the optimization of air injection FE analysis. In the air injection analysis, a node-based inbound seepage flow condition is prescribed at discrete locations corresponding to the location of the center of voids detected in the CT scans. An illustration of air injection analysis is shown in Figure 10 for a simple 2D two-ply laminate. Air is injected at constant rate at the center of the laminate within the pore-pressure cohesive layer, which results in progression of fluid-driven fracture as shown in Figure 10. The air fracture process is related to the cohesive properties of the pore-pressure interface, which can be physically interpreted as a representation of the tackiness of the raw prepreg material. More details are provided in reference [13].

13

Figure 10. Comparison of a) CT data and b) corresponding 2D visualization mesh generated using mesh morphing A flowchart diagram for creation of the Digital Twin of the curved-beam specimen using air injection analysis is illustrated in Figure 11. Three-dimensional void attributes extracted from a CT slice image are used as inputs for Abaqus CAE pre-processing of a 2D FE model representative of the corresponding cross section of the curved-beam specimen. Inputs include void centers, void thickness, void length, ply number and nominal void-free cross section thickness. First, a generic FE mesh is generated using the nominal beam thickness. The FE model is segmented and zerothickness of pore-pressure cohesive elements are inserted for all ply interfaces where a void has been detected. Nodes at void center locations and cohesive elements spanning the void length are selected and referenced in respective node and element sets. Cohesive elements assigned to void sets are defined as initially broken by initializing the state variable SDV1 described in Section 2.3 to zero prior to FE analysis. Air injection rates at all discrete void center locations are set to an initial value and air injection analysis is carried out. At the end of the air injection analysis, void attributes are extracted from the FEM and a feedback loop can be used to optimize the injection rates at discrete void locations such as simulated void attributes match the CT data.

14

Figure 11. Flowchart diagram for creation of the Digital Twin This work was limited to validation of the Digital Twin concept, and the optimization loop was not used. It was found that a single optimization variable Qglobal could be used to define all injection rates Q[i] at discrete void centers using Equation 9 to simplify the process with acceptable accuracy. This approximation relies on the physical assumption that air pressure within the voids, and consequently air injection rates, should be proportional to the voids area A[i], which are measured inputs. 𝑄[𝒊] = 𝐴[𝑖] × 𝑄𝑔𝑙𝑜𝑏𝑎𝑙

[9]

In this work, a trial-and-error approach was used to determine global optimization variable Qglobal such as correlation between the Digital Twin and the physical model is visually acceptable, as illustrated in Figure 11. It worth noting that automatic optimization could be easily implemented since the process of creating the FE model and running the air injection analysis was automated through Python scripting within Abaqus CAE interface. Nevertheless, it was interesting to observe that some voids would merge or redistribute during the simulation, such as the dynamic equilibrium between pore fluid pressure, cohesive forces and tangential flow is satisfied. This can be interpreted as a “physically-based smoothing” of the input data provided by the void extraction algorithm which might be affected by measurement noise, discretization parameters and other uncertainties. It is also worth noting that creation of the Digital Twin illustrated in this work relies on prior knowledge of the material constitutive behavior and tackiness properties. In this preliminary work, several of these constitutive properties were assumed based on reasonable engineering 15

considerations, but not explicitly measured or taken from references that might exist in the literature. The reader might refer to [13] for some details on preliminary material parameters that were also used in this work. 5. Validation of the Concept 5.1 FE model A 2D FE Digital Twin was created from a selected cross-sectional CT image of an uncured curved-beam specimen after manual layup. This model was used for simulation of the debulking process with the objectives to improve our understanding of the mechanisms driving air removal and evaluate the ability to use the model to select optimum debulking parameters for mitigating air entrapment during layup. It is worth noting that the 2D FE model is not truly representative of the physical distribution of initial bulk in the specimen. For demonstration of the concept, one CT slice image was selected and the FE simulation considered plane strain approximation. Under these assumptions, the FE model is representative of a specimen with a very large width dimension and uniform void attributes along the width. In the physical specimen, void distribution is threedimensional; therefore 2D FE results presented in this work are not comparable with experimental results after debulking. The 2D model is obtained from air injection analysis where air-filled pockets within broken pore-pressure cohesive elements develop into the final distribution illustrated in Figure 12, which representative of bulk distribution in the 2D slice of the physical model. For visualization purposes, the deformation scale in Figure 12 has been multiplied by two.

Figure 12. FE Digital Twin with initial bulk distribution used for simulation of debulking The FE model shown in Figure 12 is used in a second analysis step representative of the debulking process. In this second step, a uniform distribution of external pressure is defined on the free open surface of the FE model, as illustrated. External pressure loading is used to represent the action of the vacuum bag, which might be assimilated to a thin membrane that can sustain a significant amount of stretch and possesses no bending stiffness. Displacements boundary conditions are defined at the sides of the FE model such as normal displacement is blocked. Pore pressure boundary conditions will also be defined at the sides of the specimen to simulate vacuum application. 16

As illustrated on Figure 12, five nodes are selected at the center of five distinct air pockets for better interpretation and comparison of the simulation results. In particular, evolution of pore pressure and maximum gap opening at these locations during debulking will be discussed in the next sections. 5.2 Mechanisms Driving Air Removal To better understand the mechanisms driving air removal during debulking, a second analysis step is run after air injection with only external pressure loading applied on the open face of the model, without imposing vacuum pore-pressure conditions on the sides of the model. Air injection analysis is carried out for 200 seconds with no external pressure applied. In the debulking step, air injection flow rates at void centers are set to zero and a ramp of external pressure is applied over 500 seconds then held constant for 2300 seconds. Figure 13 shows the FE results with the initial deformed shape after injection (t=200s) and the final deformed shape at the end of pressure hold (t=3000s). Figure 13 also shows the evolution of internal air pressure within the five monitored void locations as illustrated in Figure 12.

Figure 13. a) Initial bulk content, b) final bulk content and c) evolution of pore pressure at five void centers during application of external pressure only As shown in Figure 13, despite an overall increase in internal void air pressure during application of external pressure, there is almost no change in the distribution of air pockets before and after application of the external pressure. This is because the external pressure is simultaneously compensated at each time increment by development of internal pore pressure in the existing air pockets and development of compressive stresses within the ply material. Since the local pressure/stress distribution remains relatively uniform at the length scale of individual voids, there is no significant change in equilibrium for the initially stable air pockets, which therefore stay entrapped in the laminate throughout the simulation.

17

In fact, the development of interlaminar compressive stresses during application of external pressure will tend to prevent further air evacuation. Under compression, an increase of pore pressure within the cohesive gap will have to fully compensate the compressive stresses before leading to any increase of cohesive forces, which might eventually result in air transport through air-driven fracture. Additionally, compressive stresses might close the cohesive gaps that are initially open or partially open, and restore or increase the cohesive strength due to the tackiness and pressure sensitive adhesive behavior of the uncured material, preventing or deferring further air-driven fracture. It is worth noting that a simple model of tackiness is included in the simulation, where the interfacial cohesive strength is fully restored if the cohesive gap is closed under compression, as presented previously in Section 2.3. These considerations and observations using the 2D Digital Twin suggest that debulking under vacuum bag consolidation is not driven by the squeezing effects of the bag, which in fact has a detrimental effect on the efficiency of debulking. Results presented further in Section 5.3 will show that in-plane air seepage is the driving phenomenon of the debulking process under vacuum bag consolidation. As vacuum is applied, pressure gradients develop between the vacuum air pressure at the sides of the specimen and the internal pore pressure within the entrapped air pockets. Such in-plane pressure gradients might lead to tangential air flow, or “air seepage”, within the small air channels existing at partially bonded ply interfaces. Air seepage due to in-plane pressure gradients might compete with the effect of interlaminar compressive stresses resulting from the action of the vacuum bag membrane on the top surface of the laminate, which might tend to close the cohesive gap and prevent tangential flow. At this point, it is interesting to re-examine governing Equations 2 and 3 introduced previously for description of the tangential flow within the pore-pressure cohesive gap. As shown in Equation 3, tangential permeability is a function of the gap opening d according to Reynold’s model. From Equations 2 and 3, it is shown that the tangential flow depends on the square of the gap opening distance. Therefore, in presence of very small air channels at ply interface, the air seepage phenomenon will be slow, which is qualitatively consistent with debulking operations typically spanning several hours. Furthermore, using Abaqus built-in models [30], the gap opening d is defined as: 𝑑 = 〈𝑡𝑐𝑢𝑟𝑟 ― 𝑡𝑜𝑟𝑖𝑔 + 𝑔𝑖𝑛𝑖𝑡〉

[10]

where tcurr and torig are the current and initial cohesive element geometrical thickness, respectively; ginit is an initial gap opening property and 〈〉 is the McCauley bracket. Equation 10 allows interpretation of the driving mechanisms for tangential air seepage. In the the user-defined cohesive model implemented within UMAT and USFLD illustrated in Figure 3, a cohesive element is defined as closed when relation 11 is satisfied. 𝑡𝑐𝑢𝑟𝑟 ≤ 𝑡𝑜𝑟𝑖𝑔

[11]

Therefore, when 𝑡𝑜𝑟𝑖𝑔 ― 𝑔𝑖𝑛𝑖𝑡 ≤ 𝑡𝑐𝑢𝑟𝑟 ≤ 𝑡𝑜𝑟𝑖𝑔 the cohesive element is closed but the gap opening d defined in Equation 10 is still non-zero and tangential air flow might occur. This provides the framework for introducing a dependency of the air seepage phenomenon on interfacial 18

pressure, or a “pressure sensitive” behavior, which is expected for the uncured prepreg material. In particular, the term 𝑡𝑐𝑢𝑟𝑟 ― 𝑡𝑜𝑟𝑖𝑔 represents the cohesive normal separation 𝛿𝑛, which is related to the cohesive normal traction 𝜏𝑛 through the cohesive law. In the bilinear traction separation law currently implemented, this relation is defined as: 𝜏𝑛 = (1 ― 𝑑)𝐾𝛿𝑛 ― 𝑑𝐾〈 ― 𝛿𝑛〉

[12]

where d is the cohesive damage variable and K the cohesive penalty stiffness. In presence of compressive traction with 𝛿𝑛 negative and d = 1, Equation 12 can be simplified as: 𝛿𝑛 =

𝜏𝑛 𝐾

<0

[13]

As shown in Equation 13, the penalty stiffness parameter K provides a simple way to control the pressure sensitive behavior of the cohesive interface that drives the air seepage, by controlling the gap opening as a function of the interfacial pressure or negative traction. In particular, from Equations 10 and 13, a critical compressive traction, or “pressure”, Pcr for which d = 0 can be defined as: 𝑃𝑐𝑟 = 𝐾𝑔𝑖𝑛𝑖𝑡

(13)

When the pressure or negative traction across the interface is equal to or larger than Pcr, the gap opening d is zero and tangential air flow or air seepage is prevented. Penalty stiffness K and initial gap ginit are parameters that could be obtained experimentally from inverse characterization. Ideally, the external pressure applied during debulking should be as small as possible and less than the critical pressure Pcr. It is worth noting that a double-vacuum-bag (DVB) design as introduced by NASA Langley Research Center in [20] provides an ingenious technological solution to prevent the detrimental effect of external pressure loading on vacuum-driven air removal during debulking. 5.3 Mitigation of Void Formation 5.3.1 Loading Scenarios The modeling concepts discussed in the previous section are further illustrated and verified by considering four different loading scenarios for the 2D curved-beam Digital Twin, including three scenarios representative of open face vacuum consolidation over a rigid tool, as illustrated in Table 2. These loading scenarios are also used to evaluate the relevance of the current model for mitigation of void formation.

19

Table 2. Pore pressure boundary conditions and external pressure loading for the loading scenarios considered Loading Scenario

Imposed Pore Pressure Boundary Condition at Specimen Sides

Open-face External Mechanical Pressure Loading

A

no external pressure

B

no external pressure

C

D

20

For all loading scenarios, an initial air injection step at constant air injection flow rate is carried on first for 200 seconds to generate the initial bulk distribution in the model, as previously discussed. In a second step, air injection rate is set to zero and the analysis is run for 7200 seconds using the pore pressure boundary conditions and external pressure loading detailed in Table 2. Prior to the analysis, a uniform distribution of pore pressure corresponding to the atmospheric pressure is assigned in all p-cohesive layers as initial conditions. All pressure profiles shown in Table 2 are normalized by the amplitude of this initial pressure. In scenario A, the pore pressure at both left and right sides of the specimen is imposed equal to the initial atmospheric pressure throughout the entire analysis. No external pressure is applied on the open face of the specimen. This scenario might correspond to a situation where the specimen is left “at rest” for 7200 seconds. During this rest period, stresses that have developed within the laminate during air injection might redistribute and contribute to some relatively slow air seepage within partially open channels. This scenario could be physically interpreted as the stabilization of the internal bulk content after the layup operation. In scenario B, no external pressure is applied similarly to scenario A. On the other hand, pore pressure conditions at the sides of the specimen differ after 200 seconds. Pore pressure is first set to atmospheric pressure during the first 200 seconds of air injection, similar to case A, but then ramped down to zero over a period of 200 seconds to simulate vacuum application. Vacuum is then held for the rest of the analysis (7000 seconds). This scenario conceptually corresponds to vacuum debulking in a vacuum chamber or using the double-bag-vacuum design [21]. In scenario C, pore pressure conditions at the side of specimens are similar to scenario B, with a ramp down to vacuum pressure over 200 seconds and then vacuum hold until the end of the anaylsis. However, in this scenario, a ramp up of external mechanical pressure equal to the pore pressure drop is simultaneously applied on the open face of the specimen. This scenario simulates vacuum consolidation in a single vacuum bag, with the combined action of vacuum at the sides of the specimen and mechanical external pressure exerted by the bag membrane on the open face. After reaching full vacuum (zero pressure), vacuum pressure and external pressure loading are held until the end of the analysis. Scenario D is similar to scenario C for the first 400 seconds of simulation that include the 200 seconds of air injection followed by 200 seconds of progressive vacuum. The difference is that after full vacuum is reached, the external pressure is rapidly increased to seven times the vacuum pressure over a period of 500 seconds, then held constant at this value for the rest of the analysis. Physically, this scenario could correspond to a situation where autoclave pressure is applied following vacuum application. 5.3.2 Results FE results are presented in Table 3, which shows the final deformed shape of the laminate and time history plots for the maximum gap opening at voids #1 to #5 for the four different loading scenarios.

21

Table 3. Final deformed shape and gap opening time history Loading scenario

Final Deformed Shape

Gap Opening Time History

A

B

C

D

22

In scenario A, atmospheric pressure is imposed at the sides of the specimen, which means that the excess pore pressure that is enclosed within the entrapped air pockets might lead to in-plane pore pressure gradients within the pore-pressure cohesive layers and air seepage. The air seepage is dependent on the local balance of cohesive forces and on the positive value assigned to the initial gap property. In particular, significant air seepage is observed for void #3, which almost disappears in the final deformation plot. It is worth noting that similar void closing does not occur for voids #2 and #5, on the other hand. In fact, these voids become slightly larger throughout the analysis due to stress redistribution and air seepage “into” the voids, or a phenomenon of “air suction”. As shown in the deformed shape and time history diagrams of the maximum void opening, most of the voids are fully removed in scenario B, which simulates debulking in a vacuum chamber or double vacuum bag. Void #3 is rapidly evacuated first, followed by void #4, void #1 and void #5 which is about 95% removed. Void #2 is still present at the end of the analysis, but has been 75% removed and data suggests that it will be eventually fully removed if vacuum is held a bit longer. It is worth noting that scenario B is the most efficient scenario for void mitigation out of the four scenarios considered. This result was expected, since vacuum chamber and double-bagvacuum debulking are well known in the industry as the best practices (and most expensive) for debulking autoclave composite parts. It is very encouraging to demonstrate that the Digital Concept proposed in this work is consistent with the practical experience and provides means for further quantification. In scenario C, which simulates the application of vacuum using a single vacuum bag, a significant reduction of the overall bulk content is also illustrated, but not as efficiently as observed in scenario B. It can observed that voids #3 and #4 are fully evacuated after the 7200 second of vacuum hold. Void #2, on the other hand, is only about 50% removed. The shapes of the time history plots for the gap opening tend to suggest that all voids could be eventually removed, if the vacuum is held for a long enough period. Results for scenario D illustrates further the competition between air seepage from vacuum pressure at the edge of the specimen and the effects of external pressure that tend to close the cohesive gaps and prevent air evacuation. In this case, pressure is quickly ramped up shortly after vacuum. This results in the partial removal of several voids in the curved area. In particular, voids #3 and #4 that were fully evacuated in scenario B where no additional external pressure was applied, are still present in the laminate for scenario D. The shape of the gap opening history plots for voids #3 and #4 suggests that they will remain entrapped in the laminate, even if the vacuum pressure is held for a longer period. 6. Concluding Remarks Developing the composite process simulation tools capable of predicting formation of critical defects, including voids at ply interfaces, and integrating them with available tools for structural diagnostics and performance prediction are essential for closing the loop on effective control of the manufacturing process to prevent defect formation in autoclave composite parts. In this work, we showed that integrating high-fidelity in situ NDI data capturing the physical distribution of entrapped air during layup into a Digital Twin of the composite part could become a powerful tool providing guidance on the selection of optimized parameters increasing efficiency of the debulking process. Debulking of aircraft autoclave composite parts is currently a time and labor intensive 23

(but necessary) process, where controlling parameters are mostly determined empirically or based on experience with no guarantee that they are optimum or applicable when new materials and new geometries are introduced. The concept was illustrated using X-Ray CT data of uncured unidirectional carbon-epoxy IM7/8552 curved-beam specimens that were CT scanned after manual layup. An algorithm was proposed for automatic extraction of the geometric attributes of the large pockets of entrapped air observed at ply interface in the specimens. This information was used to create a simplified twodimensional finite element (2D FE) model for Abaqus FEM representative of the initial entrapped air, or bulk, distribution. This model is referred to in this work as the Digital Twin. For demonstration of the concept, a 2D plane strain approximation was used in the simplified FE simulation subsequently to be converted into a fully three-dimentional formulation including nonuniform bulk distribution through the specimen width. The FE model relies on pore-pressure cohesive layer modeling and uses simulation concepts originally developed for numerical analysis of fluid-driven fracture in rock mechanics. In particular, the Digital Twin model was created by injecting air at discrete locations within the porepressure cohesive layers, which results in the formation of air pockets between layers of ply material assumed to be impermeable to air flow. Parameters for air injection might be optimized such as the Digital Twin reflects the initial bulk content of the physical part with the highest fidelity. The model includes implementation of a user-defined cohesive constitutive model representative of the tackiness of the uncured prepreg. Tackiness plays a critical role in the formation of voids during layup of prepreg composites, as it is responsible for air entrapment at ply interface. Preliminary results were presented for a simplified 2D Digital Twin model used for simulation of air removal during debulking. The model allowed identification of air seepage at ply interface as the primary mechanism for air removal during vacuum consolidation. It was showed that external pressure loading exerted on the open-face on the part in single bag vacuum debulking might compete with the air seepage phenomenon and be detrimental to air removal. Different loading scenarios, including different conditions for vacuum application, were studied and showed that vacuum chamber debulking, or double-vacuum-bag design, would lead to the most efficient air removal process. These observations are well-known in the industry, gained from more than five decades of trial-and-error and refinement of best practices for manufacturing composite parts. However, to the best of the authors’ knowledge, this is the first time that a predictive numerical methodology with a potential for high-fidelity quantification and mitigation of void formation during the early stages of manufacturing is presented. There are, of course, many standing challenges before the Digital Twin concept can be effectively used for real-time decision-making and optimization of the manufacturing process. One technical challenge is the extension of the current two-dimensional FE model to three-dimensions, such as the physical distribution of entrapped air is fully captured. Another challenge is the development of accurate constitutive models and characterization of several key properties that were not available in this preliminary work, including improving the constitutive model of the uncured prepreg material and modeling tackiness more accurately. The expected dependency of material behavior and tackiness properties on conditioning and environmental variables, as well as visco-elastic behavior, were not considered in this work and are worth mentioning. The ability to predict evolution of the defects throughout the entire manufacturing process, including defects geometry in the final cured part, is also highly desirable for optimization of curing cycles.

24

Continuous improvements in high-fidelity NDI and computing hardware, combined with the injection of science into practical engineering analysis, will contribute towards realization of the Digital Twin paradigm. 7. Acknowledgments This work was partly funded by Numerical Technology Company LLC as part of the “Certification Modelling for Composites with Voids and Wrinkles for Engines and Structures” effort, which was performed under Contract FA8650-16-M-5048 to the Air force Research Laboratory (AFRL). The authors also thank Mr. Brian Shonkwiler at the University of Texas at Arlington for his assistance with specimen fabrication and CT scanning. 8. References 1. Makeev, A., Seon, G., Nikishkov, Y., Nguyen, D., Mathews, P. and Robeson, M., 2019. Analysis Methods for Improving Confidence in Material Qualification for Laminated Composites. Journal of the American Helicopter Society, 64(1), pp.1-13. 2. de Almeida, S.F.M. and Neto, Z.D.S.N., 1994. Effect of void content on the strength of composite laminates. Composite structures, 28(2), pp.139-148. 3. Wisnom, M.R., Reynolds, T. and Gwilliam, N., 1996. Reduction in interlaminar shear strength by discrete and distributed voids. Composites Science and Technology, 56(1), pp.93-101. 4. Chambers, A.R., Earl, J.S., Squires, C.A. and Suhot, M.A., 2006. The effect of voids on the flexural fatigue performance of unidirectional carbon fibre composites developed for wind turbine applications. International journal of fatigue, 28(10), pp.1389-1398. 5. Seon, G., Makeev, A., Nikishkov, Y. and Lee, E., 2013. Effects of defects on interlaminar tensile fatigue behavior of carbon/epoxy composites. Composites Science and Technology, 89, pp.194-201. 6. Nikishkov, Y., Seon, G. and Makeev, A., 2014. Structural analysis of composites with porosity defects based on X-ray computed tomography. Journal of Composite Materials, 48(17), pp.2131-2144. 7. Lambert, J., Chambers, A.R., Sinclair, I. and Spearing, S.M., 2012. 3D damage characterisation and the role of voids in the fatigue of wind turbine blade materials. Composites Science and Technology, 72(2), pp.337-343. 8. Nikishkov, Y., Bostaph, E. and Makeev, A., 2015, May. Nondestructive inspection of composite structures based on limited angle X-ray computed tomography. In American Helicopter Society 71th Annual Forum (pp. 5-7). 9. Tuegel, E.J., Ingraffea, A.R., Eason, T.G. and Spottswood, S.M., 2011. Reengineering aircraft structural life prediction using a digital twin. International Journal of Aerospace Engineering, 2011. 10. Glaessgen, E. and Stargel, D., 2012, April. The digital twin paradigm for future NASA and US Air Force vehicles. In 53rd AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics 25

and Materials Conference 20th AIAA/ASME/AHS Adaptive Structures Conference 14th AIAA (p. 1818). 11. Hernández, S., Sket, F., González, C. and LLorca, J., 2013. Optimization of curing cycle in carbon fiber-reinforced laminates: void distribution and mechanical properties. Composites Science and Technology, 85, pp.73-82. 12. Nikishkov, Y., Airoldi, L. and Makeev, A., 2013. Measurement of voids in composites by Xray Computed Tomography. Composites Science and Technology, 89, pp.89-97. 13. Seon, G., Makeev, A., Nikishkov, Y., and Fergusson, L. “Defect Formation in Contoured Composite Laminates During Vacuum Consolidation”, accepted for publication in the Journal of the American Helicopter Society, 2019 14. Kardos, J.L., Duduković, M.P. and Dave, R., 1986. Void growth and resin transport during processing of thermosetting—Matrix composites. In Epoxy resins and composites IV (pp. 101123). Springer, Berlin, Heidelberg. 15. Boey, F.Y.C. and Lye, S.W., 1992. Void reduction in autoclave processing of thermoset composites: Part 1: High pressure effects on void reduction. Composites, 23(4), pp.261-265. 16. Grunenfelder, L.K. and Nutt, S.R., 2010. Void formation in composite prepregs–effect of dissolved moisture. Composites Science and Technology, 70(16), pp.2304-2309. 17. Ledru, Y., Bernhart, G., Piquet, R., Schmidt, F. and Michel, L., 2010. Coupled viscomechanical and diffusion void growth modelling during composite curing. Composites Science and Technology, 70(15), pp.2139-2145. 18. Levy, A., Kratz, J. and Hubert, P., 2015. Air evacuation during vacuum bag only prepreg processing of honeycomb sandwich structures: In-plane air extraction prior to cure. Composites Part A: Applied Science and Manufacturing, 68, pp.365-376. 19. Gangloff Jr, J.J., Simacek, P., Sinha, S. and Advani, S.G., 2014. A process model for the compaction and saturation of partially impregnated thermoset prepreg tapes. Composites Part A: Applied Science and Manufacturing, 64, pp.234-244. 20. https://store.acpsales.com/products/3491/vacuum-bagging-pleats 21. Hou, T.H. and Jensen, B.J., 2004. Evaluation of Double-Vacuum-Bag Process for Composite Fabrication. NASA Technical Report. 22. Centea, T. and Hubert, P., 2011. Measuring the impregnation of an out-of-autoclave prepreg by micro-CT. Composites Science and Technology, 71(5), pp.593-599. 23. Torres, J.J., Simmons, M., Sket, F. and González, C., 2019. An analysis of void formation mechanisms in out-of-autoclave prepregs by means of X-ray computed tomography. Composites Part A: Applied Science and Manufacturing, 117, pp.230-242. 24. Larson, N.M. and Zok, F.W., 2018. Insights from in-situ X-ray computed tomography during axial impregnation of unidirectional fiber beds. Composites Part A: Applied Science and Manufacturing, 107, pp.124-134. 25. Carrier, B. and Granet, S., 2012. Numerical modeling of hydraulic fracture problem in permeable medium using cohesive zone model. Engineering fracture mechanics, 79, pp.312328. 26

26. Yao, Y., Liu, L. and Keer, L.M., 2015. Pore pressure cohesive zone modeling of hydraulic fracture in quasi-brittle rocks. Mechanics of Materials, 83, pp.17-29. 27. Nguyen, V.P., Lian, H., Rabczuk, T. and Bordas, S., 2017. Modelling hydraulic fractures in porous media using flow cohesive interface elements. Engineering geology, 225, pp.68-82. 28. Spencer, A.J.M. ed., 1984. Continuum theory of the mechanics of fibre-reinforced composites (Vol. 282). New York:: Springer. 29. Belnoue, J.H., Nixon-Pearson, O.J., Ivanov, D. and Hallett, S.R., 2016. A novel hyperviscoelastic model for consolidation of toughened prepregs under processing conditions. Mechanics of Materials, 97, pp.118-134. 30. ABAQUS, v2016. User's Manual, ABAQUS Inc, Pawtucket, RI, USA, 2016 31. Turon, A., Camanho, P.P., Costa, J. and Dávila, C.G., 2006. A damage model for the simulation of delamination in advanced composites under variable-mode loading. Mechanics of materials, 38(11), pp.1072-1089. 32. Seon, G., Nikishkov, Y., Makeev, A. and Shonkwiler, B., 2016. Mesh morphing methodology for strength predictions in composites. Composite Structures, 140, pp.612-620.

27

Highlights     

Introduction of a Digital Twin concept for mitigating void formation in autoclave composites Algorithm for automatic quantification of bulk content in uncured carbon/epoxy curved beam specimens using in situ X-ray CT data Pore-pressure cohesive zone modeling for FE representation of entrapped air pockets at ply interface in resin-saturated uncured prepregs Automatic transfer of CT information into FE models for creation of a Digital Twin reflecting bulk content in the physical part 2D FE simulation of debulking in unidirectional curved beam specimens and verification of the concept

28