Experimental and numerical evaluation of conduction welded thermoplastic composite

The capability of joining two thermoplastic composite parts by welding is a key technology to reduce the weight and cost of assembled parts and enables high volume manufacturing of future aeronautical structures made of thermoplastic composite materials. However, there is not much experimental understanding of the mechanisms involving welded joint failure


Introduction
The use of thermoplastic composite materials is gaining momentum in the transportation sector due to their improved mechanical properties, 'unlimited' shelf life and offers a number of advantages that can benefit cost-efficient and high-volume manufacturing.One of the main manufacturing techniques that enables this is thermoplastic welding.Joining of two parts is achieved by locally melting the material through application of heat and pressure.As a consequence of this local heating, very short processing times can be achieved (seconds to minutes) compared to typical manufacturing processes such as autoclave-, ovenand press-consolidation.The most established joining techniques for continuous fiber reinforced thermoplastic composites are resistance, induction and ultrasonic welding [1].The welding technique investigated in this work is conduction welding (Patent [2]), which is currently under development at Fokker/GKN Aerospace [3,4].Conduction welding is based on heating the surface of the part through an induction heated tool and generating the weld bath by heat conduction through the laminate.The benefit of this technique is that it does not require addition of welding specific materials such as energy directors or conductive strips.Furthermore, the technique is more suitable and scalable for welding of large parts and allows for absorbing of manufacturing tolerances like gaps.However, introducing these fastener-free joints also comes with new challenges.The strength of these highly loaded joints relies on the performance of the thermoplastic matrix which may be influenced by the welding process.It also needs to be investigated if the joint strength is influenced by the intra-and interply failure mechanism of the composite material.This may make it difficult to predict the strength of the welded joints without the availability of advanced predictive tools.
A popular approach to model failure of joints and interfaces is the Cohesive Zone Model (CZM).The main advantage of this method is the ability to simulate both the onset and propagation of damage without the need for an initial flaw.This methodology may provide to be a convenient method to evaluate the strength of adhesively bonded joints [5][6][7] and delaminations [8], but due to the complex failure behavior of welded thermoplastic composites [9,10] both inter-and intralaminar failure behavior has to be taken into account, as also observed for composite bonded joints [11].https://doi.org/10.1016/j.compstruct.2021.114964Received 7 July 2021; Received in revised form 7 October 2021; Accepted 3 November 2021 Modeling of fracture and damage in composites has been achieved through various approaches that model the crack in either a discrete or smeared manner.For the discrete approach, the eXtended Finite Element Method (X-FEM) [12] and similar enriched finite elementbased approaches such as the discrete cohesive crack approach [13] and the floating-node method [14] have been employed to simulate both intra-and interlaminar damage with high accuracy but at a high cost of computational efficiency.A methodology that provides more computational efficiency is Continuum Damage Mechanics (CDM).In continuum approaches, the effect of material degradation due to damage is modeled in a smeared manner through reducing the apparent stiffness of the material (strain-softening) inside the fracture process zone.Although strain-softening is not a real physical phenomenon at the micro-scale of the material, it provides an appealing framework to simulate the effect of damage at the meso-scale [15].Continuum damage models are already used to study the failure behavior of composite plain and open-hole strength coupons [15] and bolted single lap shear joints [16], however their application to welded composite joints is rather limited in literature.
The aim of this paper is to evaluate the strength and failure behavior of conduction welded joints and to develop a validated numerical methodology to support the design of thermoplastic welded joints [17].Single lap shear specimens are designed and manufactured by means of a conduction welding robot and are tested at the Fokker/GKN Aerospace R&D facilities.For the development of the numerical methodology both a simplified and a high-fidelity approach are proposed.The validity of the modeling strategy is assessed by comparing the predictions with the experimental results.The numerical approach is then used to study the most significant material parameters and model assumptions that may influence the strength and failure modes of these welded joints.

Manufacturing and testing of welded specimens
The section describes the design, manufacturing and testing of the conduction welded thermoplastic composite single lap shear joints.The properties of the thermoplastic composite material are reported first followed by the design of the welding rig and single lap shear specimen design.The welded joints are then manufactured using a conduction welding robot and the weld quality is inspected by means of non destructive inspection.Finally, static tensile tests are performed until failure and the results are reported and discussed.

