Microsoft Word - 2180 M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 277 Fracture and Structural Integrity: ten years of ‘Frattura ed Integrità Strutturale’ A numerical model based on ALE formulation to predict crack propagation in sandwich structures Marco Francesco Funari, Fabrizio Greco, Paolo Lonetti University of Calabria, Italy marcofrancesco.funari@unical.it; https://orcid.org/0000-0001-9928-3036 fabrizio.greco@unical.it, https://orcid.org/0000-0001-9423-4964 paolo.lonetti@unical.it, https://orcid.org/0000-0003-0678-6860 Saverio Spadea University of Dundee, UK s.spadea@dundee.ac.uk, https://orcid.org/0000-0003-3525-5217 ABSTRACT. A numerical model to predict crack propagation phenomena in sandwich structures is proposed. The model incorporates shear deformable beams to simulate high performance external skins and a 2D elastic domain to model the internal core. Crack propagation is predicted in both core and external skin-to-core interfaces by means of a numerical strategy based on an Arbitrary Lagrangian–Eulerian (ALE) formulation. Debonding phenomena are simulated by weak based connections, in which moving interfacial elements with damage constitutive laws are able to reproduce the crack evolution. Crack growth in the core is analyzed through a moving mesh approach, where a proper fracture criterion and mesh refitting procedure are introduced to predict crack tip front direction and displacement. The moving mesh technique, combined with a multilayer formulation, ensures a significant reduction of the computational costs. The accuracy of the proposed approach is verified through comparisons with experimental and numerical results. Simulations in a dynamic framework are developed to identify the influence of inertial effects on debonding phenomena arising when different core typologies are employed. Crack propagation in the core of sandwich structures is also analyzed on the basis of fracture parameters experimentally determined on commercially available foams. KEYWORDS. Moving Mesh Method; Crack Propagation; Sandwich Structures; ALE; Finite Element Method; Debonding Mechanisms. Citation: Funari, MF., Greco, F., Lonetti, P., Spadea, S., A numerical model based on ALE formulation to predict crack propagation in sandwich structures, Frattura ed Integrità Strutturale, 47 (2019) 277-293. Received: 25.08.2018 Accepted: 05.10.2018 Published: 01.01.2019 Copyright: © 2019 This is an open access article under the terms of the CC-BY 4.0, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited. http://www.gruppofrattura.it/VA/47/2180.mp4 M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 278 INTRODUCTION omposites materials are widely utilized in several applications ranging from aerospace to civil engineering fields [1, 2]. Sandwich structures are a particular class of composites consisting of two thin face sheets made of stiff and strong materials such as metal or fiber reinforced composites bonded to a thick and deformable core with low density [3]. They are able to ensure a good resistance under bending/shear loading, offering a great variety of lightweight structural systems. Unfortunately, sandwich panels are affected by both macroscopic and microscopic damage phenomena, mainly produced by the heterogeneity of the layered systems, which reduce the integrity of the composite structure, leading to catastrophic failure mechanisms [4]. From physical and mathematical viewpoints, two main issues are demanding a detailed understanding of the mechanical behavior of sandwich panels: the propagation of internal macro-cracks in the core [5] and the delamination at face/core interfaces [6]. These problems have been addressed by mean of different numerical approaches, mostly developed in the framework of the Finite Element Method (FEM) due to its versatility to model complex structures. Specific modelling techniques are required to predict crack tip motion of internal material discontinuities. Interface elements based on Cohesive Zone Model (CZM) or Linear Elastic Fracture Mechanics (LEFM) are frequently used to predict crack tip evolution. Discrete or distributed interface elements can be easily incorporated into FEMs, by introducing constitutive traction forces between adherent internal surfaces [7]. These methodologies are frequently used in sandwich structures to predict the crack evolution at the core/skin interfaces, since the crack motion is expressed as a function of a linear positional variable coinciding typically with the interface coordinate. However, such modeling is affected by numerical problems due to mesh dependence, computing inefficiency, and sensitivity to the element aspect ratio. These issues may be partially addressed by adopting a very fine discretization at the crack tip front, but numerical complexity remains, due to the high number of computational points requested. Sandwich structures are also affected also by macro-cracks in the core. Quite complex scenarios are observed in presence of kinking phenomena of the crack, starting from the interfaces. In these cases, the crack growth requires more advanced numerical modeling techniques, since it needs to be expressed both in terms of angle of propagation and tip displacement. The use of CZM is quite cumbersome, since a very large number of interface elements need to be introduced, at least along the path where the crack growth is expected. Alternatively, crack propagation in 2D continuum elements can be achieved by using adaptive mesh refinement methods, in which the element boundaries are coincident with moving internal discontinuities [8]. However, an accurate description in terms of field variable interpolation is needed to describe the updated set of nodes at each mesh adaptation step and the corresponding quadrature points. Moreover, computational errors are introduced due to the projection procedures required by the re-meshing process [9]. Another possibility is to tackle the problem through the Boundary Element Method (BEM). In this case, only the structures boundaries (and not the internal domains) are represented by means of a proper mesh discretization [10]. Although such hypothesis reduces the computational costs required to generate new elements, computational complexity due to the need to define singular integrals remains. Previous formulations are classified in the literature as geometrical representation approaches [11], since an explicit definition of the cracked surfaces is required by the numerical models in order to evaluate the fracture variables and the subsequent crack propagation. Among the formulations in which an implicit crack definition is achieved, Extended Finite Element Method (XFEM) is currently used with success in many practical applications. The basic idea is to use nonconforming elements to model macro-cracks by enriching shape functions of the mesh elements by discontinuity properties. However, a further extension is required to predict fracture variables for nonlinear problems, especially in presence of frictional effects [12]. Moreover, the methodology needs a different number of kinematic variables for each node and thus the total number of mesh points may vary with the crack growth. Others methodologies based on Discrete Element Method (DEM) [13] or MeshFree Methods (MFMs) [14] have been formulated in the last decade, providing valid alternatives to study such problems. Methods based on Moving Mesh technique (MM) provide a feasible and sensible way to predict crack growth mechanisms in continuum media. Early studies were developed in [15], where MM was employed to predict energy release rate by using a virtual crack extension. The Arbitrary Lagrangian-Eulerian Formulation (ALE) was only recently implemented in Fracture Mechanics in [16, 17]. A generalization of the ALE approach was also proposed in the framework of weak based moving cohesive forces, where the interlaminar debonding phenomena are predicted without modifying the formulation of the structural problem [18]. The extension to a generalized crack path is quite rare in the literature, especially in those cases in which a proper prediction of crack angle variability is required. To the authors’ knowledge, only in [19] a generalized mesh refitting procedure applied to a continuum media is developed. The main goal of this paper is to generalize the numerical implementation proposed in [20-22] in the framework of sandwich structures, with the purpose to describe delamination phenomena along the interfaces and macrocracks evolution in the C M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 279 core. The proposed strategy explores the possibility to combine different crack growth phenomena. In detail, interface debonding phenomena at core/skin interfaces are simulated by using a weak-based moving interface strategy. Moreover, the evolution of internal cracks in the core is predicted by using a strong-based ALE strategy, in which the governing equations are expressed in the moving ALE coordinates by means of transformation rules between moving and fixed referential configurations and proper crack growth functions to identify the crack tip motion. The model is implemented numerically by using a finite element approximation. This allows the development of a series of numerical results demonstrating the effectiveness of the method to simulate both the interfacial debonding and the crack propagation in the sandwich core. Experimental tests are also performed in order to analyze the fracture parameter of a commercially available foam that are commonly used as core material in real applications. The actual material properties experimentally obtained are employed to inform the numerical simulations. THEORETICAL FORMULATION he proposed model is formulated in the framework of a two-dimensional idealization of a sandwich structure. It consists of an internal core, modelled by means of a plane stress formulation, and two external skins, following a Timoshenko beam kinematic. The formulation is able to predict crack growth of material discontinuities, which may affect the skin/core interfaces and the core. At interface level, this is achieved by the use of moving interface elements, which ensure an accurate description of the fracture variables in the process zone in terms cohesive Traction Separation Laws (TSL). For the core region, a generalization of the above interface model is proposed, since the crack tip may evolve in the two-dimensional domain. A synoptic representation of the model is reported in Fig. 1, where, without loss of generality, an initial crack length is assumed in both interface and core. Figure 1: Schematic representation of the sandwich structure: interfacial and core macrocracks. Interface/core debonding phenomena ALE strategy is implemented in the interface regions to accurately describe the evolution of debonding phenomena. In particular, moving cohesive interface elements simulate the traction forces produced by the evolution of material discontinuities. Cohesive constitutive laws are parametrized in terms of a moving interface coordinate system in order to simulate the motion of the process zone acting at the skin/core interfaces. From the mathematical point of view, the parametrization of positional variables is expressed as a function of two configurations, i.e. Referential Configuration (RC) or Moving Configuration (MC). The relationship between RC and MC coordinates of each particle is descripted by the mapping operator  (Fig.2), as follows:  , :X t with RC MC    (1) where  and X are the referential and moving coordinates, respectively. The mesh motion is expressed, by introducing the regularization or rezoning equations, defined in terms of the mesh displacement function X of the computational nodes. T M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 280 Without loss of generality, the following Laplace-based equations developed either in Statics (S) or Dynamics (D) are considered [23]:        2 3 , 1 , 12 2, , , , 0 ( ) 0 ( ) j j j jt t X X t S X X t D t                                 (2) where j, with j=(L, R), represents the index referred to the Left (L) or Right (R) process zones, whereas jX corresponds to the mesh displacement function. Previous equations should be completed by boundary and initial conditions. Homogeneous conditions are required to prescribe mesh displacements at left and right ends, in such a way to verify the coincidence in terms of positional variables:    0 0, 0 0 X X       0, 0 X L X L    (3) Figure 2: Moving and referential coordinate systems for the skin/core interfaces. In addition, the prescribed mesh motion at debonding length is defined by introducing a crack growth function. During the crack tip advance, the crack function is always equal to zero, whereas it is negative when the crack tip is stationary. A classical crack growth criterion, defined in terms of the Energy Release Rate (ERR) mode components normalized on the critical value, is assumed:      2 2 1 r r I II f IC IIC G X G X g X G G                   with   0fg X  (4) M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 281 where r is the parameter used to describe fracture in different material and  ,IC IICG G are the total area under the traction separation law. Moreover,  ,I IIG G are the individual ERRs, which correspond to the integral function of the TSL defined in terms of the maximum opening/axial stresses and displacements. For the sake of simplicity, classical bilinear relationships of the TSL are assumed in the present paper, whilst the generalization to more complex cases can be implemented straightforward. Boundary and initial conditions at the crack tip front are verified by enforcing a rigid displacement of the process zone on a specific length, namely , at the extremities of the computational nodes. This choice ensures that the NL involved in the debonding mechanisms are constrained to a small portion containing the process zone, reducing the total complexities of the model. Consequently, the following boundary conditions should be considered in the analysis:    1 j jj jX X X X    with   0 jj fg X  and   0 j fg X   (5) where + or – refer to the Right or Left debonding crack fronts, respectively, and  , j j X X  correspond to the coordinates of the extremities of the debonding region (Fig.2). It is worth noting that the conditions in Eqs.(5) are prescribed by means of a simple procedure, which consists, at first, to predict the values of the fracture function at the extremities of the debonding region and, subsequently, to enforce that each step of the crack growth corresponds to a null value of the fracture energy. Therefore, by using a linear approximation function along the debonding region, the current nominal crack tip displacements can be expressed by means of the following relationships:       0 0, 0 j fj jj j f fj jj j f f g X toll X g X g g X toll g X            (6) Figure 3: Moving and referential coordinate systems for core domain. Core debonding crack growth The evolution of preexisting cracks in the core is simulated by the generalization of the formulation developed in previous subsection to a two-dimensional domain. Two configurations are introduced to describe the mesh motion defined as referential or material ones. The latter is modified by the geometry variations produced by the crack advance, whereas the M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 282 former is basically fixed or at least re-meshed in those cases in which large distortions occur. The mapping between the two configurations is governed by the following expressions:  M RX X      1 R MX X     (7) where : R MC C   with  X Y     is assumed to be invertible with continuous inverse,  ,M RX X   are material and referential positional vector functions with  T M M MX X Y  and  T R R RX X Y  . In addition, the transformation between the RC and the MC is also expressed by gradient operator defined as a function of the Jacobian matrix, on the basis of Eq.(7):    1 M RJ       X Y R R X Y R R X Y J X Y                (8) According to the ALE methodology, rezoning or regularization equations are required to modify the mesh position, reducing mesh distortions. In addition, on the external boundaries, boundary and initial conditions should be enforced to constraint the mesh motion to do not exceed the structural domain. Therefore, according to a mesh regularization method based on Laplace’s method, the following relationships hold: 2 20 0 0, X Y UX on S             (9) where 2  is the nabla operator and SU is the external boundary of the core system. Additional equations are required to modify the mesh positions of moving boundaries, which describe preexisting internal cracks. However, in order to quantify the crack advance in terms of displacement and angle of propagation, a crack growth criterion should be introduced. Similarly to what proposed to simulate the core-to-skin interfaces debonding, a fracture function based on local-Griffith approach criterion defined in terms of the ERR is assumed. Moreover, the Maximum Energy release rate criterion is used [24] for the crack angle prediction. In particular, the crack tip boundary conditions are defined in terms of incremental displacements by means of the following equations:     2 1 cos , sin , with tan and T T T T J X Y arc J                       (10) where  is the angle of crack propagation  ,T TX Y   are the crack tip components along (XR, YR) coordinate system, (J1, J2) are the components of the J integral along normal (e1) and tangential (e2) directions. Moreover, T is the incremental displacement of the tip, which is determined by solving the following constrained optimization problem: 0, 0 0T F T F T Ff f f        (11) NUMERICAL IMPLEMENTATION he model is implemented numerically by using a finite element approximation. In particular, the structural formulation of the skins follows a Timoshenko beam model, whereas the core is based on a 2D plane-stress description. Skins and core are connected through internal forces arising from the TSL of the moving interface elements. It is worth noting that the governing equations of the skins are not affected by the ALE interface formulation, since it contributes as weak contribution in the structural formulation. Contrarily, in the case of the core region, governing equations should be modified by changing the variables and their derivatives from the MC to the RC, by using T M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 283 transformation rules defined by Eq.s (7)-(8). In the present paper, more emphasis is given to the numeric procedure implemented to simulate the crack growth in the core, since the one related to the core/skin interfaces can be recovered in previous author’s works [21, 22, 25]. From the numerical point of view, the crack growth is achieved by evaluating the interfacial variables on a small region containing the crack tip. In the proposed model, the J-integral methodology is used to evaluate on a contour enclosing the tip, in which a high number of elements is introduced to ensure accuracy in the prediction of the fracture variables. During the crack growth, such region is moved rigidly, avoiding mesh distortions in the crack tip region. Moreover, the boundaries affected by the crack advance are enforced to have the same displacements of the structure, leaving unaffected the shape of the cracked faces. This task is achieved by introducing two nodes close to the region adjoining the crack tip, which are stretched as far as the angle variation predicted by Eq.(10). 3 is lower than a fixed angle tolerance value. Once tolerance condition is satisfied, a new definition of the computational nodes, driven by the new value of crack propagation angle, is required. During the crack growth, mesh movements of the computational nodes produce distortions in the grid points, which are eliminated by the use of a re-meshing algorithm. This is able to reconstruct a new regular mesh discretization, transferring the nodal variables from the distorted to the new computational points. The procedure is recalled by means of a mesh quality parameter, which controls the allowable distortion in each element. Governing equations as well as the steps referred to above are formulated by using a customized FE subroutine in the framework of COMSOL Multiphysics software [26]. The algorithm was developed by means of script files implemented in a MATLAB® environment, which manage the parameters and the results required by the iterative procedure. In particular, the following steps are implemented in the user subroutine: 0) Read the input data: geometry, material, mesh discretization and load configuration; 1) Evaluate fracture variables and crack angle of propagation (Eq.(10).3); 2) Check crack growth conditions (Eq.(11)); 3) Solve incremental Structural and ALE problem; 4) Check tolerance conditions for the angle variation or mesh quality; 5) If the angle variation or mesh quality tolerance are satisfied, proceed to the crack growth, else employ re-meshing algorithm. The proposed procedure is quite general and can be solved either in a static or dynamic frameworks, taking into account the time dependent effects produced by the inertial characteristics of the structure and the boundary motion involved by debonding phenomena. Since the governing equations are essentially nonlinear, an incremental-iterative procedure has been adopted to evaluate the current solution. FRACTURE PROPERTIES OF FOAMS n this section, the fracture toughness of a commercially available core foam under mode I is evaluated by three point bending tests on Semi-Circular Bending (SCB) specimens. Mode II is analysed on Asymmetric Semi-Circular Bend (ASCB) specimens. Ayatollahi et al. [27] proved that the ASCB loading scheme is able to generate all range of mixed fracture modes in fragile construction materials, including a pure mode II fracture condition. The method was used by Marsavina et al. [28, 29] for PVC foams. The experimental determination of fracture properties of the core is instrumental to the implementation of the numerical simulations that follow in the paper. Experimental Setup The experimental tests were performed using a Tinius Olsen testing machine equipped with a 5 kN load cell and an appropriate bending rig. All specimens were cut from 20 mm Divinycell H100 panels in the two main directions using a Denford CNC router with a 0.1 mm resolution equipped with a 3 mm drill bit. Although the numerical cutting procedure included the mid-span sample notching, the natural crack was initiated by sliding a fresh razor across the notch root. The samples were uniformly sprayed with several coats of white paint. Subsequently a stream of black paint was applied as to obtain a random b/w speckle pattern, the ideal reference system for Digital Image Correlation (DIC) purpose [30]. Fig. 4 shows the two load conditions taken in consideration. Three samples were tested according to each scheme, accounting for a total of six specimens. All samples had identical geometrical properties: the radius, R, of the semi-circular samples was 80 mm, whereas the length of the notch, a, was 40 mm. I M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 284 The SCB loading scheme consisted of a point load centered on the top edge of the half-disc, which was spanned at 120 mm  1 2 60 mmS S  . In the ASCB loading scheme the left supports was moved to a 5.5 mm  2S distance from the specimen mid span, with the second supports remaining at a 60 mm  1S distance from midspan. All tests were performed in displacement control, at 2 mm/min rate. Load and displacement data were continuously recorded using the testing machine built-in sensors and data acquisition system. Additionally, fixed focal lens camera shots were taken every 5 second to digitally monitor the specimen change of configuration under the applied loads. (a) (b) Figure 4: Experimental Setup: (a) SCB specimen; (b) ASCB specimen. Results Fig 5(a) and (b) show the average force vs. mid-span vertical displacement relations, exhibited by the SCB and ASCB specimens, respectively, as well as the xx strain maps retrieved by means of the DIC processing. The black curves refer to the specimen displacement measured at the top point A, using the machine built-in transducer. These values are affected by the loading-pin penetration, as evidenced by the heavy non-linear behavior of the curve. An additional measure of the experimental displacement at the crack mouth (point B) was retrieved from the correlation of the digital images collected, each corresponding to a known load step. The DIC processing was performed using the freeware Matlab script Ncorr, developed by Blaber et al. [31]. The red curve obtained in this manner is a more realistic representation of the flexural deflection of the specimen, being only affected by minor local compressions occurring at the specimen supporting pins. The average values of maximum forces detected in the experimental tests, as detailed in Tab. 1, were employed for the calculation of the Stress Intensity Factors (SIFs). The SIFs for mode I, 1cK , and mode II 2cK , were computed according to the following expression, proposed by Ayatollah et al. [27]:  max 1 2, , , 1, 2 2ic i i P K a Y a R S R S R i Rt   (12) where maxP is the average value of the maximum load experimentally measured either on SCB or ASCB specimens, and the functions  1 2, ,iY a R S R S R are given by Marsavina et al. [28, 29]:                     3 2 1 2 2 2 2 5 4 3 2 2 2 2 2 2 2 2 6.235 15.069 17.229 1.062 1.884 7.309 5.037 2.77 5.075 1.983 Y S R S R S R S R Y S R S R S R S R S R S R           (13) Computed values of SIFs with reference to the two set of tests are listed in Tab. 1. M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 285 (a) (b) Figure 5: Load – midspan vertical displacement relationship and DIC strain maps: (a) SCB specimen; (b) ASCB specimen. Divinycell H100 1cK [MPa m0.5] 2cK [MPa m0.5] 1 maxP [N] 2 maxP [N] 0.222±0.001 0.116±0.002 334±1 628±5 Table 1: Experimental results in terms SIFs and maximum load reached. SIMULATION OF THE INTERFACIAL DEBONDING n this section, the capability of the proposed model to predict the debonding failure mechanism in sandwich panels at the skin/core interface is verified by means of comparisons with numerical and experimental results available in literature. Structures featuring aluminum or glass/polyester face-sheets and different type of foam cores are taken in consideration. In detail, the static behavior of interfacial crack propagation at the upper interface between face-sheet and core is investigated on classical Double Cantilever Beam (DCB) and Mixed Mode Bending (MMB) loading schemes. Subsequently, the dynamic behavior is analyzed with the aim to identify the influence of inertial effects produced by different levels of loading rate and by different core typologies. Double Cantilever Beam As a first case study, the DCB loading scheme studied by Odessa et al. [32] is investigated. Fig. 6(a) illustrates the geometry of the specimen, as well as loading and boundary conditions, whereas mechanical properties assumed for the skins, core and interfaces, are summarized in Tab.2., according to data reported in [32]. Fig. 6(b) shows the dimensionless values of the opening force  cF G B and nominal crack tip position  TX L in function of the opening displacement, with reference to two different core thickness configurations, i.e. 15 20 mmch   . The opening force values obtained are compared to both experimental [33] and numerical [32] results retrieved from literature, showing excellent agreement and demonstrating the capability of the model to efficiently simulate a mode I debonding phenomenon. It is also worth noting how an increment in the core thickness does not produce significant variations in the loading curve, due to the high deformability of the foam compared to the skin layers. I M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 286 (a) (b) Figure 6: DCB specimen: (a) loading scheme; (b) opening force and nominal crack tip position in function of the opening displacement, comparison to experimental [33] and numerical data [32]. Face-sheet (aluminium) s 11E [MPa] s 12G [MPa] s s [kg/mc] 70 103 26 103 0.33 2700 Core-PMI AIRES R 90.400 c 11E [MPa] c 12G [MPa] c c [kg/mc] 420 220 0.25 400 Interface properties cG [N mm-1] 0 [mm] 0.550 0.12 Table 2: DCB test: material and interface properties (a) (b) Figure 7: MMB specimen: (a) loading scheme; (b) opening force and nominal crack tip position in function of the opening displacement, comparison to experimental data [3]. M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 287 Mixed Mode Bending In this subsection, the MMB loading scheme shown in Fig. 7(a), experimentally analyzed by Carlosson and Kardomateas [3] and Quispitupa et al. [34], is taken into consideration. The values of material properties of the medium-density foam Divinycell H100, as reported by the referred authors and assumed in current analysis, are reported in Tab. 2. A preliminary step of the analysis involves the calibration of the cohesive model parameters, which is based on the experimental force / opening displacement relation observed in [3, 34] for two different values of the lever arm c (i.e. 40, 50 mm). Fig. 7(b) shows the best fitting curves obtained, both corresponding to a value of the separation energy for unit area, cG , and characteristic slip parameter, 0 , equal to 0.800 N/mm and 0.12 mm, respectively. Based on such interface properties, the investigation is therefore extended to a dynamic context, with the aim to investigate the influence of the loading rate and the inertial effects produced by different typologies of core. The displacement rate law adopted in the analysis is a linear ramp curve until the instant t0, with t0 assumed half of first period of vibration (T1) of the structure, followed by a constant rate value v0. A parametric study is performed according to the following values of v0: 1) v0= 1 m s-1; 2) v0= 5 m s-1; 3) v0= 10 m s-1. Face-sheet glass/polyester s 11E [MPa] s 12G [MPa] s s [kg/mc] 16.4 103 2.7 103 0.17 1500 Core-Divinycell H100 c 11E [MPa] c 12G [MPa] c c [kg/mc] 135 35 0.32 100 Interface properties cG [N mm-1] 0 [mm] 0.800 0.12 Table 3: MMB test: mechanical and interface properties. c 11E [MPa] c 12G [MPa] c [kg/mc] Divinycell H35 40 12 38 Divinycell H100 135 35 100 Divinycell H250 400 97 250 Table 4: Mechanical properties DIVINYCELL H. (a) (b) Figure 8: MMB test: (a) influence of the loading rate on the loading-displacement curve; (b) influence of the loading rate on the nominal crack tip speed. M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 288 (a) (b) Figure 9: MMB test: (a) influence of the core typology in terms on the loading-displacement curve; (b) Influence of the core typology on the nominal crack tip speed. Fig. 8(a) depicts the first set of results obtained in terms of dimensionless resistance curves. Whereas static and dynamic solutions are substantially overlapped for low value of v0, a remarkable increment of the peak load and a slightly unstable post-peak behavior can be observed when the value of v0 is increased. Furthermore, the relationship between nominal crack tip speed and applied displacement is shown in Fig. 8(b). The crack speed tends to be larger in the initiation phase while subsequently tends to decrease according to a damped oscillation, during the process of delamination. The second set of results pertains the influence of the core mechanical properties (on the interfacial debonding under the dynamic loading rate described above and adopting a constant value of v0 equal to 10 ms-1. As shown in Tab. 4, the previously analyzed core (Divinycell H100) is used as reference material in comparison to a high- density core (Divinycell H250) and a low-density one (Divinycell H35). Fig. 9(a) shows that peak loads observed in different cores are generally similar, whereas marked differences in terms of initial stiffness are visible. In particular, the higher density specimen (H250) produces a slight increment in terms of stiffness and almost same peak load as the medium-density one. On the other hand, the lower density core (H35) can guarantee relevant deformation capacity during the crack propagation, whilst its performance in terms of stiffness and peak load are less remarkable. Finally, Fig. 9(b) presents an investigation on the nominal crack tip speed in function of the applied displacement. The results show that an increment in the core density does not produce relevant amplifications in the nominal crack speed. It is also evident that more deformable cores are associated with a delayed crack onset, and consequently a higher fracture toughness. CRACK PROPAGATION IN A SANDWICH STRUCTURE CORE n this set of numerical results aims to validate the capability of model to describe the crack propagation in 2D solids. With this purpose, a sandwich structure fixed at bottom skin with an initial horizontal crack is taken into consideration as reference case study (Fig. 10). Although the proposed model is able to describe both the skin/core interface delamination and the crack propagation in the core, the aim of these simulations is specifically to describe the crack path occurred in the core of the sandwich structures and verify its effects on the interfacial traction forces. Geometry, load conditions, and material properties The reference specimen has a 150 mm length and width equal to 35 mm. The core and skin thickness are 75 mm and 2 mm, respectively. The initial crack has a length to 25 mm and it is located at the core mid-depth. Two different loading conditions are investigated: I M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 289 1) a uniformly distributed opening force applied along the upper skin of the panel (Fig. 10(a)). 2) a point opening force applied at the cracked edge of the panel (Fig. 10(b)); The idealized sandwich structure presents glass/polyester skins and a Divinycell H100 foam core, therefore the material and interface properties are the same as the ones listed in Tab. 3 with reference to the MMB test. The fracture parameters 1CK and 2CK experimentally determined Divinycell H100 foams, as detailed in the experimental section of this paper, are employed in this analysis. (a) (b) (c) Figure 10: (a) uniformly distributed opening force applied along the upper skin of the panel; (b) point opening force applied at the cracked edge of the panel (b); detail of the mesh discretization in the core (c). Modelling The 2D model involves a discretization of each face sheet by means of shear deformable beam elements, whereas the foam core is modelled by plane stress triangular elements. Both skin and core elements features cubic interpolation functions and a linear elastic constitutive behavior of the materials. ALE equations are introduced at each skin/core interfaces, in view of possible debonding phenomena, and in the core, to predict the propagation of the pre-existing crack length. Linear interpolation functions are assumed for the ALE variables. The numerical discretization used for the face sheet is assumed to be generally uniform with an element length D equal to 5mm . A coarse discretization is adopted at the interface, in absence of any pre-existing interfacial defects. With reference to the core modelling, a high number of computational points is condensed around the crack tip region, whereas a coarse discretization is adopted in the remaining domain (Fig. 10(c)). The analysis is developed in displacement control mode, with the aim to ensure a stable crack propagation. (a) (b) Figure 11: Load-displacement response: (a) uniformly distributed opening force; (b) point opening force. M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 290 Numerical Results Figs. 11 (a) and (b) show the load / opening displacement curves numerically obtained with reference to both the schemes here considered. As shown, the constitutive behavior shown is initially linear and stable. Once the crack function criterion is satisfied, the curves show a sudden descending trend. A clearly different trend can be observed in the two descending branches: the point loading causes a snap-back, a phenomenon that is not evidenced when a uniformly distributed opening load is applied to the structure. In Figs. 12 and 13, the Von Mises stress maps at four subsequent load steps are shown together with crack propagation pattern for both loading configurations. Although in both cases a mixed mode fracture arises, the crack propagation paths exhibited by the two specimens are visibly distinguishable. Whereas the point load activates a mode II dominated fracture process I, the distributed load tends to cause a mode I dominated fracture mode. In sandwich structures, loads are typically applied to the face-sheets and these are transmitted to the core by the adhesive interface. The core behavior and its sensitivity to crack propagation is, therefore, strongly affected by the interfacial stress, which are, in turn, related to the external loads, the materials stiffness’s, and the cohesive law. In this view, the interfacial traction forces detected at the upper skin-to-core interface are shown in Figs. 14 and 15 for both load configurations analyzed. Fig. 14(a) depicts the distribution of the interfacial forces corresponding to the peak load of the linear elastic branch (A) in the distributed load case. The distribution of stresses is quite uniform, which suggest the cohesive interface is not threatened by any possible debonding phenomenon. The propagation of the crack in the core is able to slightly affect the distribution of the interfacial traction forces (Figs. 14(b), (c), and (d)), whilst it does not substantially modify the quality of interface response. (a) (b) (c) (d) Figure 12: Uniformly distributed opening force: Von Mises contour plots and crack evolution for different loading steps. (a) (b) (c) (d) Figure 13: Point opening force: Von Mises contour plots and crack evolution for different loading steps. A similar analysis is illustrated in Fig. 15 (a-d) with reference to the point opening force scheme. As shown in Fig. 15(a), high values of the interfacial traction arises at the sandwich panel edge, in correspondence of the applied load. The onset condition could have been easily activated at the interface region but it is prevented, in this case, by the crack propagation in the core. This phenomenon tremendously affects the interfacial stresses, whose peak tends to move congruently with the crack tip position. However, a more even distribution of stresses suggests that an interface debonding is unluckily to occur at this stage. The results show how the interfacial stresses distribution can affect the crack propagation in the core and highlights the usefulness of a numerical model able to couple the two effects in an effective manner. M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 291 Figure 14: Uniformly distributed opening force: interfacial stresses across the upper cohesive interfaces at different value of core’s crack tip positions. Figure 15: point opening force: interfacial stresses across the upper cohesive interfaces at different value of core’s crack tip positions. M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 292 CONCLUSIONS he proposed model is developed with the purpose to study the behavior of sandwich structures affected by debonding phenomena and macrocrack in the core region. The numerical model is inspired by the previous works of the authors, performed in the framework of the layered structures and here generalized to sandwich structures. In this work, the numerical model is based on a 2D plane stress formulation to simulate the internal core, whereas the face- sheets follow a one dimensional model based on Timoshenko beam kinematics. In order to describe the delamination process, the proposed approach combines ALE formulation with a CZM. Moreover, crack growth in the core is predicted by a 2D moving mesh approach, in which a proper fracture criterion and mesh refitting procedure are introduced to predict crack tip front direction and displacement. Compared to existing formulations available in literature, this model presents lower computational complexities in the governing equations. In particular, the combination between CZM and ALE formulations, gives the possibility to introduce nonlinear interface elements in a small region containing the crack tip front, whereas in the remaining a coarse discretization is considered. In the model here proposed, macro-cracks at core/skin interface or in the core are decoupled during the crack evolution, therefore future developments are required to predict interaction phenomena, such as coalescence or initiation mechanisms. However, comparisons with numerical and experimental results, developed in both static and dynamic frameworks, show the capability of the proposed model to predict interfacial debonding phenomena at core/skin interfaces or the evolution of internal macrocracks in the core. ACKNOWLEDGEMENTS he authors acknowledge Frank McCauley and Jeremy Moore from Diab Limited UK for kindly providing the foam core material that was employed in the experimental section of this work. This work is also supported by Italian Ministry of University and Research (P.R.I.N. National Grant 2015, B86J16002300001) and SRPe (Scottish Research Partnership in Engineering) under the PECRE 2017/18 Scheme. REFERENCES [1] Spadea, S., Orr J., Ivanova K. (2017). Bend-strength of novel filament wound shear reinforcement. Composite Structures, 176, pp. 244-253. [2] Ascione, L., Berardi V. P., Giordano A., Spadea S. (2015). Pre-buckling imperfection sensitivity of pultruded FRP profiles. Composites Part B: Engineering, 72, pp. 206-212. [3] Carlsson, L. A., Kardomateas G. A. (2011). Structural and Failure Mechanics of Sandwich Composites. York SDHLN, Springer Dordrecht Heidelberg London New York, Springer. [4] Mortas, N., Reis P. N. B., Ferreira J. A. M. (2014). Impact response of balsa core sandwiches. Frattura ed Integrità Strutturale, 30, pp. 403-408. [5] Morada, G., Vadean A., Boukhili R. (2017). Failure mechanisms of a sandwich beam with an ATH/epoxy core under static and dynamic three-point bending. Composite Structures, 176, pp. 281-293. [6] Tumino, D., Ingrassia T., Nigrelli V., Pitarresi G., Urso Miano V. (2014). Mechanical behavior of a sandwich with corrugated GRP core: Numerical modeling and experimental validation. Frattura ed Integrita Strutturale, 30, pp. 317- 326. [7] Rabinovitch, O. (2008). Cohesive interface modeling of debonding failure in FRP strengthened beams. Journal of Engineering Mechanics-Asce, 134(7), pp. 578-588. [8] Judt, P. O., Ricoeur A. (2015). Crack path predictions and experiments in plane structures considering anisotropic properties and material interfaces. Frattura ed Integrita Strutturale, 9(34), pp. 208-215. [9] Rashid, M. M. (1998). The arbitrary local mesh replacement method: An alternative to remeshing for crack propagation analysis. Computer Methods in Applied Mechanics and Engineering, 154(1), pp. 133-150. [10] Mi, Y., Aliabadi M. H. (1994). Three-dimensional crack growth simulation using BEM. Computers & Structures, 52(5), pp. 871-878. [11] Ingraffea, A. R. (2007). Computational Fracture Mechanics. In Encyclopedia of Computational Mechanics , E. Stein, R. Borst and T. J. Hughes, John Wiley & Sons, Ltd. T T M.F. Funari et alii, Frattura ed Integrità Strutturale, 47 (2019) 277-293; DOI: 10.3221/IGF-ESIS.47.21 293 [12] Menouillard, T., Belytschko T. (2010). Smoothed nodal forces for improved dynamic crack propagation modeling in XFEM. International Journal for Numerical Methods in Engineering, 84(1), pp. 47-72. [13] Ciantia, M. O., Arroyo M., Calvetti F., Gens A. (2016). A numerical investigation of the incremental behavior of crushable granular soils. International Journal for Numerical and Analytical Methods in Geomechanics, 40(13), pp. 1773-1798. [14] Daxini, S. D., Prajapati J. M. (2014). A review on recent contribution of meshfree methods to structure and fracture mechanics applications. The Scientific World Journal, 2014. [15] Hyun Moo, K., Hae Sung L., Un Yong J. (1995). An incremental formulation of the moving-grid finite element method for the prediction of dynamic crack propagation. Nuclear Engineering and Design, 158(2), pp. 295-309. [16] Lonetti, P. (2010). Dynamic propagation phenomena of multiple delaminations in composite structures. Computational Materials Science, 48(3), pp. 563-575. [17] Bruno, D., Greco F., Lonetti P. (2013). A fracture-ALE formulation to predict dynamic debonding in FRP strengthened concrete beams. Composites Part B: Engineering, 46, pp. 46-60. [18] Greco, F., Leonetti L., Lonetti P. (2015). A novel approach based on ALE and delamination fracture mechanics for multilayered composite beams. Composites Part B: Engineering, 78, pp. 447-458. [19] Sudhakar, Y., Wall W. A. (2017). Mesh refitting approach: a simple method to model mixed-mode crack propagation in nonlinear elastic solids. Advanced Modeling and Simulation in Engineering Sciences, 4(1), pp. 1-23. [20] Funari, M. F., Lonetti P. (2017). Initiation and evolution of debonding phenomena in layered structures. Theoretical and Applied Fracture Mechanics, 92, pp. 133-145. [21] Funari, M. F., Greco F., Lonetti P. (2018). Sandwich panels under interfacial debonding mechanisms. Composite Structures, 203, pp. 310-320. [22] Funari, M. F., Greco F., Lonetti P. (2016). A moving interface finite element formulation for layered structures. Composites Part B: Engineering, 96, pp. 325-337. [23] Funari, M. F., Greco F., Lonetti P. (2017). Dynamic debonding in layered structures: A coupled ALE-cohesive approach. Frattura ed Integrita Strutturale, 11(41), pp. 524-535. [24] Ma, L., Korsunsky A. M. (2005). On The Use Of Vector J-Integral In Crack Growth Criteria For Brittle Solids. International Journal of Fracture, 133(4), pp. L39-L46. [25] Funari, M. F., Greco F., Lonetti P., Luciano R., Penna R. (2018). An interface approach based on moving mesh and cohesive modeling in Z-pinned composite laminates. Composites Part B: Engineering, 135, pp. 207-217. [26] COMSOL. (2015). COMSOL Multiphysics® v. 5.2, COMSOL AB, Stockholm, Sweden. [27] Ayatollahi, M. R., Aliha M. R. M., Saghafi H. (2011). An improved semi-circular bend specimen for investigating mixed mode brittle fracture. Engineering Fracture Mechanics, 78(1), pp. 110-123. [28] Marsavina, L., Constantinescu D. M., Linul E., Apostol D. A., Voiconi T., Sadowski T. (2014). Refinements on fracture toughness of PUR foams. Engineering Fracture Mechanics, 129, pp. 54-66. [29] Marsavina, L., Linul E., Voiconi T., Constantinescu D. M., Apostol D. A. (2015). On the crack path under mixed mode loading on PUR foams. Frattura ed Integrita Strutturale, 9(34), pp. 387-396. [30] Spadea, S., Orr J., Nanni A., Yang Y. (2017). Wound FRP Shear Reinforcement for Concrete Structures. Journal of Composites for Construction, 21(5). [31] Blaber, J., Adair B., Antoniou A. (2015). Ncorr: Open-Source 2D Digital Image Correlation Matlab Software. Experimental Mechanics, 55(6), pp. 1105-1122. [32] Odessa, I., Frostig Y., Rabinovitch O. (2017 (in press)). Modeling of interfacial debonding propagation in sandwich panels. International Journal of Solids and Structures. [33] Prasad, S., Carlsson L. A. (1994). Debonding and Crack Kinking in Foam Core Sandwich Beams .2. Experimental Investigation. Engineering Fracture Mechanics, 47(6), pp. 825-841. [34] Quispitupa, A., Berggreen C., Carlsson L. A. (2011). Face/core interface fracture characterization of mixed mode bending sandwich specimens. Fatigue and Fracture of Engineering Materials and Structures, 34(11), pp. 839-853. << /ASCII85EncodePages false /AllowTransparency false /AutoPositionEPSFiles true /AutoRotatePages /None /Binding /Left /CalGrayProfile (Dot Gain 20%) /CalRGBProfile (sRGB IEC61966-2.1) /CalCMYKProfile (U.S. Web Coated \050SWOP\051 v2) /sRGBProfile (sRGB IEC61966-2.1) /CannotEmbedFontPolicy /Error /CompatibilityLevel 1.4 /CompressObjects /Tags /CompressPages true /ConvertImagesToIndexed true /PassThroughJPEGImages true /CreateJobTicket false /DefaultRenderingIntent /Default /DetectBlends true /DetectCurves 0.0000 /ColorConversionStrategy /CMYK /DoThumbnails false /EmbedAllFonts true /EmbedOpenType false /ParseICCProfilesInComments true /EmbedJobOptions true /DSCReportingLevel 0 /EmitDSCWarnings false /EndPage -1 /ImageMemory 1048576 /LockDistillerParams false /MaxSubsetPct 100 /Optimize true /OPM 1 /ParseDSCComments true /ParseDSCCommentsForDocInfo true /PreserveCopyPage true /PreserveDICMYKValues true /PreserveEPSInfo true /PreserveFlatness true /PreserveHalftoneInfo false /PreserveOPIComments true /PreserveOverprintSettings true /StartPage 1 /SubsetFonts true /TransferFunctionInfo /Apply /UCRandBGInfo /Preserve /UsePrologue false /ColorSettingsFile () /AlwaysEmbed [ true ] /NeverEmbed [ true ] /AntiAliasColorImages false /CropColorImages true /ColorImageMinResolution 300 /ColorImageMinResolutionPolicy /OK /DownsampleColorImages true /ColorImageDownsampleType /Bicubic /ColorImageResolution 300 /ColorImageDepth -1 /ColorImageMinDownsampleDepth 1 /ColorImageDownsampleThreshold 1.50000 /EncodeColorImages true /ColorImageFilter /DCTEncode /AutoFilterColorImages true /ColorImageAutoFilterStrategy /JPEG /ColorACSImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /ColorImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /JPEG2000ColorACSImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /JPEG2000ColorImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /AntiAliasGrayImages false /CropGrayImages true /GrayImageMinResolution 300 /GrayImageMinResolutionPolicy /OK /DownsampleGrayImages true /GrayImageDownsampleType /Bicubic /GrayImageResolution 300 /GrayImageDepth -1 /GrayImageMinDownsampleDepth 2 /GrayImageDownsampleThreshold 1.50000 /EncodeGrayImages true /GrayImageFilter /DCTEncode /AutoFilterGrayImages true /GrayImageAutoFilterStrategy /JPEG /GrayACSImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /GrayImageDict << /QFactor 0.15 /HSamples [1 1 1 1] /VSamples [1 1 1 1] >> /JPEG2000GrayACSImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /JPEG2000GrayImageDict << /TileWidth 256 /TileHeight 256 /Quality 30 >> /AntiAliasMonoImages false /CropMonoImages true /MonoImageMinResolution 1200 /MonoImageMinResolutionPolicy /OK /DownsampleMonoImages true /MonoImageDownsampleType /Bicubic /MonoImageResolution 1200 /MonoImageDepth -1 /MonoImageDownsampleThreshold 1.50000 /EncodeMonoImages true /MonoImageFilter /CCITTFaxEncode /MonoImageDict << /K -1 >> /AllowPSXObjects false /CheckCompliance [ /None ] /PDFX1aCheck false /PDFX3Check false /PDFXCompliantPDFOnly false /PDFXNoTrimBoxError true /PDFXTrimBoxToMediaBoxOffset [ 0.00000 0.00000 0.00000 0.00000 ] /PDFXSetBleedBoxToMediaBox true /PDFXBleedBoxToTrimBoxOffset [ 0.00000 0.00000 0.00000 0.00000 ] /PDFXOutputIntentProfile () /PDFXOutputConditionIdentifier () /PDFXOutputCondition () /PDFXRegistryName () /PDFXTrapped /False /CreateJDFFile false /Description << /ARA /BGR /CHS /CHT /CZE /DAN /DEU /ESP /ETI /FRA /GRE /HEB /HRV (Za stvaranje Adobe PDF dokumenata najpogodnijih za visokokvalitetni ispis prije tiskanja koristite ove postavke. Stvoreni PDF dokumenti mogu se otvoriti Acrobat i Adobe Reader 5.0 i kasnijim verzijama.) /HUN /ITA /JPN /KOR /LTH /LVI /NLD (Gebruik deze instellingen om Adobe PDF-documenten te maken die zijn geoptimaliseerd voor prepress-afdrukken van hoge kwaliteit. De gemaakte PDF-documenten kunnen worden geopend met Acrobat en Adobe Reader 5.0 en hoger.) /NOR /POL /PTB /RUM /RUS /SKY /SLV /SUO /SVE /TUR /UKR /ENU (Use these settings to create Adobe PDF documents best suited for high-quality prepress printing. Created PDF documents can be opened with Acrobat and Adobe Reader 5.0 and later.) >> /Namespace [ (Adobe) (Common) (1.0) ] /OtherNamespaces [ << /AsReaderSpreads false /CropImagesToFrames true /ErrorControl /WarnAndContinue /FlattenerIgnoreSpreadOverrides false /IncludeGuidesGrids false /IncludeNonPrinting false /IncludeSlug false /Namespace [ (Adobe) (InDesign) (4.0) ] /OmitPlacedBitmaps false /OmitPlacedEPS false /OmitPlacedPDF false /SimulateOverprint /Legacy >> << /AddBleedMarks false /AddColorBars false /AddCropMarks false /AddPageInfo false /AddRegMarks false /ConvertColors /ConvertToCMYK /DestinationProfileName () /DestinationProfileSelector /DocumentCMYK /Downsample16BitImages true /FlattenerPreset << /PresetSelector /MediumResolution >> /FormElements false /GenerateStructure false /IncludeBookmarks false /IncludeHyperlinks false /IncludeInteractive false /IncludeLayers false /IncludeProfiles false /MultimediaHandling /UseObjectSettings /Namespace [ (Adobe) (CreativeSuite) (2.0) ] /PDFXOutputIntentProfileSelector /DocumentCMYK /PreserveEditing true /UntaggedCMYKHandling /LeaveUntagged /UntaggedRGBHandling /UseDocumentProfile /UseDocumentBleed false >> ] >> setdistillerparams << /HWResolution [2400 2400] /PageSize [612.000 792.000] >> setpagedevice