Energy minimising multi-crack growth in linear-elastic materials using the extended finite element method


In computational fracture mechanics as applied, for example, to damage tolerance assessment, it has been common practice to determine the onset and the direction of fracture propagation by post-processing the solution of the linear elastostatics problem for a given fracture configuration. Under mixed-mode loading conditions, the classic criteria used in fracture mechanics (FM) rely on the local crack tip field to determine if a crack should advance or not and in which direction. Commonly used criteria, such as maximum hoop-stress, principle of local symmetry and maximum strain energy density, are ad hoc and although they are relatively simple to apply and give reasonably good accuracy for typical fracture mechanics problems, they lack a convincing physical basis.


The aim is to develop a robust and thermodynamically consistent approach to the modelling of arbitrary crack growth under general loading conditions, e.g. body loading, thermal loading, and/or crack surface pressure loading; then, to apply this model to simulate the process of silicon wafer splitting.


We investigate multiple crack evolution under quasi-static conditions in an isotropic linear elastic solid based on the principle of minimum total energy (i.e. the sum of the potential and the fracture surface energies). The present technique enables the minimisation of the total energy functional with respect to all crack increments (taking into account their relative interactions). The crack extension directions are optimized (in the minimum energy sense) when the energy release rates associated with their virtual rotations vanish. The 2D continuum is discretised using the extended finite element method (XFEM) to enable mesh independent modelling of multiple fractures with arbitrary intersections.


An explicit crack representation is used to precisely capture crack intersections as they arise during their evolution. The energy release rates (and their rates) associated with a set of crack tip extensions at a given time are computed algebraically within the XFEM framework, based on the stiffness derivative approach, which involves geometrically perturbing mesh patches around the crack tips.


Although XFEM facilitates mesh independent modelling of cracks, the enrichment needs to be updated after every time-step. In the current implementation this is achieved by means of a systematic book-keeping of the element enrichment data, addition and removal of enrichment only where necessary, and a consistent updating of the global system of equations. Consequently, moderate computational times are obtained, even in our Matlab implementation. In the benchmark problems that we solve the largest computational expense by far is in the solution to the linear system of equations rather than in the assembly, updating or post-processing of the numerical solution. Significant speed-ups can be achieved for fracture problems involving many cracks that evolve over many time steps.

This slideshow requires JavaScript.


Upon comparing the fracture path solutions by different criteria (namely: the minimum energy criterion and the maximum hoop-stress criterion) it can be observed in the test cases that have been attempted that both criteria lead to very similar fracture path solutions. Interestingly, the  solutions appear to convergence with mesh refinement towards each other but from opposite directions.

In other words, it is found that the converged fracture path lies in between the fracture paths obtained by each criterion on coarser meshes. This indicates that the respective criteria tend to over/under estimate the crack deflections for a given discretisation. The consistency of the observations in the numerical solutions leads one to propose a possible improvement to the numerical solution. The modified crack growth direction is obtained by averaging the solutions obtained by the two criteria.

This slideshow requires JavaScript.

Consequently, significant improvements in accuracy and convergence rate of fracture paths is obtained, especially on coarse meshes. Finally, the proposed method are applied to simulate the silicon wafer splitting process. A parametric study is undertaken to quantify the influence of Smart-Cut process parameters on the post-split fracture surface roughness. Good agreement between the numerical and the experimental results can be obtained.

This slideshow requires JavaScript.