Thermoplastic composite material
The material used in this research is the Solvay (formerly Cytec) thermoplastic polymer prepreg which consists of a fast crystalizing thermoplastic matrix of poly(ether ketone ketone) commonly referred to as PEKK-FC, reinforced with the continuous unidirectional AS4D fiber with a nominal ply thickness of 0.138 mm.The coupons for material properties characterization are manufactured by means of autoclave consolidation and tested by Fokker/GKN Aerospace at room temperature ambient conditions and are reported in Table 1.The Young's modulus is given in both tensile and compressive direction.In the analysis, the longitudinal tensile modulus (loading direction) is used, while the average modulus is used for the transverse direction.For the matrix fracture toughness both the initiation and propagation values are measured.The longitudinal fracture data are not measured but are assumed to be the same as for the AS4 fiber reported in [15].The longitudinal strength and fracture toughness ratios are assumed to be similar to [18] and are validated with internal Fokker/GKN Aerospace experimental open-hole coupon data.
From literature it is known that the manufacturing process may significantly influence the material properties of thermoplastic composites.Parameters such as the local temperature profile and gradients, heat-up speed, constant time and temperature at the welded interface, may affect the melting conditions [19] or flow [20] of the polymer and the cooldown speed may significantly influence the final crystallinity of the material.Furthermore, thermal residual stress may be introduced during the welding process due to both the high temperature and shrinkage of the material during cystallization.Sacchetti et al. [21] showed an increase of factor 2.5 on the fracture toughness by changing the cooling rates during processing.The main reason can be explained by changes in crystallinity of the semi-crystalline polymer that can affect the ductility and fracture properties of the thermoplastic composite material.However, it is shown in [19] that for the material used in this work (PEKK-FC), even at high cooldown rates (−60 • C/min) the material still achieves high levels of crystallinity (Xc = 27%).This is however, still significantly lower compared to crystallinity (Xc = 36%) measured with very slow cooldown rates (−0.5 • C/min), so some influence on the local material properties due to changes in crystallinity is expected.Grouve et al. [22] investigated the sensitivity and effect of material variability of the C/PEKK-FC material during induction heating.It was found that laminates showed a inhomogeneous fiber distribution and resin rich regions were identified.This variability may have a significant effect on temperature evolution in the welded joint.Furthermore, it is also known [23] that excessive resin or a thick interface may also increase the fracture toughness, as this may allow for more plastic deformation at the crack tip.
Accounting for these influences in the analysis would require detailed knowledge on the local thermal behavior and characterization of the in-situ material properties, which is out-of-scope for this work.However, the developed numerical methodology is used to study the possible influence of changes in material properties on the joint strength and failure behavior.This is discussed in Section 4.4.

Design of a welding setup and welded single lap shear specimens
The rig shown in Fig. 1a is designed to support the welding of two laminates by means of a conduction welding robot.The laminates are welded at three positions and specimens are machined from the center section of each weld as shown in Fig. 1b in order to create a single lap shear specimen.The specimen design (Fig. 1c) follows ASTM standard D3165 [24], with overlap length of 30 mm and nominal width of 25.4 mm.The length of each composite laminate is 75 mm and the height of the weld at the interface is approximately 18 mm.This dimension is a typical value measured from the weld size of the tested samples.The layup of the composite laminates consists of 16 plies ([−45,0,45,90,−45,90,45,0]  ) for each laminate which is rotated to create a +45∕ − 45 interface at the weld.The 0-degree ply is oriented in the loading direction of the specimen.

