- Email: [email protected]
Effect of sawing induced micro-crack orientations on fracture properties of silicon wafers
Accepted Manuscript Effect of Sawing Induced Micro-Crack Orientations on Fracture Properties of Silicon Wafers Amin S. Azar, Børge Holme, Øyvind Nielsen PII: DOI: Reference:
S0013-7944(16)00017-5 http://dx.doi.org/10.1016/j.engfracmech.2016.01.014 EFM 5030
To appear in:
Engineering Fracture Mechanics
Received Date: Revised Date: Accepted Date:
3 August 2015 10 January 2016 11 January 2016
Please cite this article as: Azar, A.S., Holme, B., Nielsen, O., Effect of Sawing Induced Micro-Crack Orientations on Fracture Properties of Silicon Wafers, Engineering Fracture Mechanics (2016), doi: http://dx.doi.org/10.1016/ j.engfracmech.2016.01.014
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Effect of Sawing Induced Micro-Crack Orientations on Fracture Properties of Silicon Wafers Amin S. Azar a,1, Børge Holme a, Øyvind Nielsen b (a) SINTEF Materials and Chemistry, Forskiningsveien 1, NO-0373 Oslo, Norway. (b) Norsun AS, Karenlyst Allé 9C, NO-0278 Oslo, Norway.
Fixed abrasive wire sawing has become the method of choice for slicing single crystal silicon ingots. During the wafering process, diamond grains penetrate into the silicon while sliding along the material at the speed of several meters per second. This leaves a significant amount of sub-surface damage with micro-cracks near the wafer surface that are oriented in different directions. The micro-crack distribution depends on the general sawing parameters, and on intricate, local details of the sawing process. Handling of sliced wafers with microcracks may provide conditions for unstable growth that leads to brittle fracture. In order to categorize the effect of these micro-cracks on the overall fracture strength of the wafers, a number of cracks with different lengths and orientations were modelled. Four-point-bending of a 150 µm thick wafer was simulated. Stress intensity factor (KC) and J-Integral (JC) approaches were used to define and sort unstable crack lengths and morphologies. It was found that the cracks oriented at 20° and 60° with reference to 100 are the most and least vulnerable flaws respectively, under the identical loading conditions. Keywords: Abrasive wire sawing, Silicon, Wafer, Subsurface Cracks, Fracture mechanics.
1. Introduction The wafering technology of silicon ingots has been investigated and developed for a relatively long period. The fixed abrasive wire sawing process has been applied since the mid-1990s, as an alternative to slurry sawing technology. In slurry sawing, a smooth steel wire is used together with slurry that contains abrasive particles. The fixed abrasive sawing technology uses a diamond grit coated stainless steel wire to reduce the kerf loss, reduce the thickness variations and increase the efficiency [1]. During the sawing process, imposed local stresses can cause subsurface micro-cracks at different angles and orientations [2]. These cracks are responsible for degraded mechanical properties of the produced wafers and can cause postproduction handling problems due to wafer brittleness [3]. The target in the wafering industry is to reduce the thickness in order to mitigate the material costs. However, silicon is a brittle material that fractures elastically, so subsurface cracks introduced in the processing stage can cause major issues in reaching this target [4]. A large number of investigations have been performed to relate the mechanical properties to the sawing induced features like subsurface cracks and saw marks [5-8]. In general, these defects may align in many spatial orientations that can be defined by two characteristic angles. Fig. 1 shows an arbitrary crack that is assumed to have a planar surface. Since a
1
Corresponding author at SINTEF, Institute of Materials and Chemistry. Address: Forskningsveien 1, 0373, Oslo, Norway. Email: [email protected]. Phone: +4796833668.
1
silicon single crystal mostly cracks along low-index crystallographic planes and directions under loading, it can be assumed that the crack surfaces are near planar [9].
Fig. 1: Spatial orientation illustration of a sawing induced crack.
As shown in Fig. 1, the sawing induced crack can be described by angle α on the sawing plane and angle β on the perpendicular plane. Wu et al. [2] employed a finite element model to benchmark the effect of the β angle on the critical crack length ratio and fracture strength of the silicon wafer under uniaxial loading. The critical crack length ratio is defined as the ratio of the critical length of an angled crack with orientation divided by the critical length of a pure mode I crack (β =0) for the same applied stress. Wu et al. showed that the critical crack length ratio increases drastically when the β angle is more than 45°. According to their conclusion, the ratio is almost constant from β=0° to 45°. This means that the orientation of subsurface cracks from 0° to 45° will result in the same fracture strength level for a given constant loading condition. Wu's work was done under uniaxial loading conditions. Yang et al. [10] conducted a comparative analysis of two major sawing technologies by modelling the silicon wafer under four-point bending. The focus of their work was on the effect of saw marks on the fracture strength of the wafer. According to Yang's conclusions, the maximum principal stress was increasing as the saw mark/defect radius increased. Their work did not consider the effect of saw mark orientation on the fracture strength of the single crystal silicon wafer. As mentioned earlier, silicon fractures along the planes that are of high atomic packing factor. Three major planes are {111}, {110} and {100}. The fracture strength of these planes can be analyzed by the critical stress intensity factor (KIC) and critical fracture energy (JIC). The stress intensity factor (SIF) describes the intensity of elastic crack-tip fields, and symbolizes the linear elastic fracture mechanics (LEFM). The J-integral, however, is proposed to characterize the intensity of elastic–plastic crack-tip fields, and represents the elastic–plastic fracture mechanics (EPFM). A detailed description of these two concepts is given in [11]. The SIF and J-integral concepts were used by many researchers to describe the fracture resistance of the crystallographic planes in the silicon single crystal. The reported KIC values were 2
mostly obtained using the micro-hardness indentation tests [12, 13] for determining the fracture toughness of brittle materials. The J-integral fracture energy can be determined by density functional theory (DFT) and molecular dynamics (MD) numerical simulations. Pérez and Gumbsch [14] have performed an ab initio study of the brittle fracture anisotropy in Silicon. They concluded that crack propagates easily on {111}. They have observed instances that propagation on {110} planes are discontinued due to lattice trapping. According to them, the propagation of crack in silicon single crystal is highly anisotropic. Ding et al. [15] considered using density functional theory to study the fracture toughness in a number of brittle materials. They calculated theoretical fracture stress and theoretical fracture toughness and concluded that the role of narrow range Poisson's ratio in determining fracture toughness is much smaller than shear modulus and surface energy, meaning that anisotropic materials with large variation in Poisson's ratio may not show direct proportionality to (γsG)½. As a result of a joint effort, Masolin et al. [16] gathered all related data from previous researches and suggested recommended KIC and JIC values for the three main atomic planes in silicon crystals. The calculated values for fracture toughness on different planes were mainly performed through experiments, however, the molecular dynamics and density functional theory approaches could report both fracture toughness and fracture energy values. The recommended values are used in the current study to benchmark the crack orientation effects. In this study, the aim is to present the effect of the α-angle on the fracture strength of the single crystal silicon wafer under four-point bending conditions by modelling the relevant phenomena, including the numerical simulation of stress intensity factors and energy release rates at the crack tip for various lengths. 2. Method In order to study the effect of crack orientation, a model set-up for four-point bending was defined in Abaqus 6.13 according to Fig 2. The wafer was square and 125 mm on each side. The average thickness of the wafer was 150 µm. The inner moving rollers were set 65 mm apart from each other, and the fixed rollers were placed at 10 mm from the edges of the wafer according to the guidelines of ASTM C1161 [17]. Four-point bending test provides a constant bending momentum and force distribution between the inner rollers that sets lower geometrical contraints in the specimen [18].
3
Fig. 2: 3D model of four-point bending set-up.
A benchmarking procedure for the α-angle orientation was defined. The cracks from α=0° to 90° with increments of 10° was considered (including the special case α=45°). The computational time for running an extensive fracture strength benchmarking on the 3D model was relatively long. Therefore, an axisymmetric 2D model was used for the calculations. Fig. 3 shows the boundary conditions for the 2D model.
Fig. 3: 2D model for four-point bending with applied boundary conditions.
It was assumed that only one crack exists, right on the symmetry plane, where maximum flexural stresses develop during bending. It was also assumed that the β angle was 0° for all simulation cases. The physical rotation of the crack for varying the α angle is not possible in the 2D model. Instead of rotating the crack, the mechanical anisotropy matrix, determined by the crystallographic orientation, was transformed to replicate the rotated properties. The following matrix is the anisotropic engineering constants of the single crystal silicon [19]: Eq. 1 Eq. 2
4
where E is the elastic modulus, G is the shear modulus and ν is the Poisson's ratio. The following rotation matrix was used to transform the engineering constants for each α angle: Eq. 3
The rotation of the α angle is restricted to the (x, y) plane considering that z direction is along 100, parallel to the sawing plane (see Fig. 1). Fig. 4 shows the relative rotation of the coordinate system that is used for transforming the engineering constants.
Fig. 4: Coordination system rotation analogy.
Based on the aforementioned rotation conditions, new anisotropic engineering constants were calculated and shown in Fig. 5 for each α angle. Note that the elastic modulus for the z direction does not change, since there is no rotation around that axis.
Fig. 5: Calculated engineering constants for rotated matrix in the xy plane. The E and G values are in GPa.
The crack feature was introduced in the vicinity of the symmetry plane as illustrated in Fig. 6a.
5
The center of the wafer in a four-point bending test experiences the largest displacement. The intention of placing the crack on the centerline was to simulate the worst case scenario. The maximum applied displacement on the moving rollers was to reach 14.5 mm in the centerline. The J-integral and stress intensity factor were calculated at the crack tip for each crack length along the displacement path. For each crack orientation (α angle), 10 different crack lengths were investigated: 1 µm, 2 µm, 3 µm, 4 µm, 5 µm, 10 µm, 15 µm, 20 µm, 25 µm and 30 µm. Fig. 6b shows the mesh quality along the thickness of the wafer. The quad structured element size at the crack tip is 0.33 microns and this size is applied for 10 micron length ahead of the crack tip. These values were kept constant for all crack lengths and orientations.
Fig. 6: (a) Crack positioning and boundary conditions and (b) Mesh elements quality at the crack tip.
3. Results and Discussion Fig. 7 shows the crack tip stress field at the maximum bending in the center of the wafer for different crack lengths. As illustrated, the stress field expands as the crack length increases. The size of the process zone at the crack tip for cracks below 5 µm is very small, though the stress intensity factor and fracture energy can vary depending on the α-angle.
Fig. 7: von Mises stress field (in Pa) at the crack tip for four different crack lengths. The rectangle shows the wafer thickness and the circle inside shows the magnified crack tip field. The α-angle is 0°.
3.1. EPFM approach In the theory of elastic-plastic fracture mechanics, J-integral is used to calculate the strain energy release rate. This corresponds to the amount of energy required to propagate the crack per unit surface area. The quasi-static analysis of fracture phenomena can be defined in 3D and 3D. In this analysis, two-dimensional geometry is considered. J-integral is defined as: 6
Eq. 4
where H is: Eq. 5 Abaqus defines the domain in terms of rings of elements surrounding the crack tip. Different “contours” (domains) are created. Each subsequent contour is defined by adding the next ring of elements that share nodes with the elements in the previous contour. q is chosen to have a magnitude of zero at the nodes on the outside of the contour and to be one (in the crack direction) at all nodes inside the contour. Middle nodes have values between 0 and 1 depending on their distance to the crack tip. Them the J-integral is calculated for these paths separately while they should result in the same value according to the theory. Fig. 8 presents the results of J-integral calculations for two α-angles. Each line, as annotated on the figure, represents the J-integral value for a given amount of bending displacement in the center of the wafer under four-point flexural loading. The horizontal lines represent the recommended JIC values from the literature, for three low-index planes. The {111} plane in silicon has the highest atomic packing factor. Therefore, the required energy for generating new surfaces by planar cracking is lower than for the other two planes. The intersections of Jintegral and JIC lines show the maximum tolerable crack length prior to brittle fracturing, under the imposed amount of bending displacement. For example, if the wafer has a crack along α = 0° located in the center of the wafer, and the wafer is bent down 5 mm, the Jintegral calculation gives a maximum tolerable crack length of about 20 µm. If the α-angle is 60° and all other conditions are equal (right-hand chart of Fig 8), then the maximum allowable crack length is around 30 µm. This means that should brittle fracture happen on a {111} plane, the crack orientation of 60° would be more resistant than 0°. Fig. 9 summarizes the effect of the α-angle on the critical crack length. It includes low, medium and high levels of bending displacement (5.0 mm, 10.0 mm and 14.5mm, respectively). The variation of critical crack length for different crack orientations is more pronounced for low of bending displacements. The critical crack length tends to flatten out at high bending displacements, meaning that the stress concentration at the crack tip will already be higher than the minimum required energy for making new surfaces on {111} planes for all α-angle values. In other words, the effect of crack orientation is more significant in brittle fracture when the bending displacement is small. For larger bending displacements, the brittle fracture is almost equally critical for all crack orientations. At low bending displacement values, the critical crack length is the longest if the crack is oriented 60° and shortest if the α-angle is 20°. This relationship remains unchanged for medium and high bending displacement, but the variation is less conspicuous.
7
Fig. 8: J-integral calculation for α = 0° and 60°.
Fig. 9: Benchmarking of the α-angle's effect on the maximum allowable crack length.
3.2. LEFM approach Fig. 10 shows the stress intensity factors for two α angle values. The plots are related to Fig. 8. However, here the lines represent the crack length rather than the bending displacement, and the x-axis represents the wafer bending displacement. The recommended KIC values are extracted from [16, 20]. The intersection of KIC and SIF lines for each crack length will give the maximum allowable bending displacement in the center of the wafer.
8
Fig. 10: Stress intensity factor (SIF) calculation for α = 0° and 60°.
Fig. 11a shows how the critical bending displacement decreases as the crack length increases. As demonstrated in the figure, the crack orientation becomes less important when the crack is large. This is in agreement with findings from the J-integral analysis. Fig. 11b depicts the benchmarking results of SIF calculations for different crack orientations. The SIF analysis, however, suggests different results compared to the J-integral analysis that is shown in Fig. 9. According to Fig. 11b, if the α-angle is 30°, the wafer can tolerate the highest bending displacement, whereas if α is 0° or 90°, the maximum allowable bending displacement is the lowest. Further, the maximum allowable bending for the same crack length does not correlate between these two analyses. For example, let's say the maximum allowable bending is assumed to be 10.0 mm. At 0°, according to Fig. 9, the critical crack length is just below 5 µm, while in Fig. 11b, it is 10 µm. The SIF calculation in this case overestimates the allowable crack size by a factor of 2, compared to the value suggested by the J-integral approach.
Fig. 11: (a) Critical bending displacement and (b) α-angle benchmarking for different crack lengths, based on SIF calculations.
Hence, there are two major differences between J-integral and SIF approaches: 9
Difference in the weakest and strongest fracture crack orientation. Difference in the maximum allowable crack lengths and bending displacements.
In order to argue the probable reasons for these differences, a wafer bending experimental setup was utilized to measure the maximum bending displacement values. Fig. 12a shows the results of the designed experiment. The bending tests were performed along two perpendicular directions: 100 and 010, which are along and perpendicular to the saw marks respectively. The bending setup was a simple manual machine that was pushing the plastic anvils using a manually controlled screw. The screw was rotated at a similar low rate for all tests. Once the wafer shattered, the bending was ceased and the total displacement was measured with reference to the initial state of the screw before bending. Each presented curve is constructed from the fracture data of one hundred wafers bent under the same conditions. The maximum bending displacement at 50 % of cumulative percentage is about 8 mm for the weak direction. This value increases to about 22 mm when bending perpendicular to the sawing direction. The weakest direction can tolerate no more than 11-12 mm of bending displacement before all wafers will break. Fig. 12b shows an optical image from a tapered cross section of a wafer. The maximum crack depths after sawing are observed to be around 4 µm. The combination of experimental maximum crack depth and maximum bending displacement is in agreement with the J-integral calculations that were presented earlier. According to Fig. 9, if the wafer can tolerate 10 mm of bending displacement, the maximum crack length to initiate fracture along {111} is just below 5 µm if the crack is oriented mainly along the saw marks (α-angle < 30°).
Fig. 12: (a) Experimental measurement of bending displacement along and perpendicular to the saw marks and (b) Optical microscope image of wafer cross section showing the damage profile.
The applicability of LEFM and EPFM was explored in the development of fracture mechanics methods. For measuring and calculating the KIC under LEFM conditions, the size of the sample with respect to the crack tip plastic region plays an important role, thus, the sample geometry should be used for defining the critical SIF [21]. The suggested sizes are standardized and being practiced in the fracture mechanics experiments and simulations. The assumptions of linear elastic fracture mechanics may not hold if the plastic zone at the crack 10
tip has the size of the same order of magnitude as the crack size [22]. Furthermore, the plastic zone may change as the applied load and the crack size are increased. In other words, the LEFM approach is not suitable for short crack lengths and thin sections, both of which apply to the current study of silicon single crystal wafers. The performed simulations are excluding the involvement of dislocation role in the fracture mechanisms. The generation of interstitial dislocation loops (Franks Partials) and shear faults are the proposed mechanisms that may assist fracture initiation at atomic level [23]. Samuels and Roberts [24] have shown that silicon single crystal show some dislocation activity above ductile to brittle transition temperature, which is about 10°C depending on the type of doping. The dislocation mobility is a function of the local stress intensity that are studied and benchmarked in this paper. The Mode I fracture toughness is correlated with the J-Integral value under monotonic loading conditions; however, this relation is only for isotropic and perfectly brittle linear elastic materials under plane strain deformation [25]. Therefore, it can be inferred that the J-integral calculations are applicable for linear elastic, nonlinear elastic and elastic-plastic materials, unlike stress intensity factor calculation, which is applicable for cracks in homogeneous, linear elastic materials [26]. In order to calculate the stress intensity factor, the crack propagation direction must be defined while it is valid only for homogeneous, isotropic linear elastic materials. Calculation of J-integral is insensitive to crack propagation path as it is a path independent method. This explains the discrepancy between the results of SIF and J-integral calculations. 4. Conclusions Based on simulation and experimental results, the following conclusions can be deduced.
The elastic-plastic fracture mechanics (EPFM) approach shows consistent results when compared to the experimental test outcome. It is suggested as a suitable tool for analyzing the effect of sawing induced cracks on the structural integrity of silicon single crystal wafers. Based on the J-integral calculations, the most critical crack orientation with respect to the 100 is around 20°. If the sawing process is optimized to avoid generation of subsurface micro-cracks that are aligned close to 20°, the wafers should tolerate more bending displacement. On the other hand, if the crack orientation is at 50-60°, the brittle fracture along {111} planes becomes less susceptible, and the wafer can tolerate the highest bending displacement. The calculation of SIF using the LEFM approach is not a suitable option for analyzing the fracture properties of thin silicon wafers. It will overestimate the maximum allowable bending displacement and maximum critical crack sizes. EPFM simulation has produced more reliable results in this regard.
11
5. Further Studies The exact morphology of the cracks and its relation to the sawing parameters are not yet well understood. We have observed that the cracking initiates at different sites and grow and coalesce in milliseconds, which lead to the brittle shattering of the loaded wafer. More experimental and simulation work is required to fully understand the fracture mechanisms in silicon single crystal wafers. This may eventually lead to optimized sawing parameters that can help increase the structural integrity and fracture strength of the wafers. 6. Acknowledgment The authors would like to gratefully acknowledge the Norwegian Research Council for cofunding this research work under the contract number of 219395/O30. 7. References [1] Hardin CW, Qu J, Shih AJ. Fixed abrasive diamond wire saw slicing of single-crystal silicon carbide wafers. Mater Manuf Process. 2004;19:355-67. [2] Wu H, Melkote SN, Danyluk S. Mechanical strength of silicon wafers cut by loose abrasive slurry and fixed abrasive diamond wire sawing. Adv Eng Mater. 2012;14:342-8. [3] Brun XF, Melkote SN. Analysis of stresses and breakage of crystalline silicon wafers during handling and transport. Sol Energ Mat Sol C. 2009;93:1238-47. [4] Cook R. Strength and sharp contact fracture of silicon. J Mater Sci. 2006;41:841-72. [5] Ferrazza F. Crystalline silicon: Manufacture and properties. Practical Handbook of Photovoltaics: Fundamentals and Applications. 2012:79. [6] Coletti G, Tool C, Geerligs L. Mechanical strength of silicon wafers and its modelling. Proc of Workshop on Crystalline Silicon Solar Cells and Modules: Materials and Processes: Citeseer; 2005. p. 117-20. [7] Coletti G, Tool C, Geerligs L. Quantifying surface damage by measuring mechanical strength of silicon wafers. Presented at the 20th European Photovoltaic Solar Energy Conference and Exhibition2005. p. 10. [8] Möller HJ. Basic Mechanisms and Models of Multi‐Wire Sawing. Adv Eng Mater. 2004;6:501-13. [9] Cramer T, Wanner A, Gumbsch P. Energy dissipation and path instabilities in dynamic fracture of silicon single crystals. Phys Rev Lett. 2000;85:788. [10] Yang C, Wu H, Melkote S, Danyluk S. Comparative analysis of fracture strength of slurry and diamond wire sawn multicrystalline silicon solar wafers. Adv Eng Mater. 2013;15:358-65. [11] Zhu X-K, Joyce JA. Review of fracture toughness (G, K, J, CTOD, CTOA) testing and standardization. Eng Fract Mech. 2012;85:1-46. [12] Ponton CB, Rawlings RD. Vickers indentation fracture toughness test Part 1 Review of literature and formulation of standardised indentation toughness equations. Mater Sci TechLond. 1989;5:865-72. [13] Quinn GD, Bradt RC. On the Vickers indentation fracture toughness test. J Am Ceram Soc. 2007;90:673-80. [14] Pérez R, Gumbsch P. An ab initio study of the cleavage anisotropy in silicon. Acta Mater. 2000;48:4517-30. [15] Ding Z, Zhou S, Zhao Y. Hardness and fracture toughness of brittle materials: a density functional theory study. Phys Rev B. 2004;70:184117. 12
[16] Masolin A, Bouchard P-O, Martini R, Bernacki M. Thermo-mechanical and fracture properties in single-crystal silicon. J Mater Sci. 2013;48:979-88. [17] Standard A. C1161, 2002c,“Standard Test Method for Flexural Strength of Advanced Ceramics at Ambient Temperature”, ASTM International, West Conshohocken, PA, 2002, DOI: 10.1520/C1161-02C. [18] Popovich V, Riemslag A, Janssen M, Bennett I, Richardson I. Characterization of Multicrystalline Silicon Solar Wafers Fracture Strength and Influencing Factors. International Journal of Material Science. 2013. [19] Hopcroft MA, Nix WD, Kenny TW. What is the Young's Modulus of Silicon? Microelectromechanical Systems, Journal of. 2010;19:229-38. [20] Sen D, Thaulow C, Schieffer SV, Cohen A, Buehler MJ. Atomistic study of crack-tip cleavage to dislocation emission transition in silicon single crystals. Phys Rev Lett. 2010;104:235502. [21] Anderson TL, Anderson T. Fracture mechanics: fundamentals and applications: CRC press; 2005. [22] Miller K. The short crack problem. Fatigue Fract Eng M. 1982;5:223-32. [23] Ebrahimi F, Kalwani L. Fracture anisotropy in silicon single crystal. Materials Science and Engineering: A. 1999;268:116-26. [24] Samuels J, Roberts S. The brittle-ductile transition in silicon. I. Experiments. Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences: The Royal Society; 1989. p. 1-23. [25] Bower AF. Applied mechanics of solids, Chapter 9: CRC press; 2009. [26] Parks D. The virtual crack extension method for nonlinear material behavior. Computer Methods in Applied Mechanics and Engineering. 1977;12:353-64.
13
List of Nomenclature KC KIC
JC JIC ν E G α β x, y, z x', y', z' Γ q n W I σ u
Calculated stress intensity factor Mode I critical stress intensity factor Calculated J-integral Mode I critical J-integral Poisson ratio Elastic modulus Shear modulus Angle between saw-marks and the crack on the sawing plane Angle between sawing plane normal vector and the crack Original coordination system axes and directions Rotated coordination system axes and directions Path contour around the crack Unit vector in the virtual crack extension direction Outward normal to Γ Elastic strain energy Propagation length Stress level Displacement
14
S0013-7944(16)00017-5 http://dx.doi.org/10.1016/j.engfracmech.2016.01.014 EFM 5030
To appear in:
Engineering Fracture Mechanics
Received Date: Revised Date: Accepted Date:
3 August 2015 10 January 2016 11 January 2016
Please cite this article as: Azar, A.S., Holme, B., Nielsen, O., Effect of Sawing Induced Micro-Crack Orientations on Fracture Properties of Silicon Wafers, Engineering Fracture Mechanics (2016), doi: http://dx.doi.org/10.1016/ j.engfracmech.2016.01.014
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
Effect of Sawing Induced Micro-Crack Orientations on Fracture Properties of Silicon Wafers Amin S. Azar a,1, Børge Holme a, Øyvind Nielsen b (a) SINTEF Materials and Chemistry, Forskiningsveien 1, NO-0373 Oslo, Norway. (b) Norsun AS, Karenlyst Allé 9C, NO-0278 Oslo, Norway.
Fixed abrasive wire sawing has become the method of choice for slicing single crystal silicon ingots. During the wafering process, diamond grains penetrate into the silicon while sliding along the material at the speed of several meters per second. This leaves a significant amount of sub-surface damage with micro-cracks near the wafer surface that are oriented in different directions. The micro-crack distribution depends on the general sawing parameters, and on intricate, local details of the sawing process. Handling of sliced wafers with microcracks may provide conditions for unstable growth that leads to brittle fracture. In order to categorize the effect of these micro-cracks on the overall fracture strength of the wafers, a number of cracks with different lengths and orientations were modelled. Four-point-bending of a 150 µm thick wafer was simulated. Stress intensity factor (KC) and J-Integral (JC) approaches were used to define and sort unstable crack lengths and morphologies. It was found that the cracks oriented at 20° and 60° with reference to 100 are the most and least vulnerable flaws respectively, under the identical loading conditions. Keywords: Abrasive wire sawing, Silicon, Wafer, Subsurface Cracks, Fracture mechanics.
1. Introduction The wafering technology of silicon ingots has been investigated and developed for a relatively long period. The fixed abrasive wire sawing process has been applied since the mid-1990s, as an alternative to slurry sawing technology. In slurry sawing, a smooth steel wire is used together with slurry that contains abrasive particles. The fixed abrasive sawing technology uses a diamond grit coated stainless steel wire to reduce the kerf loss, reduce the thickness variations and increase the efficiency [1]. During the sawing process, imposed local stresses can cause subsurface micro-cracks at different angles and orientations [2]. These cracks are responsible for degraded mechanical properties of the produced wafers and can cause postproduction handling problems due to wafer brittleness [3]. The target in the wafering industry is to reduce the thickness in order to mitigate the material costs. However, silicon is a brittle material that fractures elastically, so subsurface cracks introduced in the processing stage can cause major issues in reaching this target [4]. A large number of investigations have been performed to relate the mechanical properties to the sawing induced features like subsurface cracks and saw marks [5-8]. In general, these defects may align in many spatial orientations that can be defined by two characteristic angles. Fig. 1 shows an arbitrary crack that is assumed to have a planar surface. Since a
1
Corresponding author at SINTEF, Institute of Materials and Chemistry. Address: Forskningsveien 1, 0373, Oslo, Norway. Email: [email protected]. Phone: +4796833668.
1
silicon single crystal mostly cracks along low-index crystallographic planes and directions under loading, it can be assumed that the crack surfaces are near planar [9].
Fig. 1: Spatial orientation illustration of a sawing induced crack.
As shown in Fig. 1, the sawing induced crack can be described by angle α on the sawing plane and angle β on the perpendicular plane. Wu et al. [2] employed a finite element model to benchmark the effect of the β angle on the critical crack length ratio and fracture strength of the silicon wafer under uniaxial loading. The critical crack length ratio is defined as the ratio of the critical length of an angled crack with orientation divided by the critical length of a pure mode I crack (β =0) for the same applied stress. Wu et al. showed that the critical crack length ratio increases drastically when the β angle is more than 45°. According to their conclusion, the ratio is almost constant from β=0° to 45°. This means that the orientation of subsurface cracks from 0° to 45° will result in the same fracture strength level for a given constant loading condition. Wu's work was done under uniaxial loading conditions. Yang et al. [10] conducted a comparative analysis of two major sawing technologies by modelling the silicon wafer under four-point bending. The focus of their work was on the effect of saw marks on the fracture strength of the wafer. According to Yang's conclusions, the maximum principal stress was increasing as the saw mark/defect radius increased. Their work did not consider the effect of saw mark orientation on the fracture strength of the single crystal silicon wafer. As mentioned earlier, silicon fractures along the planes that are of high atomic packing factor. Three major planes are {111}, {110} and {100}. The fracture strength of these planes can be analyzed by the critical stress intensity factor (KIC) and critical fracture energy (JIC). The stress intensity factor (SIF) describes the intensity of elastic crack-tip fields, and symbolizes the linear elastic fracture mechanics (LEFM). The J-integral, however, is proposed to characterize the intensity of elastic–plastic crack-tip fields, and represents the elastic–plastic fracture mechanics (EPFM). A detailed description of these two concepts is given in [11]. The SIF and J-integral concepts were used by many researchers to describe the fracture resistance of the crystallographic planes in the silicon single crystal. The reported KIC values were 2
mostly obtained using the micro-hardness indentation tests [12, 13] for determining the fracture toughness of brittle materials. The J-integral fracture energy can be determined by density functional theory (DFT) and molecular dynamics (MD) numerical simulations. Pérez and Gumbsch [14] have performed an ab initio study of the brittle fracture anisotropy in Silicon. They concluded that crack propagates easily on {111}. They have observed instances that propagation on {110} planes are discontinued due to lattice trapping. According to them, the propagation of crack in silicon single crystal is highly anisotropic. Ding et al. [15] considered using density functional theory to study the fracture toughness in a number of brittle materials. They calculated theoretical fracture stress and theoretical fracture toughness and concluded that the role of narrow range Poisson's ratio in determining fracture toughness is much smaller than shear modulus and surface energy, meaning that anisotropic materials with large variation in Poisson's ratio may not show direct proportionality to (γsG)½. As a result of a joint effort, Masolin et al. [16] gathered all related data from previous researches and suggested recommended KIC and JIC values for the three main atomic planes in silicon crystals. The calculated values for fracture toughness on different planes were mainly performed through experiments, however, the molecular dynamics and density functional theory approaches could report both fracture toughness and fracture energy values. The recommended values are used in the current study to benchmark the crack orientation effects. In this study, the aim is to present the effect of the α-angle on the fracture strength of the single crystal silicon wafer under four-point bending conditions by modelling the relevant phenomena, including the numerical simulation of stress intensity factors and energy release rates at the crack tip for various lengths. 2. Method In order to study the effect of crack orientation, a model set-up for four-point bending was defined in Abaqus 6.13 according to Fig 2. The wafer was square and 125 mm on each side. The average thickness of the wafer was 150 µm. The inner moving rollers were set 65 mm apart from each other, and the fixed rollers were placed at 10 mm from the edges of the wafer according to the guidelines of ASTM C1161 [17]. Four-point bending test provides a constant bending momentum and force distribution between the inner rollers that sets lower geometrical contraints in the specimen [18].
3
Fig. 2: 3D model of four-point bending set-up.
A benchmarking procedure for the α-angle orientation was defined. The cracks from α=0° to 90° with increments of 10° was considered (including the special case α=45°). The computational time for running an extensive fracture strength benchmarking on the 3D model was relatively long. Therefore, an axisymmetric 2D model was used for the calculations. Fig. 3 shows the boundary conditions for the 2D model.
Fig. 3: 2D model for four-point bending with applied boundary conditions.
It was assumed that only one crack exists, right on the symmetry plane, where maximum flexural stresses develop during bending. It was also assumed that the β angle was 0° for all simulation cases. The physical rotation of the crack for varying the α angle is not possible in the 2D model. Instead of rotating the crack, the mechanical anisotropy matrix, determined by the crystallographic orientation, was transformed to replicate the rotated properties. The following matrix is the anisotropic engineering constants of the single crystal silicon [19]: Eq. 1 Eq. 2
4
where E is the elastic modulus, G is the shear modulus and ν is the Poisson's ratio. The following rotation matrix was used to transform the engineering constants for each α angle: Eq. 3
The rotation of the α angle is restricted to the (x, y) plane considering that z direction is along 100, parallel to the sawing plane (see Fig. 1). Fig. 4 shows the relative rotation of the coordinate system that is used for transforming the engineering constants.
Fig. 4: Coordination system rotation analogy.
Based on the aforementioned rotation conditions, new anisotropic engineering constants were calculated and shown in Fig. 5 for each α angle. Note that the elastic modulus for the z direction does not change, since there is no rotation around that axis.
Fig. 5: Calculated engineering constants for rotated matrix in the xy plane. The E and G values are in GPa.
The crack feature was introduced in the vicinity of the symmetry plane as illustrated in Fig. 6a.
5
The center of the wafer in a four-point bending test experiences the largest displacement. The intention of placing the crack on the centerline was to simulate the worst case scenario. The maximum applied displacement on the moving rollers was to reach 14.5 mm in the centerline. The J-integral and stress intensity factor were calculated at the crack tip for each crack length along the displacement path. For each crack orientation (α angle), 10 different crack lengths were investigated: 1 µm, 2 µm, 3 µm, 4 µm, 5 µm, 10 µm, 15 µm, 20 µm, 25 µm and 30 µm. Fig. 6b shows the mesh quality along the thickness of the wafer. The quad structured element size at the crack tip is 0.33 microns and this size is applied for 10 micron length ahead of the crack tip. These values were kept constant for all crack lengths and orientations.
Fig. 6: (a) Crack positioning and boundary conditions and (b) Mesh elements quality at the crack tip.
3. Results and Discussion Fig. 7 shows the crack tip stress field at the maximum bending in the center of the wafer for different crack lengths. As illustrated, the stress field expands as the crack length increases. The size of the process zone at the crack tip for cracks below 5 µm is very small, though the stress intensity factor and fracture energy can vary depending on the α-angle.
Fig. 7: von Mises stress field (in Pa) at the crack tip for four different crack lengths. The rectangle shows the wafer thickness and the circle inside shows the magnified crack tip field. The α-angle is 0°.
3.1. EPFM approach In the theory of elastic-plastic fracture mechanics, J-integral is used to calculate the strain energy release rate. This corresponds to the amount of energy required to propagate the crack per unit surface area. The quasi-static analysis of fracture phenomena can be defined in 3D and 3D. In this analysis, two-dimensional geometry is considered. J-integral is defined as: 6
Eq. 4
where H is: Eq. 5 Abaqus defines the domain in terms of rings of elements surrounding the crack tip. Different “contours” (domains) are created. Each subsequent contour is defined by adding the next ring of elements that share nodes with the elements in the previous contour. q is chosen to have a magnitude of zero at the nodes on the outside of the contour and to be one (in the crack direction) at all nodes inside the contour. Middle nodes have values between 0 and 1 depending on their distance to the crack tip. Them the J-integral is calculated for these paths separately while they should result in the same value according to the theory. Fig. 8 presents the results of J-integral calculations for two α-angles. Each line, as annotated on the figure, represents the J-integral value for a given amount of bending displacement in the center of the wafer under four-point flexural loading. The horizontal lines represent the recommended JIC values from the literature, for three low-index planes. The {111} plane in silicon has the highest atomic packing factor. Therefore, the required energy for generating new surfaces by planar cracking is lower than for the other two planes. The intersections of Jintegral and JIC lines show the maximum tolerable crack length prior to brittle fracturing, under the imposed amount of bending displacement. For example, if the wafer has a crack along α = 0° located in the center of the wafer, and the wafer is bent down 5 mm, the Jintegral calculation gives a maximum tolerable crack length of about 20 µm. If the α-angle is 60° and all other conditions are equal (right-hand chart of Fig 8), then the maximum allowable crack length is around 30 µm. This means that should brittle fracture happen on a {111} plane, the crack orientation of 60° would be more resistant than 0°. Fig. 9 summarizes the effect of the α-angle on the critical crack length. It includes low, medium and high levels of bending displacement (5.0 mm, 10.0 mm and 14.5mm, respectively). The variation of critical crack length for different crack orientations is more pronounced for low of bending displacements. The critical crack length tends to flatten out at high bending displacements, meaning that the stress concentration at the crack tip will already be higher than the minimum required energy for making new surfaces on {111} planes for all α-angle values. In other words, the effect of crack orientation is more significant in brittle fracture when the bending displacement is small. For larger bending displacements, the brittle fracture is almost equally critical for all crack orientations. At low bending displacement values, the critical crack length is the longest if the crack is oriented 60° and shortest if the α-angle is 20°. This relationship remains unchanged for medium and high bending displacement, but the variation is less conspicuous.
7
Fig. 8: J-integral calculation for α = 0° and 60°.
Fig. 9: Benchmarking of the α-angle's effect on the maximum allowable crack length.
3.2. LEFM approach Fig. 10 shows the stress intensity factors for two α angle values. The plots are related to Fig. 8. However, here the lines represent the crack length rather than the bending displacement, and the x-axis represents the wafer bending displacement. The recommended KIC values are extracted from [16, 20]. The intersection of KIC and SIF lines for each crack length will give the maximum allowable bending displacement in the center of the wafer.
8
Fig. 10: Stress intensity factor (SIF) calculation for α = 0° and 60°.
Fig. 11a shows how the critical bending displacement decreases as the crack length increases. As demonstrated in the figure, the crack orientation becomes less important when the crack is large. This is in agreement with findings from the J-integral analysis. Fig. 11b depicts the benchmarking results of SIF calculations for different crack orientations. The SIF analysis, however, suggests different results compared to the J-integral analysis that is shown in Fig. 9. According to Fig. 11b, if the α-angle is 30°, the wafer can tolerate the highest bending displacement, whereas if α is 0° or 90°, the maximum allowable bending displacement is the lowest. Further, the maximum allowable bending for the same crack length does not correlate between these two analyses. For example, let's say the maximum allowable bending is assumed to be 10.0 mm. At 0°, according to Fig. 9, the critical crack length is just below 5 µm, while in Fig. 11b, it is 10 µm. The SIF calculation in this case overestimates the allowable crack size by a factor of 2, compared to the value suggested by the J-integral approach.
Fig. 11: (a) Critical bending displacement and (b) α-angle benchmarking for different crack lengths, based on SIF calculations.
Hence, there are two major differences between J-integral and SIF approaches: 9
Difference in the weakest and strongest fracture crack orientation. Difference in the maximum allowable crack lengths and bending displacements.
In order to argue the probable reasons for these differences, a wafer bending experimental setup was utilized to measure the maximum bending displacement values. Fig. 12a shows the results of the designed experiment. The bending tests were performed along two perpendicular directions: 100 and 010, which are along and perpendicular to the saw marks respectively. The bending setup was a simple manual machine that was pushing the plastic anvils using a manually controlled screw. The screw was rotated at a similar low rate for all tests. Once the wafer shattered, the bending was ceased and the total displacement was measured with reference to the initial state of the screw before bending. Each presented curve is constructed from the fracture data of one hundred wafers bent under the same conditions. The maximum bending displacement at 50 % of cumulative percentage is about 8 mm for the weak direction. This value increases to about 22 mm when bending perpendicular to the sawing direction. The weakest direction can tolerate no more than 11-12 mm of bending displacement before all wafers will break. Fig. 12b shows an optical image from a tapered cross section of a wafer. The maximum crack depths after sawing are observed to be around 4 µm. The combination of experimental maximum crack depth and maximum bending displacement is in agreement with the J-integral calculations that were presented earlier. According to Fig. 9, if the wafer can tolerate 10 mm of bending displacement, the maximum crack length to initiate fracture along {111} is just below 5 µm if the crack is oriented mainly along the saw marks (α-angle < 30°).
Fig. 12: (a) Experimental measurement of bending displacement along and perpendicular to the saw marks and (b) Optical microscope image of wafer cross section showing the damage profile.
The applicability of LEFM and EPFM was explored in the development of fracture mechanics methods. For measuring and calculating the KIC under LEFM conditions, the size of the sample with respect to the crack tip plastic region plays an important role, thus, the sample geometry should be used for defining the critical SIF [21]. The suggested sizes are standardized and being practiced in the fracture mechanics experiments and simulations. The assumptions of linear elastic fracture mechanics may not hold if the plastic zone at the crack 10
tip has the size of the same order of magnitude as the crack size [22]. Furthermore, the plastic zone may change as the applied load and the crack size are increased. In other words, the LEFM approach is not suitable for short crack lengths and thin sections, both of which apply to the current study of silicon single crystal wafers. The performed simulations are excluding the involvement of dislocation role in the fracture mechanisms. The generation of interstitial dislocation loops (Franks Partials) and shear faults are the proposed mechanisms that may assist fracture initiation at atomic level [23]. Samuels and Roberts [24] have shown that silicon single crystal show some dislocation activity above ductile to brittle transition temperature, which is about 10°C depending on the type of doping. The dislocation mobility is a function of the local stress intensity that are studied and benchmarked in this paper. The Mode I fracture toughness is correlated with the J-Integral value under monotonic loading conditions; however, this relation is only for isotropic and perfectly brittle linear elastic materials under plane strain deformation [25]. Therefore, it can be inferred that the J-integral calculations are applicable for linear elastic, nonlinear elastic and elastic-plastic materials, unlike stress intensity factor calculation, which is applicable for cracks in homogeneous, linear elastic materials [26]. In order to calculate the stress intensity factor, the crack propagation direction must be defined while it is valid only for homogeneous, isotropic linear elastic materials. Calculation of J-integral is insensitive to crack propagation path as it is a path independent method. This explains the discrepancy between the results of SIF and J-integral calculations. 4. Conclusions Based on simulation and experimental results, the following conclusions can be deduced.
The elastic-plastic fracture mechanics (EPFM) approach shows consistent results when compared to the experimental test outcome. It is suggested as a suitable tool for analyzing the effect of sawing induced cracks on the structural integrity of silicon single crystal wafers. Based on the J-integral calculations, the most critical crack orientation with respect to the 100 is around 20°. If the sawing process is optimized to avoid generation of subsurface micro-cracks that are aligned close to 20°, the wafers should tolerate more bending displacement. On the other hand, if the crack orientation is at 50-60°, the brittle fracture along {111} planes becomes less susceptible, and the wafer can tolerate the highest bending displacement. The calculation of SIF using the LEFM approach is not a suitable option for analyzing the fracture properties of thin silicon wafers. It will overestimate the maximum allowable bending displacement and maximum critical crack sizes. EPFM simulation has produced more reliable results in this regard.
11
5. Further Studies The exact morphology of the cracks and its relation to the sawing parameters are not yet well understood. We have observed that the cracking initiates at different sites and grow and coalesce in milliseconds, which lead to the brittle shattering of the loaded wafer. More experimental and simulation work is required to fully understand the fracture mechanisms in silicon single crystal wafers. This may eventually lead to optimized sawing parameters that can help increase the structural integrity and fracture strength of the wafers. 6. Acknowledgment The authors would like to gratefully acknowledge the Norwegian Research Council for cofunding this research work under the contract number of 219395/O30. 7. References [1] Hardin CW, Qu J, Shih AJ. Fixed abrasive diamond wire saw slicing of single-crystal silicon carbide wafers. Mater Manuf Process. 2004;19:355-67. [2] Wu H, Melkote SN, Danyluk S. Mechanical strength of silicon wafers cut by loose abrasive slurry and fixed abrasive diamond wire sawing. Adv Eng Mater. 2012;14:342-8. [3] Brun XF, Melkote SN. Analysis of stresses and breakage of crystalline silicon wafers during handling and transport. Sol Energ Mat Sol C. 2009;93:1238-47. [4] Cook R. Strength and sharp contact fracture of silicon. J Mater Sci. 2006;41:841-72. [5] Ferrazza F. Crystalline silicon: Manufacture and properties. Practical Handbook of Photovoltaics: Fundamentals and Applications. 2012:79. [6] Coletti G, Tool C, Geerligs L. Mechanical strength of silicon wafers and its modelling. Proc of Workshop on Crystalline Silicon Solar Cells and Modules: Materials and Processes: Citeseer; 2005. p. 117-20. [7] Coletti G, Tool C, Geerligs L. Quantifying surface damage by measuring mechanical strength of silicon wafers. Presented at the 20th European Photovoltaic Solar Energy Conference and Exhibition2005. p. 10. [8] Möller HJ. Basic Mechanisms and Models of Multi‐Wire Sawing. Adv Eng Mater. 2004;6:501-13. [9] Cramer T, Wanner A, Gumbsch P. Energy dissipation and path instabilities in dynamic fracture of silicon single crystals. Phys Rev Lett. 2000;85:788. [10] Yang C, Wu H, Melkote S, Danyluk S. Comparative analysis of fracture strength of slurry and diamond wire sawn multicrystalline silicon solar wafers. Adv Eng Mater. 2013;15:358-65. [11] Zhu X-K, Joyce JA. Review of fracture toughness (G, K, J, CTOD, CTOA) testing and standardization. Eng Fract Mech. 2012;85:1-46. [12] Ponton CB, Rawlings RD. Vickers indentation fracture toughness test Part 1 Review of literature and formulation of standardised indentation toughness equations. Mater Sci TechLond. 1989;5:865-72. [13] Quinn GD, Bradt RC. On the Vickers indentation fracture toughness test. J Am Ceram Soc. 2007;90:673-80. [14] Pérez R, Gumbsch P. An ab initio study of the cleavage anisotropy in silicon. Acta Mater. 2000;48:4517-30. [15] Ding Z, Zhou S, Zhao Y. Hardness and fracture toughness of brittle materials: a density functional theory study. Phys Rev B. 2004;70:184117. 12
[16] Masolin A, Bouchard P-O, Martini R, Bernacki M. Thermo-mechanical and fracture properties in single-crystal silicon. J Mater Sci. 2013;48:979-88. [17] Standard A. C1161, 2002c,“Standard Test Method for Flexural Strength of Advanced Ceramics at Ambient Temperature”, ASTM International, West Conshohocken, PA, 2002, DOI: 10.1520/C1161-02C. [18] Popovich V, Riemslag A, Janssen M, Bennett I, Richardson I. Characterization of Multicrystalline Silicon Solar Wafers Fracture Strength and Influencing Factors. International Journal of Material Science. 2013. [19] Hopcroft MA, Nix WD, Kenny TW. What is the Young's Modulus of Silicon? Microelectromechanical Systems, Journal of. 2010;19:229-38. [20] Sen D, Thaulow C, Schieffer SV, Cohen A, Buehler MJ. Atomistic study of crack-tip cleavage to dislocation emission transition in silicon single crystals. Phys Rev Lett. 2010;104:235502. [21] Anderson TL, Anderson T. Fracture mechanics: fundamentals and applications: CRC press; 2005. [22] Miller K. The short crack problem. Fatigue Fract Eng M. 1982;5:223-32. [23] Ebrahimi F, Kalwani L. Fracture anisotropy in silicon single crystal. Materials Science and Engineering: A. 1999;268:116-26. [24] Samuels J, Roberts S. The brittle-ductile transition in silicon. I. Experiments. Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences: The Royal Society; 1989. p. 1-23. [25] Bower AF. Applied mechanics of solids, Chapter 9: CRC press; 2009. [26] Parks D. The virtual crack extension method for nonlinear material behavior. Computer Methods in Applied Mechanics and Engineering. 1977;12:353-64.
13
List of Nomenclature KC KIC
JC JIC ν E G α β x, y, z x', y', z' Γ q n W I σ u
Calculated stress intensity factor Mode I critical stress intensity factor Calculated J-integral Mode I critical J-integral Poisson ratio Elastic modulus Shear modulus Angle between saw-marks and the crack on the sawing plane Angle between sawing plane normal vector and the crack Original coordination system axes and directions Rotated coordination system axes and directions Path contour around the crack Unit vector in the virtual crack extension direction Outward normal to Γ Elastic strain energy Propagation length Stress level Displacement
14