A review of experimental and theoretical fracture characterization of bi-material bonded joints

High performance structures require the use of different materials to meet their demanding requirements. Especially fibre reinforced polymer composites are nowadays often bonded to metals in order to take the most advantage of the materials properties and to minimize their disadvantages. However, the interface in such bi- material assemblies often represents the weakest point and thus has to be carefully addressed to ensure structural integrity. This review paper presents an overview of the research on bi-material interface crack problems over the past 30 years. Three categories of the research are discussed: mechanical testing, crack driving force and mode partitioning. The literature reveals that the key element to the fracture analysis of the bi-material interface crack is how to perform the mode partitioning. The proposed theories for mode partitioning by many researchers are meaningful yet underdeveloped and need further experimental validation.


Introduction
There is a renewed interest in the subject of bi-material interface fracture owing to the necessity of understanding the failure of adhesively bonded composites-to-metal assemblies which are increasingly employed in engineering structures. The fibre reinforced polymers (FRP) are attractive for high-performance structural applications. Adhesive bonding of FRP composites-to-metal gives the opportunity to provide hybrid structures an excellent performance in traditionally metal-based industries. Several examples can be found in different industries. In the aerospace industry, the fibre metal laminate design concept, where thin metal sheets are bonded together with FRP, has succeeded commercially with the application of Glare as fuselage skin material on the jumbo airliner Airbus A380 [1][2][3]; Titanium strips are bonded to the edges of the aero engine fan blades that are made of carbon fibre reinforced polymers (CFRP) in order to improve the integrity of the blades without a significant weight penalty [4]; Metal skins-to-CFRP stiffeners are being investigated in order to design leading edge structures with the hybrid laminar flow control technology for reducing the drag and improving the fuel efficiency of airplanes [5][6][7][8]. In the maritime industry, steel is being partly replaced by composites to construct hybrid ship hulls to achieve a better stealth performance and weight reduction leading to significant research on steel-to FRP bonded joints [9][10][11][12]. In the automotive industry, bonded hybrid FRP/steel structures have been explored to be used as crash absorbers [13].
A major concern about the adhesive bonding of composites-to-metal is the failure at the interface: a bi-material interface failure problem. Strictly speaking, a bi-material interface failure is not only the failure at the apparent composite/metal interface, but also the delamination between two plies of different fibre orientations in a composite laminate. The difference in the stiffness between the two delaminated plies due to fibre orientations results in a crack interface between two asymmetrical sections. Thus it should be regarded as a bi-material interface crack problem [14,15].
It is imperative to understand the fracture behaviour of such bimaterial interface failure to ensure structural integrity. The interface crack problem has attracted attention for several decades. Far earlier than any experimental investigations on bi-material interface failure, the theoretical analysis of interface fracture mechanics had started since 1959. Williams published the first paper on the interface crack between two dissimilar homogeneous and isotropic planes [16]. Oscillatory stress singularity in the vicinity of the crack tip and unrealistic crack surface interpenetration were found based on his theoretical work. Following his work, Erdogan [17], England [18] and Rice [19] analytically derived stress intensity factors (SIF) for the interface crack between two half planes in spite of the unrealistic crack surface overlap, while many others employed numerical approaches [20,21]. To overcome the unrealistic crack surface interpenetration, many theoretical investigations adopted the modification that a small portion of two crack surfaces behind the crack tip could be in contact [22][23][24][25][26][27][28][29]. Nevertheless, the prediction models with and without the modification provide stresses that agree well except in the small contact zone very near the crack tip. Other researchers also took the effect of crack tip plasticity into consideration [30][31][32][33]. Dynamic crack growth along bi-material interfaces were experimentally investigated to obtain the near tip fields and propose failure criteria by Rosakis and co-workers [34][35][36][37].
The knowledge developed by analysing the interface crack between two half planes has been extended to the interface problem in layered structures. Since late 1980s, papers on the interface fracture mechanics in layered bi-material systems started to be published, including theoretical analysis [38,39] and experimental work [40][41][42][43][44]. In order to avoid the complexity of calculating the stress singularity at the bi-material interface, strain energy release rate (SERR), G, is regularly chosen as the crack driving force [45]. SERR can be calculated by two means: Irwin's crack closure integral [46] and the compliance method [47]. J-integral has also been applied as the crack driving force to account for the failure behaviour of the ductile adhesive in bonded bi-material joints [48,49]. It is worth noting that most bi-material interface problems are inherently mixed mode problems as a result of loading asymmetry and material property asymmetry across the interface [50]. The challenge of using the SERR or J-integral as the crack driving force, however, is that they provide an overall quantity without detailed information regarding the percentage of mode I and mode II, i. e., mixed mode ratio. Mode partitioning has to be performed for this case since the fracture behaviour in bonded joints is mode mixity dependent [40,47,51].
In the context of describing mixed mode fracture problems, two terminologies are frequently used in literature: mixed mode ratio and mode mixity. It is deemed beneficial to provide somewhat clarification on these two terms. Mixed mode ratio is normally referring to the ratio of the mode II SERR G II , to the total SERR G total , i.e., G II / G total . The total SERR, G total , is the sum of the mode I SERR, G I , and mode II SERR, G II . Mode mixity or mode mixity phase angle is given by arctan( or arctan(K II /K I ) with K I and K II being the mode I and mode II SIFs.The objective of this review is to provide an overview of the available characterization methods of fracture behaviour in adhesively bonded bimaterial joints. The layered nature of most composites and the popularity of testing composites with beam-like specimens motivate this review paper to have a strong focus on the interface problem in beam-like joints with dissimilar materials bonded together. The key to the bimaterial interface crack problem is to calculate the crack driving force such as the SIF, SERR or J-integral for a given bonded assembly under a prescribed loading case. Over the last 3-4 decades, many bi-material joint designs together with calculation methods of the crack driving force are proposed. Based on the fracture modes, the characterization methods are categorized into three groups in this review paper: • Mode I bi-material interface problem. The purpose is to characterize the opening mode of an interface crack. • Mode II bi-material interface problem. The purpose is to characterize the in-plane sliding mode of an interface crack. • Mixed mode bi-material interface problem. The purpose is to characterize the state of an interface crack subjected to the opening mode and in-plane sliding mode.
As it has been explained earlier, most bi-material interface cracks are inherently mixed mode. Mode I and mode II bi-material interface cracks can only be achieved by deliberately designing a bonded joint under a specifically prescribed loading case. Regarding the mode partitioning methods, no consensus has been reached by researchers on which method is the most appropriate. The popular methods to partition the overall SERR into two components of pure modes are summarized and reviewed.

Mode I bi-material interface testing
Mode I loading is the severest loading case for delamination in composites and interface cracking in bonded assemblies. Consequently, the bi-material interface crack subjected to mode I loading has attracted extensive attention from researchers, leading to numerous publications that will be discussed hereafter.
The double cantilever beam (DCB) specimen configuration has been widely applied to obtain the fracture toughness or investigate the fatigue behaviour under mode I loading of FRP composites and adhesively bonded joints. Two standards [52,53] are available for the characterization of mode I fracture toughness of unidirectional fibre reinforced plastic composites. The DCB configuration has also been modified to serve other purposes, such as the Tapered DCB (TDCB) specimen, which is adopted to obtain the SERR as a function of the contour [54][55][56], and the Asymmetric DCB (ADCB) specimen, which is adopted to obtain mixed mode failure behaviour [57][58][59].
The DCB specimen configuration has been extended to test the interface crack between two adhesively bonded dissimilar materials. In Fig. 1 a general DCB configuration is illustrated. In open literature two types of bi-material DCB configurations have been proposed to characterize the bi-material interface crack subjected to mode I loading: 1) the first one with the two arms having the same flexural rigidity [60][61][62][63][64], and 2) the second one with the two arms having the same longitudinal strain distribution at the faying surfaces [14,65,66].
The mode I testing of the bi-material interface crack is initially based on the DCB specimen designed for cracks between two identical arms with the same flexural stiffness which bend symmetrically under a pair of opening loads. The same line of reasoning has been followed by many researchers for testing bi-material interfaces, resulting in the bi-material DCB specimen with two arms that have the same flexural stiffness, which can be mathematically expressed as [60][61][62][63][64]: where E and I denote the longitudinal Young's modulus and the moment of inertia respectively. The subscripts 1 and 2 denote the two different adherends respectively, in accordance with Fig. 1. A simple manipulation of Eq. (1) gives the following expression: where h 1 and h 2 denote the thickness of the upper and lower adherends respectively. Many researchers have claimed that the shear loading which induces a mode II fracture component can be eliminated by using the DCB configuration designed according to Eq. (2) [60][61][62][63][64]. Mode II fracture component was found in the numerical results, in contrast to what they believed [61,62]. It was therefore proven that this phenomenological design method of matching the bending stiffness of both arms for mode I fracture in bonded bi-material joints fails to guarantee a pure mode I loading condition.
In the second type of bi-material DCB configuration, the joint design is such that the two arms have the same longitudinal strain distribution at the faying surfaces, i.e., ε 1 = ε 2 , under the loading condition shown in Fig. 1. This is proposed to remove the in-plane sliding associated with the mode II fracture component [65]. The following equation describes the joint configuration that meets this design criterion.
Ouyang et al. [14,66] reach the same conclusion that Eq. (3) should be used to design bi-material DCB joints to achieve the pure mode I fracture. However, they proposed the following equation as the joint design criterion, based on the plane strain condition: where v xz and v zx are the Poisson's ratios where x and zdenote the length and width directions of the DCB specimen, respectively. It is clear that Eq. (4) can be simplified into Eq. (3) when the Poisson's ratio of the two arms is the same or slightly different [14,66].
Research on bi-material DCB joints design configuration following Eq. (3) is not as extensive as the one following the same flexural rigidity. However, the limited research of the former agrees that the mode II fracture component is sufficiently suppressed when Eq. (3) is followed in designing the bi-material DCB joints. The experimental work by Ouyang [14] shows that the crack tip tangential slip is significantly suppressed. The fractographic results by Wang et al. [65] demonstrate the absence of fracture surface features related to the mixed mode fracture process, i.e., the absence of the mode II fracture component. FEA results by Zambelis et al. [67] also show that the mode II SERR is 0 in their bi-material DCB joints.
Apart from the above two DCB design configurations, several researchers investigated the mode I dominant fracture behaviour of bimaterial interface cracking using ADCB specimens [51,68,69] [54,[70][71][72][73][74][75]. The total SERR was calculated in all these studies. While some researchers overlooked the mode II component in their analysis [68,69], the rest did not conduct any analysis pertaining to the mode II SERR component [51,54,[70][71][72][73][74][75]. They overlooked the mode II fracture component as the researchers adopted the Williams theory [68,69]. Williams assumed that pure mode I loading requires two opposite moments of the same value applied to the two arms [47]. In other words, the adoption of the Williams theory leads to neglecting the mode II SERR component for any general bi-material DCB specimens as long as they are loaded as shown in Fig. 1.

Mode II bi-material interface testing
The significance of mode II interface crack problem is justified by the fact that most bonded joints are designed to transmit loads in the form of shear. Extensive research has been carried out to investigate the delamination propagation in laminated composites and interface crack growth between two bonded adherends, leading to few standard test procedures for the three-point bending end-notched flexure (ENF) and end-loaded split (ELS) test configurations [76,77]. The ENF test configuration, proposed by Russel [78,79], is illustrated in Fig. 2 while the ELS test configuration is schematically represented in Fig. 3. The advantage and disadvantages of using ENF, ELS and some other test configurations for mode II fracture testing are summarized in the literature [80,81]. The mode II fracture tests on bi-material joints are summarized in the following.
In contrary to the extensive research on the mode I interface crack problem summarized in the preceding subsection, less work on the mode II interface crack between two dissimilar materials can be found in open literature.
The main concern about using the ENF coupon is how the bi-material ENF coupon can be designed to achieve a pure mode II loading at the crack tip. To the authors' best knowledge, four different joint designs dedicated to obtaining mode II fracture toughness in bi-material assemblies can be found from open literature.
Ning et al. [82] used the standard ENF configuration that comprises two adherends of the same thickness to investigate the mode II fracture toughness of CFRP/Al interface. The configuration was adopted to examine the toughening effects of vapour-grown carbon fibre interleaves. However, it is questionable that this standard ENF configuration provides pure mode II of the bi-material interface. Ning et al. [82] did not provide any information on the mode mixity of the chosen ENF specimen.
Ouyang and Li [83] proposed a bi-material ENF specimen configuration that can be described by the previously presented Eq. (4) for measuring the pure mode II fracture toughness. As has been explained in the previous subsection, the joint configuration can also be described with Eq. (3) for two adherends with very close Poisson's ratios. The development of the bi-material ENF specimen is based on their previous work [14,66] on developing a bi-material DCB specimen for pure mode I fracture toughness. According to Ouyang and Li [14,66,83], Eq. (4) is a decoupling condition such that the mode I loading does not cause tangential deformation components related to mode II fracture behaviour. Besides this, the mode II loading does not cause normal deformation components related to mode I fracture behaviour. They used this specimen configuration to study the influences of shear strength and mode II fracture toughness on the critical load level of ENF testing.
Other researchers intended to design the bi-material ENF configuration as described by Eq. (3) for obtaining mode II fracture toughness experimentally [84][85][86][87]. Inconsistently, they finally designed the  bi-material ENF specimen by matching the stiffness of upper and lower arms as depicted by Eq. (1) and Eq. (2). They did not conduct detailed analysis to prove that mode I fracture component was not presence in their joint design. However, the bi-material DCB specimen designed by matching the stiffness shows the coupling effect, i.e., the mode II fracture component is found in the specimen subjected to mode I loading. It is therefore unlikely that the bi-material ENF specimen designed according to Eq. (2) achieves pure mode II mode at the crack tip.
To overcome the complexity in detecting the location of the crack tip in the ENF testing, Qiao and Wang [88][89][90][91] proposed a tapered ENF configuration for examining the mode II fracture of bonded interfaces between wood and E-glass FRP. Based on the elastic foundation theory [90], Qiao and Wang developed a theory to design the contour of tapered adherends, such that the compliance of the specimen remains constant with crack growth. An illustration of the designed configuration of the tapered wood-composite ENF specimen is given in Fig. 4.
Qiao and Wang [88,89] conducted experimental work to examine the wood-composite fracture toughness under mode II loading. The tested specimen configuration is illustrated in Fig. 4. They claimed that the crack at the FRP/wood interface was under mode II shear loading since the stiffness of the upper adherend comprising red maple and FPR composite was very close to that of the lower adherend. They did not further examine the mode mixity at the crack tip of the proposed ENF specimen subjected to 3-point bending. One should realize that the proposed bi-material ENF specimen is designed by matching the stiffness of two adherends, in spite of the effort to achieve a constant compliance. This design can also be described with Eq. (2). In another paper from Qiao et al. [92], asymmetric four point bending fracture specimens, similar to the hybrid ENF specimen, were adopted for characterizing the fracture toughness of the wood-FRP interface. They considered the interface deformations owing to interfacial stresses and shear deformation by modelling the adherends as shear deformable beams. Liu and Qiao [93] studied the fracture toughness of a GFRP-concrete interface using 4-point asymmetric ENF specimens. The GFRP substrate was stiffened by bonding an additional aluminium strip to achieve a mixed mode fracture in such specimen. Opening of the cracked portion as shown in Fig. 5b has been observed for the 4-point asymmetric ENF specimen in their FEM results. They obtained three mixed mode ratios by changing the stiffness of the hybrid Al-GFRP adherend.
Another piece of significant work in this regard is from Sundararaman and Davidson [94]. The authors discussed the possibility of having the mode I fracture component in a bi-material ENF test. To eliminate the effects of the crack plane contact outside of the right support pin (see Fig. 5) on the load transfer between the two arms and thus mode mixity, the authors proposed the uneven ENF configuration that is illustrated in Fig. 5. They argued that a classic plate theory approach would predict an identical cross-sectional rotation for the two arms at the crack tip, whereas the more rigid upper arm dominated the rotation and the less rigid lower arm would rotate less, from a physical point of view as illustrated in Fig. 5b [94]. The relative rotations of two arms could lead to opening at the crack tip, resulting in a mode I fracture component. Based on a plate theory approach, they enforced two torsional springs for the two arms at the crack tip to account for the effects of shear deformation and the portion of loads carried by respective arms on the mode mixity. However, it is not clear from this paper how the specimen should be designed to achieve a pure mode II fracture behaviour. It is noteworthy that no literature with regard to experimentally obtaining mode II fracture in a bi-material joints discussed so far considers the effects of residual thermal stresses on the mode mixity. When considering the effects of residual stresses due to difference in the thermal expansion coefficients, the relative position of the two arms, i. e., which arm is placed in the upper position, will also affect the mode mixity.
Based on the literature survey on the fracture testing of bi-material interface cracks subjected to mode II loading, the mode I fracture behavior can be induced in such tests by the local opening at the crack tip. Ouyang and Li suggested obtaining pure mode II fracture toughness by subjecting bi-material joints designed according to Eq. (4) to mode II loading [83]. Nevertheless, rigorous proof that the crack tip opening does not occur in such joints under mode II loading is missing. It is worth highlighting that the effect of thermal residual stresses is not considered in the joint design to achieve pure mode II fracture. More research is needed in this regard.

Mixed mode bi-material interface testing
It is crucial to examine the mixed mode failure of adhesively bonded joints. Adhesively bonded joints are normally subjected to mixed mode loading and it is well known that the fracture toughness, G total = G I + G II Fig. 3. Illustration of a generic ELS configuration for mode II fracture tests. Fig. 4. Tapered ENF specimen for mode II wood-composite interface fracture. This figure is adapted from Ref. [90].
, is a function of the mode I, mode II fracture components and their ratio. Even though the majority of research results show that mode II loading can increase the total fracture toughness [40,[95][96][97], some contradictory research results can also be found in the literature [59].
Several test setups can be adopted for mixed mode testing. For studying the mixed mode failure in adhesively bonded joints with two adherends of the same material and symmetric geometry, several specimen configurations and loading scenarios have been extensively used, including the mixed mode bending (MMB) test, the ADCB test, the single leg bending (SLB) test etc. [4,59]. These test setups have been adopted for studying the mixed mode failure in adhesively bonded bi-material joints as well. The key is how the mode I and mode II fracture components (G I , G II ) can be extracted since the geometrical and stiffness asymmetries of the two arms in bi-material specimens influence G I and G II values. Besides, the mismatch of the thermal expansions of two arms contributes to the mode mixity when the bonding is performed at an elevated temperature and residual thermal stresses are generated [98][99][100]. The mode mixity introduced by the asymmetries and thermal expansion differences jeopardizes the validity of pure mode I loading of traditionally adopted DCB specimens, and pure mode II loading of ENF specimen, which has been explicitly shown in the preceding subsections. Different test setups for mixed mode testing of bi-material joints will be discussed hereafter. Fig. 6 depicts a general MMB test setup, which is wildly used for mixed mode fracture testing. The MMB test was first proposed by Reeder and Crews [101,102] and has been standardized for specimens with symmetric arms [103]. Technically, all mode mixities can be achieved with this test setup by changing the loading moment arm. For this reason, the MMB test is widely employed for characterizing crack growth in bi-material joints.
Shahverdi and Vassilopoulos et al. [104][105][106][107] have employed the MMB specimens to study the mixed mode fracture in pultruded glass FRP joints. Specimens with crack paths that were not in the symmetric plane of the bonded joints were investigated. They used a visual extensometer technique to monitor the crack length. More than 20 pairs of black circular dots were evenly located above and below the crack path and their relative displacements were continuously monitored. The crack reached the location of a pair of dots when their relative displacements exceeds 2% of the original displacement. To separate the SERR components for mode I and mode II, they extended the Williams "Global method" following the same assumptions (further described in Section 4.2). A simple yet closed form analytical mode partitioning method was therefore developed.
Arouche and Wang et al. [108] also employed the MMB test to study the crack propagation at bi-material CFRP-to-Steel bonded interface. The specimen was designed according to Eq. (3) (matching longitudinal strains) so that a pair of opposite moments could only generate pure mode I loading at the crack tip. The same assumption for obtaining pure mode II of the Williams method was also followed: pure mode II is obtained by having the same curvature in the two arms. However, the total SERR calculated with the proposed method was less than that obtained with FEM method. This is attributed to the fact that the analytical method in the study neglects the shear energy as a result of transverse loads in the MMB setup.
Bi-material DCB loaded with uneven bending moments (DCB-UBM) is another test configuration used. Similar to the MMB test, the DCB-UBM test offers the mode mixity ranging from pure mode I to pure mode II by varying the applied moments on the two arms [95]. It also permits stable crack growth. To apply different moments on the two arms, a special test fixture is required. Rask and Sørensen [49]   developed a curvature based J-integral approach to characterize the total fracture energy. This method circumvents the inconvenience of determining material properties and the J-integral is obtained by determining the curvatures of the two arms with measured strain distribution. The mode partitioning method of Suo and Hutchinson [38] was used for the partitioning of the mode I and mode II. Dimitri et al. [109] has developed an analytical mode partitioning method for the DCB-UBM configuration by assuming the bondline as an elastic interface. This model is computationally expensive since many differential equations have to be solved in order to derive the energy release rate components.
The single leg bending (SLB) test shown in Fig. 7 is another commonly applied experimental specimen for the mixed mode fracture test of achieve. For symmetric joints, this configuration can achieve a mixed mode ratio of G II /G total = 0.429 [4], while da Silva [59] reported this value to be around 0.465 or mode mixity of around 41 • . Davidson and Sundararaman [110] has demonstrated the applicability of the SLB for testing bi-material interfaces. They pointed out that the SLB specimen essentially provides the same mode mixity range as the notched four point bending test shown in Fig. 9. The mode mixity increases with modulus E 1 /E 2 ratio for the SLB test. Eventually, they performed a parametric study to investigate the range of mode mixity for different thickness ratios of the two arms. Their study indicates that a reasonably wide range of mode mixity could be achieved with this specimen.
The mixed mode end loaded split (ELS) test, shown in Fig. 8, is also popular for mixed mode testing [111,112]. It is noteworthy that the ELS test can also be used for mode II fracture tests as illustrated in Fig. 3. However, the loading is different for the mixed mode ELS test from the loading for the mode II ELS test. The mode mixity of the ELS is roughly independent of the crack length. It has been reported that the symmetric mixed mode ELS specimen also yields a mixed mode ratio of G II / G total = 0.429 [4]. The mixed mode ratio can be altered when changing the relative mechanical properties and thicknesses of the two arms, i.e., in the case of bi-material or asymmetric joints. The mixed mode ratio is almost independent of the crack length. Guo and Weitsman [113] performed a parametric study on the effects of the thickness ratio h 1 / h 2 on the mixed mode ratio. They found that a range of 0.3 < G II / G I < 4 can be achieved with this configuration by changing the thickness ratio. They also found that the mixed mode ratio is sensitive to the adhesive thickness. Qiao and Wang [114] studied the deformations at the crack tip of bi-material ELS specimens using their "interface compliance" method.
Finally, the notched four point bending specimen is another specimen configuration that has been adopted for analysis of mixed mode failure of bi-material joints [110,[115][116][117]. Fig. 9 shows the specimen configuration. Essentially, this specimen setup is equivalent to the SLB specimen [110].
As can be seen, the specimen configurations that have been well established for examining mixed mode failure behaviour in adhesively bonded joints with geometrical and material symmetries have also been extended for testing bi-material joints. The specimen design and the mode partitioning for bi-material joints, however, are different from those for conventional joints. The material asymmetry and geometric asymmetry of bi-material joints require accurate and reliable mode partitioning methods which can be used to guide the design of such test specimens and analyze the results afterwards. The theoretical work pertaining to the crack driving force and mode partitioning is presented in subsequent sections.

Crack driving force
The crack driving force, either SIF, K, SERR, G or J-integral is always pursued to predict the onset and propagation of delamination and disbond growth. Depending on which crack driving force is used for the characterization of a bi-material interface crack, two approaches can be identified: the local approach using K and the global approach using the concept of SERR G or J-integral. K and G are essentially equivalent as elastic fracture parameters and the use of either of the two parameters provides the same information. The difficulties in calculating Kfor inhomogeneous composites have made G a preferable parameter for characterizing damage growth in layered composites. The validity of K and G is confined to crack problems where the extent of the inelastic processes is much smaller than the relevant geometric length scales. Jintegral has to be used when the validity conditions of K or G are not met.
The derivation of these three forms of crack driving force for bimaterial bonded joints will be concisely summarized hereafter.

Stress intensity factors
In the domain of Linear Elastic Fracture Mechanics (LEFM), the SIF concept is frequently applied to characterize the crack growth behaviour. It particularly describes the singular stress and displacement fields in front of the crack tip. Based on the Westergaard method [118], Irwin in 1957 introduced the concept of SIF with a focus on the near-tip stress and displacement fields [119]. This method was soon elaborated by other researchers and regarded as the driving force for static/fatigue crack growth in isotropic homogeneous materials [120]. The calculation methods of SIFs, including opening and shearing modes, for cracks in isotropic homogeneous metallic materials are comprehensively summarized in handbooks [121][122][123]. The derivations of the SIFs, crack tip stress and displacement fields will not be discussed here. It is highlighted that the crack tip stress and displacement fields possess the singularity of r − 1/2 with rbeing the radial distance from the crack tip. This singularity is different from the singularity of a bi-material interface crack.
As has been discussed in the introduction section, Erdogan [17], England [18] and Rice [19] analytically derived stress intensity factors (SIF) for the bi-material interface crack problem between two isotropic homogenous half planes. It is noted that the derived crack tip stress and displacement fields of a bi-material interface crack possess the singularity of type r − 1/2+iε with ε being the oscillatory index. As a result, the calculated near-tip stress and displacement fields show an oscillatory behavior [124]. The mode I and mode II SIFs for bi-material joints will be explained in detail in Section 4.1 in this review. One has to keep in mind that the calculation of K is difficult due to the material asymmetry across the interface crack in bi-material joints.

Strain energy release rate
Unlike the near tip stress field method discussed for K, the concept of energy release rate associates the crack growth with the energy variation of the system containing a growing crack. The total SERR, G total , is defined as the decrease in total potential energy, Π, during unit crack extension: The total SERR, G total , is an overall quantity which can be decomposed into a mode I component, G I , and a mode II component, G II . The process of decomposition is the so called mode partitioning which will be addressed in Section 4. Nevertheless, the calculation of the total SERR is essential in the fracture analysis of a bi-material interface fracture and is addressed briefly in this section. The total SERR can be calculated by two means: Irwin's crack closure integral [46] and the compliance method [47]. The former method lays the cornerstone for the numerical VCCT method which has been extensively applied to calculate the total SERR and perform mode partitioning for adhesively bonded joints. For the essential knowledge of the VCCT and its applications in fracture analysis of composite structures, adhesively bonded joints and particularly bi-material joints, one is referred to the literature [125][126][127][128][129][130][131][132]. This section mainly focuses on the compliance method. Particularly the material and geometry asymmetries of the bi-material interface crack should be taken into consideration when using the compliance method.
The Irwin-Kies method [133] is used to calculate the SERR in the compliance method: dC da (6) with Cbeing the compliance of the specimen containing a crack of length a. P denotes the applied load and B is the width of the specimen. It is noteworthy that different methods exist for the calculation of the compliance.
The most intuitive method is to experimentally derive the function C = f(a) and then differentiate the function with respect to the crack lengtha. Usually the function takes the following two forms: where C i=0,1,2,3 are constants that can be derived by curve fitting either Eq. (7) or Eq. (8) to experimental compliance data. In Eq. (7), the exponent nis usually 3 [4]. This method requires the continuous monitoring and recording of loadP, displacementδand crack lengtha during testing. It is not difficult to discover that the highest exponent in both function forms being 3 is in accordance with theoretically obtained compliance functions of a beam theory. Theoretical analysis using the beam or plate theory is also extensively applied to derive the functionC = f(a). The easiest way is to apply the simple beam theory to calculate the deflection of the cantilever beam of lengtha. It is worth noting that the compliance of each adherend of the bi-material joint has to be calculated in order to obtain the compliance of the whole system [67]. The shear deformation is also considered in some cases, but it is argued that the shear modulus of a laminated composite adherend is much less than the longitudinal modulus and thus the shear deformation cannot be neglected. In this regard, the Timoshenko beam theory has been extensively applied.
In order to take into account the influence of adhesive on the compliance of the beam-like specimen containing a crack, the deflection of a beam resting on an elastic foundation is derived and the compliance function can then be easily obtained [134,135].

J-integral
The LEFM has enjoyed great success when the plastic deformation in front of the crack tip is much smaller in comparison with other characteristic dimensions. The applicability of LEFM is violated when the plasticity has to be taken into account. Rice [136] developed the path independent integral scheme to calculate the energy release rate taking plasticity ahead of the crack tip into consideration. The J-integral is expressed mathematically as following for a two-dimensional case shown in Fig. 10a: where Wis the strain energy per unit volume, T the traction vector and Γan arbitrary contour starting from the bottom of the crack and ending  at the upper end of the crack. uis the displacement in the x direction. J = 0 for the case the contour does not contain a crack face [137]. The J-integral is equivalent to SERR when considering only elasticity [126]. For a beam-like joint, the contour is typically chosen such that dy = 0 or T = 0for some segments of the contour to facilitate the analytical derivation of Eq. (9) [126]. The typical contour is schematically shown in Fig. 10b. Cui has derived the J-integral for DCB specimens under mode I loading [138], while the J-integral for ENF specimens under mode II loading can be found [139,140].

Mode partitioning
Inherently an interface crack is subjected to a mixed mode loading condition. It has been proven by extensive experimental work that the fracture toughness is a function of mode mixity [40,59,124,141,142]. A general practice is to extrapolate a limited number of mixed mode results to obtain the total fracture energy of any mode combinations: the critical envelop method. As discussed in the preceding sections, the total SERR is relatively easy to calculate. On the other hand, it is troublesome to partition the total SERR into its mode I and mode II components. In this section, the analytical mode partitioning methods are summarized.

Complex stress intensity factor method
The pioneering work pertaining to the analysis of an interface crack between two isotropic but dissimilar materials was carried out by Rice, Hutchinson and Suo et al. [19,38,41,50,123,143,144]. They used the concept of complex SIF. The nomenclature of this type of analysis is provided in Fig. 11. The overall solution is provided here.
Rice [143] and Hutchinson et al. [144] proposed the complex SIF concept for the bi-material interface crack problem, which is expressed as: It is noteworthy that K 1 and K 2 are the real and imaginary parts of the complex SIF. They are not mode I and mode II SIFs respectively. Based on the complex SIF concept, the shear stress component, σ xy , and normal stress component, σ yy , at the distance r in front of the crack tip can be asymptotically expressed as: (11) where ε is the oscillatory index, and it is related to the Dunder's parameter [38,145]. The relative displacements behind the crack tip are given by with c 1 and c 2 being two compliance parameters given in Ref. [38]. Based on the near tip stress and displacement fields given by Eq. (11) (12), a relation between the total SERR, G, and the complex SIF, K, is established: After the total SERR is derived for the given loading case, the complex SIF can be analytically solved for, expect one angular quantity ω. The determination of ω is carried out by integral equations.
Strictly speaking, this approach is applicable for interface cracks between two elastic homogenous layers only. The determined complex SIF gives the near tip stress and displacement fields for this type of crack problems. As it has been pointed out by many researchers, the determined crack opening has unrealistic crack surface penetration immediately after the crack tip using this method.
The applicability of this approach for interface cracks with inhomogeneous layered composites involved is highly debated. The near tip stress and displacement fields may be not accurate for this type of crack problems.

The Williams partitioning method
As shown in Fig. 12, Williams [47] has made two assumptions for the partition of the mixed mode: 1) pure mode I loading requires two opposite moments of the same value applied to the two arms, β = 1 in Fig. 12; 2) pure mode II is obtained when the curvature of the two arms is the same, resulting in ψ = (1 − ξ/ξ) 3 with ξ = h 1 /(h 1 + h 2 ). Following the assumptions, the applied moments can be decomposed into two moments corresponding to pure mode I and pure mode II loading cases: (14) Substitute the above equation into the equation for G total , and one should be able to decompose the total SERR into the following where E and I are the equivalent bending flexural modulus and moment of inertial of the bonded portion of the beam. Shahverdi and Vassilopoulos et al. [104,105] followed the two  Fig. 11. Nomenclature for an interface crack between two isotropic materials.
assumptions in the Williams method, but argued that the curvature of the beam relies more on the flexural stiffness than the thickness of each arm. Consequently, ψis expressed by Eq. (16). Then G I and G II can be obtained by following the same procedure as in the Williams method.
Schapery and Davidson [141] has proven the assumption made in the Williams method for obtaining pure mode I SERR not accurate. It is only reasonable for the symmetric case [146]. It has been shown in section 2.1, the interface crack between two bonded dissimilar materials under mode I loading is a mixed mode problem, unless the joint is designed deliberately such that the longitudinal strain distributions of the faying surfaces [65].
To overcome the identified deficiency of the Williams method, Arouche and Wang et al. [108] proposed another partitioning method by modifying the condition for mode I SERR. This approach was motivated by the strain-based design of bi-material DCB for pure mode I loading [65]. As illustrated in Fig. 12, in order to have the same strain distribution at the faying surfaces of the two arms for the pure mode I SERR, β is expressed as

The Davidson method
Davidson and his co-workers have published a series of research papers on the crack tip element (CTE) method for mode partitioning of crack/delamination at the bi-material interface [94,110,111,141,[147][148][149][150][151]. The crack-tip element approach was first proposed by Schapery and Davidson in 1990 [141]. In Fig. 13, the nomenclature and loadings on a crack-tip element are depicted. In the context of classic plate theory, concentrated loads and moments act at the orthotropic axes of three sections: cracked upper arm, cracked lower arm and perfectly jointed beam. The crack tip is assumed to be represented by a concentrated crack tip shear force N c and momentM c , as shown in the close-up in Fig. 13, and no tractions act along the crack prolongation ahead of the crack tip.
The total SERR is derived by using the concept of Irwin's VCCT, which is given by Eq. (18). The calculation of N c and M c is essential and comprises mainly two steps: 1) consider the uncracked portion as one laminate and calculate the force and moment resultants in the upper and lower lamina by employing the classic plate theory; 2) calculate N c and M c by employing the equilibrium of the upper arm with the externally appliedN 1 ,M 1 and calculated force and moment resultants known. In Eq. (18), subscript 1 and 2 denote the upper and lower arms respectively, and u and θrespectively are the displacement and rotation at the crack tip due to N c and M c for a crack increment of Δa. They can be calculated by using the classic plate theory.
Rearrange Eq. (18), one can get the following equation: where c 1 ,c 2 and c 12 are known parameters related to the material properties and geometry dimensions. Based on a mathematical observation, Eq. (19) can be decomposed into two SERR components corresponding to mode I and mode II in the following manner [141]: The decomposition can be finished by deriving k i=1,2,3,4 parameters. By matching the coefficients forN c , M c and N c M c in Eq. (19) and Eq. (20), three equations can be obtained for the four unknowns. It is apparent that one more condition is required to resolve for the fourk i parameters. An FEM analysis of a simple loading case was proposed by Schapery and Davidson [141] to gain the additional equation.
In 1995, Davidson et al. [147] correlated the CTE approach and the complex SIF method for the mode partitioning. They expressed the complex SIF as a function of the internal N c and M c instead of external loads as developed by Suo [38]. Similar to the method proposed by Suo, the "mode-mix parameter" Ω cannot be calculated based on the CTE approach.
Davidson et al. [147] suggested that an FEM analysis of the CTE with one special loading case needs to be performed to obtainΩ. In addition, the calculation of the complex SIF in this method involves the choice of a characteristic dimension Lon which the SIF is based. Discussions on how to choose a special loading case for an FEM analysis and a value for Lto the simplify the calculation are provided in the literature [147]. This method has been extended for 3D problems [148,149].
It is evident that the CTE approach requires the choice of a characteristic dimension and FEM analysis, making it troublesome for analyzing the failure of bi-material interfaces. The validity of SIF in laminated FRP structure is also questionable, especially for cases where the fracture process zone is higher than the lamina thickness.

The Qiao & Wang method
Qiao and Wang and their co-workers extended the respective methods proposed by Suo [38] and Davidson [147] to include the contributions of shear deformation and near crack tip deformations to the total SERR [92,93,114,146,[152][153][154][155][156][157]. They categorized the CTE methods developed by others and themselves into three models, see Fig. 14. The ones proposed by Suo, Davidson and Williams, in which the uncrack portion of the CTE is modelled as a classic laminate where the cross-section at the crack tip remain one plane, are the "rigid joint model". The model developed by Qiao and Wang in which the shear deformation is considered is called the "semi-rigid joint model". A further step was taken by Qiao and Wang to model the crack tip deformations as a result of interfacial stresses, leading to the so-called "flexible joint model".
The semi-rigid joint model is to improve the rigid joint model. In the rigid joint model, the shear deformation cannot be considered as a result of the nature of the classical beam/plate theory. As a result, the SERR of the interface cracks in beam like coupons is underestimated, which is more pronounced for FRPs [158]. This shortcoming was overcame by modelling the two adherends bounding the crack plane as first-order shear deformable plates by Qiao and Wang. The kinematics of the Mindlin-Reissner plate theory, the equilibrium condition of the upper and lower plates and the displacement continuity equations along the interface were used to solve for the axial force, transverse shear force and bending moment distributions in the uncracked portion of the CTE. The force equilibrium conditions of the upper plate in the immediate vicinity of the crack tip were adopted to derive a concentrated normal force N c and a shear force Q c at the crack tip [146,153].
For obtaining the SERR of the CTE, either J-Integral [153] or Irwin's crack closure integral [146] has been employed. And the mode I and mode II SERRs are expressed as eventually: with δ N and δ Q mathematically determined and denoting the compliances of the CTE underN c and Q c respectively. It is noteworthy that the displacement continuity at the bond interface is described as Eq. (24), no restraint on the relative rotation angles of the two plate at the crack tip. The upper and lower plate can have different rotations at the crack tip as shown in Fig. 14b. Whereas the cross-section at the crack tip in the rigid joint model remain a plane and thus the upper and lower plates have the same rotation, as shown in Fig. 14a.
In the "flexible joint model", not only the relative rotations of the two plates, but also the deformations at the crack tip due to interfacial stresses are relaxed [114]. The displacement continuities along the bond interface are expressed as Eq. (25). In comparison to Eq. (24), this joint model includes the deformations induced the normal stresses and shear stresses at the interface. C si and C ni are the compliance parameters in their notations. i = 1, 2 denoting the upper and lower plates. Estimations of the two compliance parameters are also provide by the authors. With these deformations included, Qiao and Wang managed to obtain the interfacial stress distributions instead of concentrated forces at the crack tip. Further improvement is achieved.
In later work from them [157,159], the thermal residual stresses and hydrothermal influence are also taken into consideration. As the approach of Davidson, Qiao and Wang extensively employed the Suo method for mode partitioning in addition to their proposed global partitioning method.
Luo and Tong [160] developed a mode partition method similar to that of the semi-rigid joint, whereas they took the effects of transition near the crack tip into the displacement continuity. They expressed the continuity in the following equation: where φ 1 = (φ 1 +φ t )/2 and φ 2 = (φ 2 +φ t )/2 and φ t is the curvature of the bonded region.

The Wang and Harvey method
Wang and Harvey have also published a series of papers on mode partitioning of bi-material interface fracture [161][162][163][164][165]. Using two orthogonal mode vectors, the authors wrote the total SERR in terms of orthogonal modes, which is given by Eq. (27). The orthogonality condition actually means that no cross product terms are generated.
In Eq. (27), α θ and α β are two mode partition coefficients respectively corresponding to two orthogonal vectors. Two types of pure modes were assumed: 1) Pure mode I was defined when crack tip relative shearing displacement was zero or pure mode II was defined when crack tip opening force was zero. 2) Pure mode I was defined when crack tip shearing force was zero or pure mode II was defined when the crack tip relative opening displacement was zero. Under either condition, α 2 θ G θ denoting mode I SERR and α 2 β G β denoting mode II SERR can be calculated. Wang and Harvey used both Euler and Timoshenko beam theories to perform the calculations. Compared to the mode I SERR and mode II SERR obtained using the Timoshenko beam theory, Euler beam theory provides a higher mode I SERR and a lower mode II SERR. They argued that the two theories acted as upper and lower bounds of the mode partition and the practical partition could be the average of the two partition theories. Equivalently, the Euler beam theory provides the mode partition similar to the Davidson method and the Timoshenko beam theory results in mode partition similar to the "semi-rigid joint" of the Qiao and Wang method. Wang and Harvey argued that the mode I SERR should be the average of the G I from the Davidson method and the G I from the Qiao and Wang method and the same for mode II SERR.
However, it is noteworthy that over certain range of loading conditions the Wang and Harvey method provides unrealistic results where either mode I or mode II SERR is negative [166].

Interfacial stresses based methods
The discussed methods above focus on the calculation of loads and displacements in the vicinity of the crack tip for the mixed mode partitioning, whereas the methods that will be discussed in the following calculate the interfacial stresses along the crack plane bounded by bimaterial adherends. The normal and shear interfacial stresses are calculated in these models to derive their contributions.
The elastic interface model has been considered widely by many researchers [57,109,158,167,168]. In this method, the problem is represented by two sub-laminates partially bonded with an elastic interface which exerts uncoupled normal and tangential stresses. A set of differential equations are deduced and need to be solved, sometimes numerically, to obtain the interfacial stresses. By modelling the two sub-laminates using the classic plate theory or first-order shear deformable plate theory, Bruno and Greco [158] derived closed-form mode partition and they argued that considering the shear effects at the crack tip was important in calculating the mode partition. Bruno and Greco [167] further refined the interface method by modelling each adherend as a multi-layers laminate. Valvo et al. [57,168] also adopted the elastic interface model for mixed-mode delamination problems, they provided a different way of solving the differential equations. Dimitri et al. [109] used the method to analyze the mixed-mode delamination in moment -loaded DCB composite specimens.
Instead of using the elastic interface model, Guo and Weitsman [113] calculated the shear stress and normal stress by adopting the weight function approach and the Timoshenko beam theory.

Discussion
Even though many mode partitioning methods have been developed, they inadequately address the effect of thermal residual stresses on the fracture analysis of bi-material interface cracks. Besides, it is also crucial to understand the limits of the proposed mode partitioning methods. These two topics are covered in this section to complement the discussion on the fracture analysis of bi-material interface cracks.

The role of thermal residual stresses
Only few scientific research work has considered the effects of residual stresses in the interface crack analysis of bi-material joints [1,64,[98][99][100]159,[169][170][171][172][173][174]. Residual thermal stresses are developed at the interface in bi-material joints made of two dissimilar materials with different thermal expansion coefficients when the service temperature differs from the manufacturing temperature. The developed residual stresses result in curved beam specimens and are released once the interface crack propagates. They contribute to the total SERR of the interface crack and mode mixity of the bi-material interface failure. It is crucial to take the effects of residual stresses into consideration.
The total SERR of the bi-material interface crack with residual stresses can be broken into three components according to the literature [99,100,169]. It can be mathematically expressed as Eq. (28). G m is the SERR resulting from the mechanical loading, G th is the SERR resulting from the thermal loading and G ad is the SERR term that indicates the work done by the mechanical forces in the residual displacement field.
G total = G m + G th + G ad (28) Consequently, the total SERR cannot be simply calculated by superimposing the two SERR items due to mechanical loading and thermal loading. There is a coupling item G ad which is dependent of how the mechanical loads and thermal loads interplay.
Another important factor is that the thermal residual stresses in each arm of the bi-material joints can result in contact of the two beams. The resultant contact loads also contribute to the total SERR by affecting the loading condition of the crack tip in addition to the externally applied loads and thermal loads [100,169,170]. Additionally, the contact loads also alter the mode mixity of the interfacial crack by altering the loads and moments in the two arms near the crack tip.
As for the mode partitioning, the semi-rigid joint in the Qiao and Wang method has been extended to account for the residual stresses effects [159,173]. A crack-tip element is typically adopted to describe the bi-material interfacial crack problem in this method. The resultant loads and moments in the arms are calculated with the shear deformable beam theory, the residual stresses are also included. The effect of the residual stresses on the mode partitioning is therefore considered in this method in the calculated crack-tip concentrated forces. One needs to be cautious when contact between two arms due to residual stresses happens, the contact forces have to be included in the calculation of the resultant loads and moments of the crack-tip element [170]. However, G ad in Eq. (28) has yet been considered in the mode partitioning methods.

The validity of the mixed mode partitioning methods
As it has been reported by many researchers, while the total SERR calculated by different mode partitioning methods yield similar values, the predicted mode mixities are different [160,161]. This difference can even be significant [175].
The validity of the Suo method is limited to cases where the fracture process zone should be much smaller than the K-dominance zone. The applicability of the Suo method for laminated composite specimens becomes questionable due to the validity of K-dominance [161][162][163][164][165]. Similarly, the Davidson method has the same problem associated with its validity. Several researchers also indicated that the Williams method has a defective assumption regarding the calculation of mode I SERR [141]. The Qiao and Wang method represents the crack tip with concentrated forces, which may misrepresent a crack with a large fracture process zone.
Experimental efforts have been made in order to indicate the applicability of the Suo method and the Williams method [176,177]. They are inclined to conclude that the Williams method is more applicable for laminated composite specimens. No further systematic experimental investigation has been carried out to compare the rest of the proposed mode partitioning methods. A predominant challenge is that the total SERR can be experimentally obtained while the mode I and mode II SERR components have to be calculated with one of the proposed methods. However, none of the theories has been proven highly accurate or has a very sound physical basis. A benchmark is lacking.

Conclusions
The difference in Young's moduli of two adherends and material asymmetries of bi-material interface cracks lead to the deficiency of employing the experimental methods and analysis tools that have been well established for crack problems without material asymmetry.
However, the bi-material interface normally represents a weak plane and is prone to failure and hence needs attention. Extensive research work has been dedicated to analysing the bi-material interface crack problems. The key to analysing bi-material interface crack problems is to partition the total SERR into its mode I and mode II components accurately and reliably.
While it is well known that the fracture toughness of an interface crack depends on the mode mixity, it is troublesome to perform the mode partitioning for the bi-material interface crack. Inherently, this type of cracks undergo mixed mode loading conditions. This poses problems in terms of experimentally determining pure mode I and mode II fracture toughness values and mode mixities for mixed mode tests. In most cases only mode I dominant or mode II dominant fracture toughness can be obtained experimentally. The mode mixity of a mixed mode test relies on how to perform the mode partitioning with one of various proposed theoretical work frames. It has to be pointed out that most experimental work has not considered the effects of residual stresses.
This review summarizes many theoretical work frames developed for the bi-material interface crack problem. No consensus has been reached on which theoretical work provides the best prediction results for bimaterial interface crack problems. The latest research work suggests that the accuracy of the prediction models depends on the damage process behaviour of the interface crack. More research in this area is needed.

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.