Conduction welding of the laminates
All laminates used in the test campaign are manufactured from two large autoclave consolidated AS4D/PEKK-FC laminates according to a typical process cycle for aeronautical parts and are cut into smaller laminates for welding.Additional PEKK-FC foil is included in each laminate to make sure sufficient resin is present at the weld interface.The welding robot, shown in Fig. 2, consists of a clamping system which applies local heat and pressure through an induction heated stamp.The two AS4D/PEKK-FC laminates are positioned and clamped by the robot after which heat and pressure are applied to locally melt the thermoplastic material.The system is heated from one side at approximately 400 • C in order to reach the typical processing temperature of 377 • C [25] at the welded interface and to ensure a complete melt [19].The maximum temperature is limited to approximately 400 • C in order to prevent polymer degradation at the laminate surface near the welding tool.This value is based on Fokker/GKN experience and can be determined by techniques such as thermogravimetric analysis (TGA).The welding temperature is measured at the weld interface by means of several thermocouples and the quality of the weld is inspected by means of Non Destructive Inspection (NDI) by using a phased array setup that inspects the width and the position of the conduction weld.The equipment consist of the Omniscan MX and a probe that generates ultrasonic soundwaves.When putting the probe on the composite laminate the sound waves travel through the composite plate and reflect back to the probe.The sound waves are dampened when a change in density is detected.This can be an indication of voids, delaminations or an area which is not welded.The different type of scans that are performed are shown in Fig. 3 on a conduction welded laminate.NDI Cscan results are shown in Fig. 4 where the colors represent the damping of the sound.The color red is 0% and blue is 100% attenuation of the ultrasonic sound waves.The phased array probe can perform different types of scans.The A-scan provides information about a local point of the weld; this point is perpendicular to the welding location.The C-scan provides information about the weld in top view.The S-scan provides information about the cross-section of the weld, perpendicular to the weld, this is a combination of the data from the A-and C-scan.
The C-scan results of the first welded laminate (Specimens 1-3, Table 2) is compared in Fig. 4a to a welded laminate with a more consistent weld quality and higher strength (Specimens 7-9, Table 2) in Fig. 4b.Fig. 4a shows some indication of delamination or variation through the thickness (blue-green region with <50% attenuation) compared to the higher quality welds in Fig. 4b that show a more consistent C-scan result above 50% attenuation.

Experimental results
In total, 21 conduction welded single lap shear joints are manufactured and tested.Static tensile tests are performed according to ASTM D3165 [24] until failure.
The tests are carried out using a ZWICK 100 kN universal test machine with self-aligning wedge grips.The cross-head speed is set to 1 mm/min and the test is automatically stopped after measuring a load drop to 75% of the highest measured force.It is noticed that the selfaligning grips do not provide clamped conditions at the grip area and some gapping is observed due to secondary bending, so the effective specimen length between the grips is unknown.Furthermore, the clamping pressure of these wedge grips increases during loading due to axial movement.This makes comparing load-displacement curves of the experimental and numerical results difficult, as the machine displacement includes both the machine compliance and the movement of the wedge grips.Therefore, the numerical load-displacement curve is only compared to the experimental failure loads.
The results are shown in Table 2.The specimens are grouped according to the laminates from which they are machined.The relation of the weld position with the specimen, is also given.Furthermore, the height of the weld and the measured weld area based on the fracture surface are also reported.The average failure load is 18841  and the minimum and maximum failure loads are 12307  and 21556 , respectively.
The failure modes of a typical lower (Specimen 3) and higher bound experimental result (Specimen 8) are shown in Fig. 5.The fracture surface of the weld is shown for both the top and bottom laminate as indicated in Fig. 5a.For Specimen 3 the failure mode extended outside of the weld (Fig. 5b), which may be related with the NDI indication at the same location as shown in Fig. 4a.It is also observed that the failure mode of Specimen 8 show clear signs of a resin rich area at the tip of the weld and clear signs of delamination propagation (light colored zone) as indicated in Fig. 5c.The failure modes are investigated and explained in more detail in Section 4 evaluating the results from the numerical analysis.

Modeling strategy
Two numerical strategies are explored to analyze the strength of the thermoplastic conduction welded single lap shear joints.The first simplified approach considers failure only at the welded joint and simplifies all damage mechanisms into a single cohesive interface (Fig. 6), while the second follows a high-fidelity approach that considers both inter-and intralaminar damage in each ply of the single lap shear joint and the welded interface (Fig. 7).This section first explains the general model parameters and geometry, followed by the damage models used in the work and a study on how to efficiently simulate single lap shear joints in ABAQUS/Explicit [26].

