# Fileset

[Direct Measurement and Verification_submission.pdf](https://mdr.nims.go.jp/filesets/a1e1fa5f-f625-4b5d-b56d-2719c5b6895c/download)

## Creator

[Yuta Tobata](https://orcid.org/0000-0003-0730-5934), [Kimiyoshi Naito](https://orcid.org/0000-0002-3334-4876), [Jonathon Tanks](https://orcid.org/0000-0002-0232-8240)

## Rights

This is the pre-peer reviewed version of the following article:Direct measurement and verification of cohesive zone model parameters for basalt/PP rods using the transverse tensile test and virtual double cantilever beam test, which has been published in final form at https://doi.org/10.1002/pc.27145. This article may be used for non-commercial purposes in accordance with Wiley Terms and Conditions for Use of Self-Archived Versions.[In Copyright](http://rightsstatements.org/vocab/InC/1.0/)

## Other metadata

[Direct measurement and verification of cohesive zone model parameters for basalt/PP rods using the transverse tensile test and virtual double cantilever beam test](https://mdr.nims.go.jp/datasets/aed10430-205a-445c-8859-3304b7ac5336)

## Fulltext

Microsoft Word - Direct Measurement and Verification_submission.docxDirect Measurement and Verification of Cohesive Zone Model Parameters for Basalt/PP Rods using the Transverse Tensile Test and Virtual Double Cantilever Beam Test Yuta Tobataa*, Kimiyoshi Naitoa,b and Jonathon Tanksa a Research Center for Structural Materials, Polymer Matrix Hybrid Composite Materials Group, National Institute for Materials Science (NIMS), 1-2-1 Sengen, Tsukuba, Ibaraki 305-0047, Japan b Tohoku University, Department of Aerospace Engineering 6-6-01, Aramaki Aza Aoba, Aoba-ku, Sendai, Miyagi, 980-8579, Japan Email: TOBATA.Yuta@nims.go.jp and yuta.tobata@community.isunet.edu  Abstract To further our understanding of failure in thermoplastic-based composite rods, the mode Ⅰ-governed cohesive zone model parameter of a basalt fiber reinforced polypropylene composite rod is evaluated. The relations between stress and strain, the maximum normal traction and final separation displacement are examined using a transverse tensile test method. The measured traction is verified by three-dimensional finite element analysis of rod-shaped double cantilever beam specimens. The numerical results show the close relation between the crack length and cube root of compliance, and yields a representative value of the mode I-governed critical fracture energy.     Graphical abstract Key word Cohesion, Finite Element Analysis (FEA), Tension, Toughness  1. Introduction Fiber reinforced polymer (FRP) composite rods have been utilized for the reinforced concrete (RC) as an alternative reinforcement to the carbon steel because of its superior specific tensile strength and erosion resistance1-6, particularly glass fibers. Meanwhile, basalt fiber has drawn attention due to having composition, properties, and cost comparable to glass. Kessler et al.7 revealed that the tensile strength and Young’s modulus of the basalt fiber reinforced polymer (BFRP) composite rods are greater than those of a glass fiber reinforced polymer (GFRP) composite. Sharma et al. also reported that Taguchi’s design of experimental method showed better wear resistance of BFRP than that of GFRP.8 Ranganathan et al. revealed that adding basalt fibers to vinyl ester reinforcing with silane treated glass fibers indicated good flexibility of the interphase between fiber and matrix during impact loading.9 Lu et al.10 conducted a comparative study of the tensile strength retention of the vinyl ester reinforced with basalt and glass fibers after immersion in alkali solution, seawater and tap water. However, the study reported that the strengths of BFRP are weaker than those of the GFRP at all the immersion conditions. Jiang et al. reported the same tendency regarding BFRP using epoxy matrix.11 The research also indicated that the degradation derives from water absorption of the fiber. Consequently, water resistant matrix can utilize the potential of a basalt fiber. It is known that PP has low water absorption rate, and superior alkali and acid resistance.12-14 In addition, Poornima reported that adding Maleic anhydride improved damping effect of the basalt fiber reinforced polypropylene composite (BF/PP), and interfacial interaction between the fiber and matrix, which can further improve its mechanical properties.15 It also indicated that chlorine sustains water absorption. Therefore, it is anticipated that BF/PP rods have a potential use for the platform of the offshore wind turbine.16  To that end, experimental and numerical data guaranteeing their reliability are essential. Tanks et al.17 extensively reported the creep and fatigue behavior of the BF/PP rod at room temperature. However, the relation between damage characteristics and strength remains unknown. An accurate fracture prediction requires the deeper understanding of the fracture mechanics. Referring to studies of RC using BFRP rods, the single edge notch bending (SENB) test 18-22 and impact test 23 were conducted for the evaluation of the fracture energy. Liang et al.22 conducted a cohesive zone model (CZM)-based finite element analysis (FEA) for the verification of the measured fracture energy, and predicted mechanical responses of the SENB specimen during the crack propagations. CZM is an effective analysis method of predicting the fracture, requiring the fracture energy and the maximum traction. 24-29 Thereby, these parameters should also be obtained for the deeper understanding of the fracture mechanics of BF/PP rods as well as the case of RC. A typical evaluation method of a fracture energy of mode-Ⅰ is the double cantilever beam (DCB) test.30-35 Many studies conducted the DCB test for the rectangular specimen of the laminate composites. On the other hand, there is no report using the test for the rod shape specimens. Since the normal stress distribution in the radial direction around the crack tip differs from that of rectangular specimens, cracks may begin to propagate from the center region in the rod, while no crack propagation can be confirmed on surface using conventional measurement techniques. Considering this different crack propagation behavior, an evaluation method for the rod shape composites had been expected as well as the case of the laminate composites. The authors previously proposed a method of evaluating the mode Ⅰ-governed critical fracture energy SEⅠC for the rod shape composites by a modified DCB test. 36 The SEⅠC is a representative fracture energy parameter. However, the other necessary parameter, the maximum normal traction, of the BF/PP rod still remains unknown. Basically, the tensile strength of the transverse direction is lower than that of the longitudinal. Proper adhesive selection possibly enables the transverse tensile load to a rod without interface delamination between jigs and rod. Because of the round shape geometry, however, it is challenging to define the representative stress and strain of the composite rod. Accurate measurement of the stress-strain relation is also expected to measure the transverse Young’s modulus for the FEA. A direct transverse tensile testing method for the composite rods needs to be constructed. Furthermore, the previous study only experimentally verified the SEⅠC by comparing to that of the neat PP because the trigger of the entire fracture seemed to derive from the crack propagation in the resin area by considering the fracture observation. The verification seems to be insufficient. A numerical verification is desirable to support the effectiveness of the modified DCB testing method.  This study involves transverse tensile testing of the BF/PP rod for measuring the maximum fracture stress (maximum normal traction). In addition, the measured maximum normal traction was also verified by a three-dimensional finite element method (FEM) model. The estimated relation between the crack length and cube root of compliance, and SEⅠC were compared to those of the previous study.   2. Material and methods 2.1. Material and specimen Heat-resistant basalt fiber (Nakagawa Sangyo, Co., Ltd., Japan) and maleic anhydride-compatibilized PP (Prime Polymer, Co., Ltd., Japan; Mitsui Chemicals, Co., Ltd., Japan) were used. BF/PP rods were fabricated at the Kanazawa Institute of Technology. The diameter of BF/PP rods was 4.4-4.5 mm. The fiber volume content calculated from the cross section was 44 %. More details are described in the previous reports.17, 36 The specimen is composed with a 20 mm length BF/PP rod and two glass fiber reinforced plastic (GFRP) jigs as shown in Fig. 1. A BF/PP rod and GFRP jigs were attached with acrylic adhesive (3M™ Scotch-Weld™ Structural Plastic Adhesive DP8005, 3M Co., Ltd., USA), which has good adhesion strength to the PP.   Fig. 1. Dimensions of specimens  2.2. Testing method Tensile load was applied with a tensile testing machine (EZ-LX, SHIMADZU Co., Ltd., Japan) under crosshead displacement control (1 mm/min) as shown in Fig. 2. The test was repeated on three samples. The total displacement including that of the GFRP jigs and BF/PP rod was measured from a video type extensometer (TRViewX, SHIMADZU Co., Ltd., Japan), by attaching two CCD markers on the GFRP jigs (25 mm away from the central axis of the rod) as shown in Fig. 1 and Fig. 2. The coordinate system is also drawn in the figure. y is a longitudinal direction of a rod. z is a loading direction. x is an orthogonal axis to the others. The strains of GFRP were measured using two strain gages (5 mm gage length, KFGS-5-120-C1-11, Kyowa Electronic Instrument Co. Ltd., Japan) and averaged. Strain distribution of the BF/PP rod was also measured by digital image correlation (DIC) method. In-situ deformation behaviour of the surface was observed with a digital microscope (VHX-6000, KEYENCE Co., Ltd., Japan) as shown in Fig. 2, and image processing was conducted with a commercial software (VIC-2D, correlated SOLUTIONS, Inc., USA). Finally, the strains were averaged.    Fig. 2. Configuration of the transverse tensile test  2.3 Cohesive zone model parameter The maximum normal traction σnt, max was measured from the maximum applied stress of the stress-strain relation of the transverse tensile test. The bilinear traction-separation law proposed by Alfano et al.27 was applied. The final separation displacement δn, sep was calculated by using the fracture energy SEIC measured in the previous study.36 SEIC would be equivalent to the energy release rate if the plastic deformation did not occur.   2.4. Finite element analysis and material model     Fig. 3 shows a typical three-dimensional symmetrical half cut material model of double cantilever beam (DCB) specimen of BF/PP rod and the coordinate system. Y-axis is a direction parallel to the longitudinal axis of BF/PP rods, Z-axis is a direction parallel to the cross section, and X-axis is an orthogonal axis to the others. Linear elastic analysis was conducted using commercial FEM software (ANSYS mechanical 2021 R1). The representative total length and diameter of a BF/PP rod are determined to be 80 mm and 4.5 mm respectively. The notch represented as a sharp crack was placed in the center of the model (Fig. 3). Four models with different notch lengths (5 mm, 10 mm, 15 mm, and 20 mm) were built. Tensile load was applied to the top of the upper beam while the bottom point of the lower beam was fixed (Fig. 3). The Y-Z cross section was supported in the X direction (frictionless in the plane of the crack). The Z-X cross section was supported so that the Y direction displacement of its centroid is zero (zero remote displacement was applied in ANSYS). 8-node brick isoparametric elements were used for the bulk geometry, and CZM was carried out using a bilinear traction-separation law followed by the procedure of ANSYS.37 Cohesive elements were applied to the interface of the notch front area using an 8-node brick isoparametric element as indicated with an arrow in Fig. 3.   Fig. 3. FEM material model.   3. Results and discussions 3.1. Stress-strain relation and CZM parameter     Strains measured from the method proposed in this study are calculated from the following equations;  𝜀!"# =∆𝑢$"$%& − ∆𝑢'()*𝐷 (1) ∆𝑢'()* = 𝜀'()* ∙ 𝐿+, (2) Δutotal is the total displacement including a BF/PP rod and GFRP jigs measured from the video extensometer. ΔuGFRP is a displacement of GFRP jigs. D is a diameter of a BF/PP rod. εGFRP is an average strain measured from strain gages. Lmk is an apparent marker distance (45.5 mm) provided by subtracting the rod diameter (4.5 mm) from the marker distance (50 mm). The stress was calculated by dividing the load by the x-y cross section area.   Fig. 4(a) shows the stress-strain relation measured from the proposed method (solid line) and DIC (plots) method. The result of DIC was obtained by averaging the strain distribution. Both methods apparently show the similar stress-strain relation in linear region below 0.5 MPa (see. magnified image of stress-strain relation in Fig. 4(a)). The distribution measured by DIC at 0.5 MPa was homogeneous as shown in Fig. 4(b). Consequently, the both strain measurement methods seem to be effective. Above the 0.8 MPa, however, the local higher strain field (region “A” depicted in Fig. 4(c)) was confirmed, and caused greater strain comparing to that measured from the proposed method. This strain gap seems to derive from the crack propagation as indicated with arrows in Figure 4(d). Generally, crack propagations disturb the accurate strain measurement of DIC. Consequently, the measured strain seems effective for the measurement of the stress-strain relation because it can be verified in the linear region by DIC. As a result, the transverse modulus ET was measured 1.1±0.04 GPa. The variance is around 10 % of the average. Therefore, the proposed testing method seems to be able to measure effective ET.            Fig. 4. Comparison of the stress-strain relation between proposed method and DIC. (a) stress-strain relation (b)-(d) strain distribution [%] measured by DIC at 0.5 MPa, 0.8 MPa and 2.6 MPa respectively.  Fig. 5(a) shows an optical microscope image taken at the maximum stress. The crack propagated in the rod, while there is no delamination in the interface between the BF/PP rod and adhesive. Fig. 5(b) depicts a surface image of the specimen after the total fracture. Fiber bundles attached on the both sides of the separated specimen. This also suggests that cracks properly propagated in the BF/PP rod during the test. Consequently, the transverse tensile test conducted in this study seems to provide the effective maximum normal traction σnt from the maximum applied tensile stress. As a result, σnt was measured 2.8±0.3 MPa from the maximum subjected stress. Therefore, it is also suggested that the proposed testing method can measure effective σnt. Applying the bilinear law 27 so that the area of stress-displacement relation is equivalent to the SEIC (1.5 kJ/m2) of the previous study 36, the maximum separation displacement δn, sep was estimated to be 1.1 mm.   Fig. 5. Optical micro scope surface images of BF/PP rod. (a) at the maximum stress (b) after the fracture.   3.2. Numerical verification     CZM parameters were verified with FEM by comparing the relation between crack length and cube root of compliance, and SEIC. The previous study applied preloading, and evaluated load-displacement relations of the reloading. Basically, the preloading decreases the accumulation of permanent strain, and causes almost linear stress-strain relation in reloading. Actually, Ogawa et al. 38 succeeded in estimating the relation between crack length and cube root of compliance of the DCB test for the polyurethane adhesives inserted between two steel plates without considering material nonlinearity. Thereby, this study also did not consider the material non-linearity effect. Geometric non-linear effect was also not considered because the straight-line shape material model (See. Fig. 3) does not seem to cause large strain and rotation deformation. The crack length was measured from the straight-line distance between the crack tip and load axis, following our previous studies.36, 38  Material properties used for the numerical analysis are listed in Table 1, which reflect the isotropic and orthotropic cases. Poisson’s ratio ν of each direction is assumed 0.3 for the simplification. Longitudinal shear modulus (GXY and GYZ) was estimated to be 0.20 GPa from the previous study.38 Transverse shear modulus GXZ was estimated to be 0.23 GPa from the following equation, considering the isotropy in the X-Y plane of the FEM model: 𝐺-. =𝐸/2(1 + 𝜈) (2) 3.2.1. Crack length vs cube root of compliance The relation between crack length and compliance of the preloading was evaluated because the relation does not change until crack propagation begins. The compliance was measured by dividing the displacement by the load where a linear correlation between the two quantities is observed. Fig. 6 presents the relation between the crack length and cube root of compliance of FEM. The experimental results of the previous report 36 were also plotted in Fig. 6 (hollow circle). The square plots present the result of the isotropic model. The solid circle plots present the result of the orthotropic model. The results of the isotropic model are close to those estimated from the beam theory, but differ from the experimental results. Our previous study speculated that the reason of the difference comes from neglecting the orthotropy.36 In the present study, the result of the orthotropic model shows good accordance with the experimental results. Therefore, the FEM model considering the orthotropy seems to reproduce the relation between the crack length and cube root of compliance.  In the previous report, the experimental relation between the crack length and cube root of compliance was expressed by the following linear equation.36  𝐶0/2 = 𝐴0𝑎 + 𝐴3 (3) The relation estimated from the FEM model was also approximated. As a result, A1 was calculated as 1.8 [×10-3 / (m2N)1/3], and A0 was calculated as 0.024 [(m/N)1/3].   Fig. 6. Comparison of the relation between crack length and cube root of compliance.  3.2.2. Mode Ⅰ fracture energy (separation energy)  The reloading was subjected as well as the case of the previous study.36 The bilinear traction-separation law 27, 36 was used by substituting the measured maximum normal traction σnt (=2.8 MPa) and quoting SEⅠC (=1.5 kJ/m2) of the previous study (See. Table 1).  Table 1 Material properties for FEM analysis Property Value Isotropic model EL [GPa] 27 ET [GPa] 27 νLT [-] 0.3 GLT [GPa] 10.4 Anisotropic model EL [GPa] 27 ET [GPa] 1.1 νLT [-] 0.3 GXY, YZ [GPa] 0.2 GXZ [GPa] 0.42 SEⅠC [kJ/m2] 1.5 σnt [MPa] 2.8 Fig. 7(a) shows the load-displacement curves of CZM. Table 2 also presents specific loads (PNL, pre, PNL, re, and Pmax) obtained from the curve. Under the preloading condition, the linear relation was confirmed at the initial stage of the loading. Subsequently, the nonlinear behavior started at PNL, Pre of each notch length specimen. After the maximum load Pmax, the load decreased with increasing displacement, and the relation was linear during unloading. During reloading, a linear load-displacement relation was observed until reaching a certain load (PNL, re in Table 2) which varied for different notch lengths. Above these loads the curves exhibited nonlinear behavior again.  To elucidate the reason of the difference of the nonlinear start load, PNL between preloading and reloading, the three-dimensional crack shape formation was investigated. In preloading, the crack shape showed a straight-line because of the initial geometry until the load reached PNL. However, in reloading, a curved profile around the crack tip was found before the crack started propagating again at PNL (Fig. 7(b), 5 mm notch length). Generally, plane strain field is dominant at the center of the rod. This can cause higher stress distribution near the center comparing to the surface. Consequently, the debonding seems to start from the center of the crack tip in preloading before the stable crack front shape is formed, and to cause the lower value of PNL in preloading.    Fig. 7. Load-displacement relation and crack propagation behavior estimated from CZM. (a) load-displacement relation. Small figure shows a definition of PNL,re and Parrst with experimental results 36 (b) normal tensile stress distribution of reloading [MPa] (5 mm notch length specimen)       SEⅠC can be obtained from the following equation.36  𝑆𝐸45 = 𝑃%!!6$72𝐷73𝐴02𝑎7 + 6𝐴07𝐴3𝑎 + 3𝐴0𝐴37: (4) Parrst is an initial unload after the cross head holding at each cycle obtained in our previous work. During the holding period, cracks propagate and cause a decrease of the load. This is why Parrst was employed for obtaining the SEⅠC instead of the Pmax. However, analysis of data from our previous report showed that PNL (or PNL, re) was nearly equivalent to the Parrst, within 10% standard deviation. Thus, although there is no theoretical justification for this at this time, for the convenience of the numerical analysis PNL, re was used to calculate SEⅠC by substituting Parrst = PNL, re in Eq. (4). The initial crack length a of each specimen in reloading was calculated from the Eq. (3) by substituting the calculated cube root of compliance, and obtained the values in Table 2.    Table 2 Parameters used for finite element analysis Notch length [mm] PNL, pre [N] PNL, re [N] Pmax [N] a , Eq. (4) [mm] SEⅠC [kJ/m2] 5 11 37 40 10 1.3 10 7.8 27 31 17 1.3 15 6.4 23 26 23 1.5 20 5.6 19 22 29 1.3  The SEⅠC obtained by FEM showed 1.4±0.09 kJ/m2. This value is apparently close to that measured in the previous study (1.5±0.5 kJ/m2).36 Therefore, we concluded that the proposed testing method can provide the effective CZM parameters for FEA. The accordance of SEⅠC also suggests the effectiveness of the modified double cantilever beam test proposed in the previous study. In this study, we substituted σnt at room temperature and fixed the crosshead speed. However, the rheological effect of polymer matrix can cause the variation of σnt. Therefore, further investigations on the dependency of CZM parameters under these test conditions are expected. In addition, this study conducted linear analysis. However, PP matrix has plasticity. Therefore, material nonlinearity effect on the evaluation of the parameters also needs to be elucidated.   4. Conclusion This study measured the mode Ⅰ-governed CZM parameter, the maximum traction, of a rod shape basalt fiber reinforced polypropylene composite by transverse tensile testing, and verified it by the three-dimensional FEM analysis. The analysis applying the measured CZM parameter apparently reproduced the relation between the crack length and cube root of compliance, and the SEⅠC to each of the previous report. Therefore, we conclude that the proposed evaluation method of CZM parameter for the composite rod is effective. The accordance of the SEⅠC also suggests the effectiveness of the modified double cantilever beam test proposed in the previous study.  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.   Acknowledgements  This research was promoted by COI program "Construction of next-generation infrastructure using innovative materials ~Realization of safe and secure society that can coexist with the Earth for centuries~ supported by Japan Science and Technology Agency (JST) Grant Number JPMJCE1315.   ORCID Yuta Tobata https://orcid.org/0000-0003-0730-5934  Reference [1] V. Dhand, G. Mittal , K.Y. Rhee, S. Park, D. Hui, Compos Part B, 73, 2015, 166-180. [2] D.K. Rajak, D.D. Pagar, P.L. Menezes, E. Linul, Polymers. 11(10), 2019, 1667.  [3] H. Zhou, B. Jia, H. Huang, Y. Mou, Mater. 13(6), 2020, 1362.  [4] O.A. Mohamed, W.A. Hawat, M. Keshawarz, Polymers. 13(9), 2021, 1402.  [5] E. Quagliarini, F. Monni, S. Lenci, F. Bondioli, Constr Build Mater. 34, 2012, 372-380 [6] Onur Aykanat, Mahmut Ali Ermeydan, Polym. Compos, 43(6), 2022, 3519-3530, [7] E. Kessler, R. Gadow, J. Straub, AIMS Mater Sci. 3(4), 2016, 1561-1576 [8] Sharma, V, Meena, ML, Kumar, M, Patnaik, Polym. Compos. 41, 2020, 3717– 3731.  [9] G, J., Ranganathan, N., Kader M, A., Polym. Compos. 43(6), 2022, 3609.  [10] C. Lu, M. Ni, T. Chu, L. He, J. Mater in Civil Eng. 32(7), 2020, 04020170 [11] Lu, Z., Jiang, M., Pan, Y., Xian, G., Yang, M., Polym. Compos. 43(4), 2022, 1961.  [12] Q.T.H. Shubhra, A. Alam, M. Quaiyyum, J. Thermoplast. Compos. Mater. 26(3), 2014, 362-391. [13] H.A. Maddah, Am J Polymer Sci. 6(1), 2016, 1-11. [14] Liang, S., Liu, J., Guo, Y., Luo, J., Liu, H., Peng, S., Polym. Eng. Sci., 62(10), 2022, 3168.  [15] Poornima, C., Uthamballi Shivanna, M., Sathyanarayana, S., Polym. Compos. 2022, 1.  [16] K. Uzawa, Y. Saito, A. Hokura, JSCE, Ser. A1 (Structural Engineering & Earthquake Engineering (SE/EE)) 73(5), 2017. [17] J. Tanks, K. Naito, H. Ueda, Polymers. 13, 2021, 3136. [18] D.P. Dias, C. Thaumaturgo, Cem Concr Compos, 27(1), 2015, 49-54. [19] Y. Zhao, X. Sun, P. Cao, Y. Ling, Z. Gao, Q. Zhan, X. Zhou, M. Diao, Hindawi Advances in Civil Engineering. 2019, 5630805. [20] M.K. Kazberuk, J. Krassowska, MATEC Web of Conferences. 265, 2019, 01008. [21] F. Shi, T.M. Pharm, H. Hao, Y. Hao, Constr Build Mater. 30, 2020, 120108. [22] N. Liang, L. Ren, S. Tian, X. Liu, Z. Zhong, Z. Deng, R. Yan, Int J Concr Struct Mater. 15, 2021, 35. [23] W. Li, J. Xu, Mater. Sci. Eng. A 505, 2009, 178-186. [24] G.I. Barenblatt, Adv. Appl. Mech. 17(1), 1962, 55-129. [25] D.S. Dugdale, J Mech Phys Solids 8(2), 1960, 100-104. [26] A. Needleman, J. appl. Mech. 54(3), 1987, 525-531 [27] G. Alfano, M.A. Crisfield, Int J Numer Methods Eng. 50, 2001, 1701-1736. [28] R.D. Borst, Eng Fract Mech. 70(14), 2003, 1743-1757. [29] C. González, J.L. Lorca, Compos Sci Tech 67, 2007, 2795-2806. [30] J.G. Williams, J. Strain Anal Eng Des 24(4), 1989, 207–14.  [31] L.A. Carlsson, D.F. Adams, R.B. Pipes, Characterization of delamination failure. In: Experimental characterization of advanced composite materials. 4th edition. New York: CRC Press. 2014, p. 250–270.  [32] C. Colombo, L. Vergani, M. Burman, Compos Struct. 94(3), 2012, 1165–74.  [33] T. Chen, C.M. Harvey, S. Wang, V.V. Silberschmidt, Int J Fract. 229(2), 2021, 269–76.  [34] R. Jones, A.J. Kinloch, J.G. Michopoulos, A.J. Brunner, N. Phan, Chapter thirteen – delamination growth in polymer – matrix fiber composites and the use of fracture mechanics data for material characterization and life predictions, aircraft sustainment and repair. Oxford: Butterworth-Heinemann, 2018, 763–97 [35] M. Sabaghi, F. Taheri-Behrooz, M. Salamat-Talab, Polym. Compos. 43(9), 2022, 6135.  [36] Y. Tobata, K. Naito, J. Tanks, Compos. Struct. 279, 2022, 114778. [37] Material model – Bilinear Behavior, 4.12.1.2., ANSYS, Mechanical APDL Theory Reference release 18.2 [38] Y. Ogawa, K. Naito, K. Harada, H. Oguma, Int J Adhes Adhes. 117, 2022, 103172