Models and parameters
The two models follow the specimen design and geometry as shown in Fig. 1.The single lap shear joint is divided into two zones and assembled from separate parts.The composite sections at the load introduction, which are far from the welded joint, are discretized with through-thickness continuum shell elements (SC8R) as shown in Fig. 6a and considers only linear-elastic material behavior.The central zone of the joint in the simplified model is also modeled using SC8R elements but a smaller mesh size is used to meet mesh size requirements [27] for cohesive zones.The central section consists of a surface with contact (Fig. 6b) and cohesive surface for the welded interface (Fig. 6c).
The central zone of the high-fidelity model, also referred to the damage zone, is shown in Fig. 7a and follows a ply-by-ply modeling strategy where each ply is discretized as a layer of reduced-integration solid elements (C3D8R).A fiber aligned mesh with an aspect ratio of three is used following the guidelines given in [15].The welded interface, as in the simplified model, consists of a cohesive surface for    the weld and a contact definition outside of the weld, as shown in Fig. 7b.
The connection between the damage zone and supports is implemented with kinematic constraints (*TIE) in order to transfer displacements and rotation across their boundaries.The boundary conditions are applied by means of a velocity amplitude profile imposed at the top and bottom surfaces of the supports.The fixed boundary condition at the bottom surface is enforced by means of conditions of zero velocity while the velocity at the top surface is ramped-up to a constant velocity until failure.The simplified model consists of approximately 10.000 elements, while the high-fidelity model consists of approximately 100.000 elements.The optimization of the analysis parameters to achieve efficient analysis times is explained in more detail in Section 3.3.

Modeling damage behavior of welded single lap shear specimens
Failure of the welded joints is considered to be similar to the failure of general composite laminates due to the nature of the joining process.The laminate is locally melted and consolidates in a similar manner as consolidating plates in an autoclave process, if the right manufacturing conditions are respected.This means that the welded interface is locally indistinguishable from the bulk material.
Therefore, the failure of the welded joint is expected to behave the same as interlaminar damage in the thermoplastic composite material and that a zero-thickness interface can be assumed, which is different from bonded joints that have a non-zero adhesive thickness.
The failure mechanisms of the welded joint can therefore be divided into the typical categories of composite failure modes, namely interlaminar and intralaminar damage.Failure of the weld and failure of the ply-to-ply interface of the laminates fall into the interlaminar category.This damage behavior affects the separation between plies, which forms a delamination that can occur under different opening modes.Intralaminar damage considers all the failure modes that occur within each ply such as fiber and matrix failure.The interlaminar model already available in ABAQUS/Explicit has been used, while the intralaminar model has been implemented in a Continuum Damage Model (CDM) through a user-defined ''VUMAT'' subroutine using a numerically explicit integration scheme.

Interlaminar damage
The interlaminar behavior of the welded joints is modeled through the general contact algorithm available in ABAQUS/Explicit [26].This method takes care of the kinematics of surface contact, cohesive and frictional behavior.The cohesive zone model describes the opening of the delamination in terms of tractions and displacements.Once damage is initiated, the model reduces the stiffness of the cohesive surface thus decreasing the traction while dissipating the fracture energy corresponding to the specific mixed-mode opening mode, as given by the Benzeggagh-Kenane criterion [28].Damage initiation is identified by means of a quadratic nominal stress failure criteria which is a function of the interlaminar strength values for each damage mode.Frictional effects are considered on the surface to include possible effects of ply friction within the delaminations.For the simplified model, this behavior is only considered at the weld interface as shown in Fig. 6c, while for the high-fidelity model it is considered for each ply interface.Outside the weld surface, only contact and frictional behavior are considered.A value of 200 000 N/mm 3 is used for the mode I penalty stiffness and Turon's equation [29] is used to calculate the shear penalty stiffness value.Friction is considered to be ply interface angle dependent and follows the approach and values from [30].

Intralaminar damage
CDM is used for modeling damage in the plies taking into account the three-dimensional stress states through the physically-based threedimensional failure criteria proposed by Catalanotti et al. [31].This effect is considered to be important for the analysis of single lap shear joints because significant out-of-plane loading may be present due to secondary bending.The implementation follows the approach defined in [30] and is originally based on the work of Maimi et al. [32] that guarantees the correct energy dissipation for each composite fracture mode.The model is also enhanced to account for large shear deformations through the use of the 2nd Piola-Kirchhoff stress and Green-Lagrange strain to accurately predict matrix cracks and ply splits as demonstrated in [33,34].Capturing large shear deformations is relevant for accurately predicting the matrix-dominated failure modes of the welded joints.The implementation following [35] is briefly explained below and the influence of some of the assumptions are discussed in Section 4.
The stress and strain in the element is based on a Lagrangian kinematic measure where the constitutive equations can be defined within a orthonormal material frame.Therefore, fiber rotation due to large shear deformation is intrinsically captured.The 2nd Piola-Kirchhoff stress  is the work conjugate to the Green-Lagrange strain , which is determined from the deformation gradient tensor : where  is the identity tensor.The 2nd Piola-Kirchhoff stress  can be determined from the material stiffness tensor  and the Green-Lagrange strain : and can be mapped to the current configuration by where  and  are the Kirchhoff and Cauchy stress.The stress has to be returned in the co-rotational basis of Abaqus (  ) through a change of basis operation by using the rotation matrix .
The strain tensor for an orthotropic ply is defined as {} = [(  )] {} + {} , where {} is the tensor of effective stresses in Voigt notation,  is the isotropic thermal load and {} are the ply thermal expansion coefficients following the normal orthotropic directions.Although the model can take into account the influences of thermal stress due to the welding process, this has not been considered in this work, as this would require detailed knowledge and prediction of the local thermal gradients and material behavior during welding conditions.
The compliance tensor, following the convention adopted in the ABAQUS/Explicit implementation, is expressed as: Softening laws are used to ensure physically correct dissipation of fracture energy for each failure mode.The exponential softening law in general form, which is used for all matrix failure modes, is expressed in Eq. (6).
where   ( = 1±, 2±, 3±, 4, 5, 6) are elastic domain thresholds, which are initially 1.0 for an undamaged material and increase after damage initiation.The physically-based three-dimensional failure criteria proposed by Catalanotti et al. [31] is used to identify initiation of damage and can also account for combined loading effects.The parameter   ensures that the dissipated energy is independent of mesh refinement [36] by relating the element characteristic length  *  ( = 1, 2, 3) for a fiber aligned mesh with the material and fracture properties in the corresponding directions.For matrix damage the equation is as follows [30]: Matrix cracking is assumed to occur under general mixed-mode (mm) conditions, with initiation   predicted by the three-dimensional failure criteria and the propagation fracture energy   by the energy-based Benzeggagh-Kenane (BK) criterion [28].
Softening in fiber direction requires taking into account both fiber pull-out and bridging and is implemented through a trilinear softening law following the implementation of CompDam [18].However, this failure mode is less relevant for the matrix dominated failure modes of the welded single lap shear joints.
Failure in out-of-plane direction is accounted for by the cohesive surfaces between the plies and can be disabled in the CDM, but nonlinear elastic-plastic response is taken into account for each shear direction (1-2, 1-3, 2-3).This is achieved using the Ramberg-Osgood equation following the implementation of CompDam [18].The engineering shear strain (  ) is defined as: where   and   are the fitting parameters obtained with experimental data from in-plane shear tests.The nonlinear shear response is assumed to be plastic until the onset of damage, after which the element is damaged following an exponential softening law that is regularized according to the mode II fracture toughness.The use of exponential softening laws allows for gradual changes in stiffness and a more robust solution, however during softening the strains in the elements may reach very high strain levels.Furthermore, highly distorted elements, which can be detected by sudden changes in volume, may influence the stability of the analysis.These elements can be detected through the determinant of the deformation gradient det() [37].The criteria for element deletion is adopted from [15], with more strict conditions during very large deformation.In summary, the

Analysis parameters for the explicit finite element method
The explicit finite element method implemented in ABAQUS is selected for the analysis, as it involves solving a highly nonlinear dynamic problem such as large displacements, non-linear material behavior including damage and complex contact interaction with damage and frictional behavior.In order to perform computational efficient analyses in ABAQUS/Explicit, the total analysis time should be as short as possible.This can be achieved by mass-scaling or through the loading velocity.Mass-scaling up to 1000x as proposed in [15] combined with high loading velocities may not be a issue for in-plane coupons such as plain and open-hole strength, but they may cause unwanted dynamic effects such as oscillations and overshoots for problems that are more sensitive to dynamic or inertia effects.This may be the case for single lap shear joints as in-plane loading results in out-of-plane deformation due to secondary bending.
The three main parameters that determine the total analysis time are (1) the stable time increment related to element size and density, (2) the ramp up time or loading amplitude and (3) the loading velocity.Although the element size is a parameter that could be adjusted, it is indirectly limited through the mesh size requirement of the damage models.This leaves the density, loading amplitude and velocity as the main parameters for the investigation.
The simplified model has been used to conduct the parameter study.The following combination of parameters is investigated: Where the baseline analysis is defined with a density scaling factor of 1, a loading amplitude of 0.01s and a velocity of 25 mm/s.The total analysis time for this case is over 2 h on 10 CPUS (Intel(R) Xeon(R) CPU E5-2640v4 @ 2.40 GHz).The predicted failure load is checked against the solution of the implicit solver and the results are nearly identical.
If is found that by using the maximum mass-scaling factor and loading velocity, severe dynamic effects occur and the failure load is over-predicted by approximately 10 percent.This means that the typical parameters for in-plane coupons [15] are not valid for single lap shear joints.
Based on the results of this study, it is chosen to select a loading velocity of 50 mm/s in combination with a ramp up time of 0.0025 s and 10 times mass-scaling.This results in an analysis time of approximately 20 min on 10 CPU which is a significant reduction in analysis time compared to the baseline of over two hours, while staying within 1% of the expected failure load.On the high-fidelity model this results in a analysis time of 8 to 10 h when using 20 CPU.

Numerical evaluation of thermoplastic conduction welded single lap shear joints
The strength and failure modes of the welded joint are evaluated using the two modeling approaches and comparisons are made against the experimental results.First, only failure at the welded interface is considered using the simplified approach followed by a discussion on the apparent fracture toughness of the welded joint.The high-fidelity model is then used to investigate failure of the surrounding plies and to study the behavior due to these additional failure modes.This is followed by a discussion on the influence of the inter-and intralaminar material properties and some limitations of the modeling approaches are identified.

Numerical evaluation using the simplified modeling approach
The single-lap shear specimens are at first simulated using the simplified approach which only accounts for damage in the welded interface.The material properties used are reported in Table 1 and the ply thickness is adjusted to 0.144 mm, based on the measured thickness of the test specimens.The mode I and mode II propagation fracture toughness values measured from unidirectional specimens are used in the cohesive model for the weld.The grip length is conservatively taken as the full specimen length as support conditions of the self-aligning wedge grips are difficult to determine.Some analyses are performed to study the influence of the gripping distance and the predicted failure loads would increase by 1%-2% due to slightly reduced secondary bending of the specimen.The numerical load-displacement curve is compared against the upper and lower bound experimental results in Fig. 8a.
The status of the weld can be represented by the cohesive damage variable (CSDMG) that identifies the area of the weld undergoing degradation.This is shown for two locations along the load-displacement curve in Fig. 8b.During damage initiation at the edge of the welded joint a Fracture Process Zone (FPZ) develops where the damaged area becomes critical, causing the load-displacement curve to suddenly drop.The predicted failure load using this approach is 11 597  which is just below the lower experimental bound equal to 12 307 .
The mismatch between the numerical and experimental results can be caused by several factors.The most likely causes are the influence of the +45∕ − 45 interface angle, the more complex failure modes, that are not taken into account by the simplified modeling approach, and the influence of the manufacturing process on the material properties as discussed in Section 2.

Apparent fracture toughness of conduction welded joints
The simplified modeling approach is used to investigate how the fracture toughness determines the upper and lower bounds of the experimental campaign.The strength of the joint is analyzed scaling the interlaminar fracture toughness of the material, which can be referred to as the apparent fracture toughness of the joint.
The result of this study is shown in Fig. 9 where the joint strength is plotted against a factored interlaminar fracture toughness.It becomes clear that the interlaminar fracture toughness has to be increased by approximately a factor three in order to achieve the upper bound experimental results.The average experimental results can be approached by using a factor of 2.5.
The apparent fracture toughness of the joint can be used in an engineering approach to take into account effects such as changes in material properties, dependence of the fracture properties on the interface mismatch angle or complex failure mechanisms that cannot be predicted with the simplified modeling approach.The benefit of this is that the effective material properties can be determined at the coupon level and applied in coarser models at the structural level to efficiently perform simulation at that scale.Compared to thermoset composites, the fracture toughness of thermoplastic composites is much higher, making it easier to meet mesh-size requirements in cohesize zone or VCCT analyses [27,38] with a coarse mesh and reduce mesh dependency.However, the validity has to be further investigated for different materials, loading modes and interface angles.

Numerical evaluation using the high-fidelity modeling approach
The high-fidelity modeling approach described in Section 3 is used to predict both the failure load and failure mode of the welded single lap shear joint.The same geometry and unidirectional material properties are used as the simplified model, but now following the plyby-ply modeling strategy as shown in Fig. 7 that allows for both intraand interlaminar damage.The numerical load-displacement curve is compared with the upper and lower bound experimental results, and also with the curve obtained using the simplified model, in Fig. 10a.
It is evident that the joint strength is highly influenced by the complex interaction of failure modes in the surrounding plies of the laminates.The predicted failure load following the high-fidelity approach is 14 962 N. Fig. 10b shows the failure modes at two points along the load-displacement curve prior to final failure.
Both the weld failure (delamination, CSDMG) and matrix damage (D2) are shown in Fig. 10b for each side of the joint indicated with the top and bottom view.The failure process predicted by the numerical model at the first point is as follows: (i) Matrix damage starts to develop in both 45 degree plies near the edges of the weld.This occurs at approximately the failure load of the simplified model.(ii) The damage appears to slow down the initiation of delamination in the weld interface, while the matrix damage starts to develop ply splits along the fiber direction in the 45 degree plies.Further along the loaddisplacement curve at the second point these ply splits (iii) continue to grow and migrate into delaminations between the first and second ply above the weld interface and (iv) will cause the joint to fail after they become critical.
The final failure mode of the numerical model is compared in Fig. 11a to the experimental failure mode in Fig. 11b and similarities can be identified.The main failure mode, that most likely caused the failure of the single lap shear joint, is the large delamination between the first ply at 45 degree and second ply at 0 degree (I) and not at the +45∕ − 45 degree welded joint interface, because a large portion of the zero degree ply (II) is visible.
Another interesting feature is the location at where the 45 degree ply failed by means of fiber failure (III) which is at an off-set from where the delamination initiation is predicted by the numerical model.The white markings on the test sample are a sign of resin rich areas at the edge of the weld which may increase the initiation strength and fracture toughness locally [23].Several ply splits of the +45∕ − 45 ply (IV) are visible on both sides of the specimen, which led to the pull-off of one of the 45 degree plies (V) due to the interaction with delamination on the 0/45 interface.

Influence of the fracture properties
The influence of the intra-and interlaminar and weld interface fracture properties on the single lap shear joint strength and failure mode is investigated by using the high-fidelity modeling approach.Changes in these properties could, as already discussed, be caused by the manufacturing process but also by assumptions in the modeling approach such as taking into account in-situ material strength or by accounting for the full cohesive law of the different fracture modes and/or the effect of interface layup angles at the welded interface versus +45∕ − 45).
Several different models, see Fig. 12, are run with different fracture toughness values for intra-and interlaminar fracture properties in order to gain insight into how the material properties influence the joint strength and corresponding failure modes.
The CS and HF models are the models considered in the previous sections and are the simplified continuum shell and high-fidelity model, respectively.The HF: weld model is the high-fidelity model with an increased fracture toughness at the welded interface.A factor of 2.5 is chosen based on the prediction of the average experimental results by using the simplified approach.This brings the joint strength to 15776 N, which is much less compared to scaling the weld fracture toughness in the simplified approach.This difference between the simplified and high-fidelity approach is explained by the change in failure mode to matrix cracking and failure of the first ply.Interestingly, this also suggests that having increased material properties at only the welded interface cannot explain the upper bound experimental results.
If only the matrix fracture toughness of the ply (HF: intra) is increased the failure mode becomes delamination dominated, as shown in Fig. 12b.Delaminations initiate (I) at the same level as the HF model and migrate through matrix damage (II) to the first ply and eventually pull-off (III) the first ply.Due to the increased matrix fracture toughness the load is increased from 14 962 N to 17 509 N.
When both, the interlaminar fracture toughness of the first ply and the weld (HF weld+inter), are modified, the strength further increases to 18 685 N and a more ply dominated failure mode, as shown in Fig. 12c, is observed.In this situation, there are no signs of initial delaminations (IV) and failure of the joint is caused by matrix failure (V) of the first ply.The strength can be even further increased if the fracture toughness of all adjacent plies and interfaces is increased (HF: all).
The under-prediction of the numerical simulations cannot be attributed only to the properties of the weld and of the plies, and require further analysis on the damage development and assumption of the numerical model.When taking a closer look at matrix failure of the ply, it is observed that the loading of the matrix cracks changes as damage progresses.Matrix damage and welded interface initiation is   mostly due to in-plane shear loading with some influence of mode I opening due to secondary bending, but as damage progresses the cracks become predominately loaded in out-of-plane shear (Fig. 13a) and element deletion is triggered when 100 percent out-of-plane shear is reached.This identifies one of the limitation of the current modeling approach.Changes in loading mode of existing in-plane cracks, changes in crack angles of partially developed cracks and possible effects of friction due to crack closure are not taken into account.The influence of these effects may effectively delay the propagation of damage and significantly increase the strength of the matrix dominated failure modes.
Although the modeling strategy is already able to limit the delamination growth due to damage and softening of the plies, it is also observed that the predicted crack angles for matrix damage are generally between 20 and 45 degrees as shown in Fig. 13b, where a 0 degree crack represents a pure in-plane matrix crack [31] or no damage initiation.At some locations near the start of the fully developed ply splits the angles are closer to in-plane failure while at most other locations they are at the lower angles, which supports the experimentally observed migration of delamination to the second interface.As the current modeling approach is using the Abaqus/Explicit built-in failure criteria, there is no control in preventing failure of the welded interface when crack migration is detected.
Another important consideration is that the bending stiffness of the laminate will influence secondary bending of the welded joint and will change the local mode mixity at the welded interface.This means that a more stiff or hard laminate will experience less secondary bending and a mode II induced failure mode.Secondary bending is higher for soft laminates and will result in a higher mode I component at the welded interface.Furthermore, since the joint strength appears to be highly influenced by the matrix dominated failure modes of not only the weld but also the composite laminate, it is evident that the interface angle and ply orientation near the weld plays an important role for the strength of the joint.These effects may have to be considered when establishing design guidelines for composite welded joints.

Conclusions
Single lap shear thermoplastic composites welded joints have been designed, manufactured and tested.Evaluation of the experimental results provided new insights into the complex failure behavior of the joint and of the interaction with the failure modes of the laminate.For the numerical analysis both a simplified and a high-fidelity approach are developed and evaluated.It is found that the simplified approach, based on the cohesive zone method, predicts a very conservative joint strength when unidirectional interlaminar fracture toughness properties are used.The apparent fracture toughness of the joints is found to be approximately 2.5 times higher.A high-fidelity modeling approach is developed to investigate the influence of the complex failure behavior of the welded joint using an improved Continuum Damage Model to accurately predict matrix failure and ply splitting.This model provides new insights in the failure behavior of the joint and is able to accurately predict the failure mode.However, predicting the upper bound experimental results is still difficult and some limitations to the numerical methodology are identified.This includes changes in modemixity during crack propagation, non-zero matrix crack angles near the welded interface and not considering frictional effects on the fracture plane.The analysis using the high-fidelity model shows that different material properties for the inter-and intralaminar failure modes have a strong effect on the joint strength and may significantly influence the failure modes.A better understanding of the material properties of the welded joint is still needed and design guidelines may need to not only consider the welded interface but also the surrounding plies.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 7 .
Fig. 7. High-fidelity single lap shear model: (a) Overview of different model features; (b) Details on welded surface showing only the lower laminate.

Fig. 8 .
Fig. 8. Simplified model result: (a) Load-displacement curve; (b) Initiation of weld failure and failure of welded interface after reaching the peak load.

Fig. 10 .
Fig. 10.High-fidelity model results: (a) Load-displacement curve compared to simplified model; (b) Predicted weld and matrix damage at two points along the load-displacement curve.

Table 2
Welded single lap shear joint experimental results.