CN103268064B - Analog calculation method of ablating silicon nitride by ultrashort pulse laser - Google Patents

Abstract

The invention discloses a simulation calculation method for ultrashort pulse laser ablation of silicon nitride, using Matlab to simulate calculation and analysis of processing parameters, the main steps include: establishing an ablation model for ultrashort pulse laser ablation of silicon nitride, The calculation parameters are initialized; the plasma density boundary conditions are defined, and the ablation threshold, depth and volume under different simulation parameters are calculated through the model; based on the calculation results and the residual height model, the laser ablation residual height is analyzed and evaluated, and given out guidance parameters. The present invention can avoid the process of obtaining processing parameters through repeated experiments, and use the simulation results to optimize the processing parameters, thereby shortening the product cycle, reducing processing costs, and improving production efficiency. The processing process has important guiding value.

Description

A Simulation Calculation Method for Ultrashort Pulse Laser Ablation of Silicon Nitride

technical field

The invention belongs to the technical field of ultrashort pulse laser micromachining, in particular to a simulation calculation method for ultrashort pulse laser ablation of silicon nitride.

Background technique

With the rapid development of advanced manufacturing technology and industrial level, research on silicon nitride materials has also attracted much attention. Because silicon nitride materials have the characteristics of high chemical stability, high thermal shock resistance, high hardness, high temperature resistance, radiation resistance, corrosion resistance, good thermohardness, and excellent optical properties. Therefore, it has been widely used in optoelectronics, machinery, atomic energy, aerospace and other industries, and as it is widely used in precision and ultra-precision fields such as microelectronics, it is particularly important to perform microprocessing on it. However, silicon nitride is a high-hard and brittle material, and it is difficult to precisely process it with existing processing methods. Therefore, it is a key problem to be solved to achieve controllable micro-removal of it.

Ultra-short pulse laser micromachining technology is a unique processing feature of femtosecond laser ultra-fast, ultra-strong, and ultra-high. It generates local micro-explosions by generating high-temperature and high-pressure plasma in a very short time, so as to achieve material processing. Microfabrication method for precise removal. Since the femtosecond laser interacts with matter, it has unique advantages such as small area of action, no thermal effect, high processing precision, and can break through the diffraction limit, which makes it possible to precisely process silicon nitride. This not only overcomes the difficulties in processing caused by its high hardness characteristics, but also provides a precision processing method for high-hard and brittle materials, thereby improving the processing quality. However, since the removal accuracy of femtosecond laser ablation of silicon nitride is not easy to control, and it is very cumbersome to find the optimal process parameters through trial and error alone, and the processing efficiency cannot be guaranteed, so there is an urgent need for a method that is both A method that can conveniently predict accuracy without losing processing efficiency. Based on this, the present invention proposes a simulation calculation method for ultrashort pulse laser ablation of silicon nitride, and the key lies in the establishment of the ablation model and the evaluation of simulation calculation results with different parameters. Regarding the theoretical model and related modeling methods of ultrashort pulse laser ablation of dielectric materials, the following literatures have reported:

American scholars: M.D.Perry, B.C.Stuart, P.S.Banks, et al.Ultrashort-pulse lasermachining of dielectric materials[J].J Appl Phys,1999,85(9):6803-6810.

French scholars: L.Sudrie, A.Couairon, M.Franco, et al.Femtosecond laser-induceddamage and filamentary propagation in fused silica[J].Phys Rev Lett,2002,89(18):186601.

American scholars: C.H.Fan, J.Sun, J.P.Longtin. Plasma absorption of femtosecond laser pulses in dielectrics[J].J Heat Transfer,2002,124(2):275-283.

Chinese scholars: Liu Qing, Cheng Guanghua, Wang Yishan, etc. Three-dimensional optical storage and mechanism of femtosecond pulses in transparent materials[J]. Acta Photonica Sinica, 2003,32(3):276-279.

Chinese scholars: Li Xiaoxi, Jia Tianqing, Feng Donghai, etc. Ablation mechanism and ultrafast kinetics of lithium fluoride under ultrashort pulse irradiation[J]. Acta Optics Sinica, 2005,25(11):1526-1530.

Through literature research and analysis, scholars at home and abroad have proposed many theoretical models describing the process of avalanche ionization and photoionization when studying ultrashort pulse laser-induced dielectric materials, and have studied the evolution law of free electron density in the conduction band and the criticality of material damage. Plasma density is described. However, most of these studies focus on the discussion of a single model and a single factor, and seldom carry out systematic multi-model comprehensive analysis; and studies on ablation threshold and ablation depth rarely fully involve the control parameters in the material removal process. The research on ablation volume focuses more on experimental analysis. Domestic scholars' research in this area is mainly aimed at the discussion of the interaction mechanism, and a large number of experimental studies have been carried out. There are still defects in theoretical modeling, and the research objects are mostly quartz glass and other materials. The simulation calculation of silicon nitride materials Analysis is relatively rare. Due to its important application background and the shortcomings of current simulation analysis, it is particularly important to find a method that can not only accurately characterize the interaction process between ultrashort pulse laser and material, but also provide guidance parameters for the actual removal process.

Contents of the invention

In order to solve the problems that the removal accuracy is difficult to control and the process parameters are difficult to determine in the process of ultrashort pulse laser ablation of silicon nitride, the present invention overcomes the shortcomings of the existing simulation calculation method, and studies the single effect of a single variable and multiple The law of joint action of two variables, the key lies in:

1. When preprocessing the avalanche ionization coefficient and photoionization coefficient, the relationship between the ionization rate and the electric field strength should be determined first, and then the relationship between the ionization rate and the energy density should be determined, so as to facilitate the subsequent simulation calculation of the ablation threshold.

2. The spatial energy distribution of laser pulses is described by Gaussian functions, and Gaussian beams can break through the processing diffraction limit, thereby improving processing accuracy.

3. The conditions when the conduction band electron density reaches the critical density should be accurately controlled to ensure that the ablation depth, volume and residual height are calculated under the premise of damage to the material.

The purpose of the present invention is to provide a simulation calculation method for ultrashort pulse laser ablation of silicon nitride, which can not only realize the control of ablation accuracy and the optimization of process parameters, but also avoid blindness caused by repeated trials and reduce processing costs ,Increase productivity.

To achieve the above object, the present invention adopts the following technical solutions:

A. Establish the ablation threshold, depth and ablation volume model of ultrashort pulse laser ablation of silicon nitride, and initialize the parameters of the model constants; at the same time, define the plasma density boundary when ablation of silicon nitride occurs in the built model The condition ρ _cr is 1.6×10 ²¹ cm ^-3 , and the avalanche ionization coefficient, photoionization coefficient and conduction band free electron density of silicon nitride are pretreated;

B. Determine the laser wavelength, set the pulse width calculation range to 10fs-10ps, and calculate the silicon nitride ablation threshold. When the electron density in the conduction band exceeds the critical plasma density, that is, when ρ _c (x,t) > ρ _cr , silicon nitride will be damaged. When ρ _c (x, t) ≤ ρ _cr , the initial variable must be reassigned calculate;

C. Set the value range of the laser energy density to 1J/cm ² ~ 8J/cm ² , and simulate the ablation depth based on the threshold result; under the assumption that the single-pulse ablation morphology is a conical volume modeling assumption, and Define the radius of the laser beam waist, and then use the simulation results of the ablation depth to obtain the ablation volume under the corresponding conditions; the volume model is

$\mathrm{<span\; class="notranslate">\; V</span>}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; \&CenterDot;</span>{{\mathrm{<span\; class="notranslate">\; \&pi;\&omega;</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; 6</span>}}\mathrm{<span\; class="notranslate">\; ln</span>}<span\; class="notranslate">\; (</span>\frac{\mathrm{<span\; class="notranslate">\; f</span>}}{{\mathrm{<span\; class="notranslate">\; f</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}}<span\; class="notranslate">\; )</span>$

In the formula: V is the ablation volume, x is the ablation depth, ω ₀ is the beam waist radius, F _th is the ablation threshold, F is the energy density;

D. Establish the relational model of scanning velocity v, line overlap rate δ and residual height Δx, define the boundary condition of the ablation residual height, and load the threshold and depth model calculation results into the residual height model; the scanning velocity v and the residual height The relationship model of height Δx is

$\mathrm{<span\; class="notranslate">\; v</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; \&Delta;\; L</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; \&xi;D</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; \&Delta;x</span>}}{\mathrm{<span\; class="notranslate">\; x</span>}}<span\; class="notranslate">\; \&times;</span>\mathrm{<span\; class="notranslate">\; 100</span>}<span\; class="notranslate">\; \%</span><span\; class="notranslate">\; \&Center\; Dot;</span>\mathrm{<span\; class="notranslate">\; D.</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; f</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>$

In the formula: ΔL is the distance between two pulses in the direction of scanning velocity, ξ is the relative residual height, D is the ablation diameter, f is the pulse frequency, x is the ablation depth; Δx is the residual height, which refers to the ablation profile peak The distance between the line and the bottom line of the contour valley; the relationship model between the line overlap rate δ and the residual height Δx is

$\mathrm{<span\; class="notranslate">\; \&delta;</span>}<span\; class="notranslate">\; =</span>\frac{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; d</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; \&Delta;d</span>}<span\; class="notranslate">\; )</span>}{\mathrm{<span\; class="notranslate">\; d</span>}}<span\; class="notranslate">\; \&times;</span>\mathrm{<span\; class="notranslate">\; 100</span>}<span\; class="notranslate">\; \%</span><span\; class="notranslate">\; =</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; \&xi;</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; \&times;</span>\mathrm{<span\; class="notranslate">\; 100</span>}<span\; class="notranslate">\; \%</span>$

In the formula: Δd is the ablation line spacing, Δd=ξd; d is the ablation line width, d=D;

E. Assign a value to the pulse frequency, and calculate the residual height under the condition that the energy density is 4.0J/cm ² -8.0J/cm ² and the scanning speed is 0-3mm/s. If Δx>Δx _max does not meet the requirements, the variable calculation should be re-input. If Δx≤Δx _max , the processing accuracy is satisfied, and the guidance parameters are output to complete the simulation calculation.

The avalanche ionization coefficient is:

$\mathrm{<span\; class="notranslate">\; \&eta;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; E.</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; v</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; \&CenterDot;</span>\mathrm{<span\; class="notranslate">\; \&alpha;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; E.</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; v</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}\mathrm{<span\; class="notranslate">\; E</span>}}{\mathrm{<span\; class="notranslate">\; \&Delta;</span>}}\mathrm{<span\; class="notranslate">\; exp</span>}<span\; class="notranslate">\; (</span><span\; class="notranslate">\; -</span>\frac{{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}{\mathrm{<span\; class="notranslate">\; E.</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>\frac{\mathrm{<span\; class="notranslate">\; E.</span>}}{{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; P</span>}}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; kT</span>}}}<span\; class="notranslate">\; )</span>$

Where: v _s is the saturation drift rate, α(E) is the Townsend coefficient, e is the electronic charge, and Δ is the forbidden band width. E _i , E _p and E _kT are electric field strengths required for carriers to overcome ionization, optical phonon and thermal scattering effects, respectively.

Described photoionization coefficient is:

In the formula: ω is the laser frequency, m=m _e m _h /( _me + m _h ) is the reduced effective mass of electron-hole pairs, _me and m _h are the effective masses of electrons and holes respectively; n=<Δ/(hω)+1> indicates the number of photons to be absorbed to excite an electron from the valence band to the conduction band; U=Δ-e ² E ² /(4mω ² ), Φ(z) is the Dawson integral.

The described ablation threshold model is:

${\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; cr</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 0</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; \&tau;</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; [</span>{\mathrm{<span\; class="notranslate">\; W</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; h</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; ]</span>\mathrm{<span\; class="notranslate">\; exp</span>}<span\; class="notranslate">\; [</span><span\; class="notranslate">\; (</span>\frac{{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; W</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}}{{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; h</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; \&tau;</span>}<span\; class="notranslate">\; ]</span><span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; W</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>}{{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; W</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; h</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}\mathrm{<span\; class="notranslate">\; exp</span>}<span\; class="notranslate">\; [</span><span\; class="notranslate">\; (</span>\frac{{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; W</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}}{{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; h</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; \&tau;</span>}<span\; class="notranslate">\; ]</span><span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; h</span>}}_{\mathrm{<span\; class="notranslate">\; the\; th</span>}}<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}\mathrm{<span\; class="notranslate">\; 0</span>}}<span\; class="notranslate">\; )</span>}$

In the formula: ρ _a0 , ρ _c0 are the initial valence band, conduction band electron density, respectively, ρ _cr is the critical electron density. H _th =η(F _th ); W _th =w _PI (F _th ), F _th is the ablation threshold, and τ is the pulse width.

The described ablation depth model is:

In the formula: x is the ablation depth, h is Planck’s constant; N=(<ε>+Δ)/hω, <ε> is the average energy of electrons in the plasma, W _F =[F/F _th ] ⁿ W _PI (F _th ), _HF = [F/F _th ] ^1/2 η(F _th ).

The described conduction band free electron density equation is:

$\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; =</span>{<span\; class="notranslate">\; (</span>\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; )</span>}_{\mathrm{<span\; class="notranslate">\; AI</span>}}<span\; class="notranslate">\; +</span>{<span\; class="notranslate">\; (</span>\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; )</span>}_{\mathrm{<span\; class="notranslate">\; P.I.</span>}}<span\; class="notranslate">\; -</span>{<span\; class="notranslate">\; (</span>\frac{{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; \&rho;</span>}}_{\mathrm{<span\; class="notranslate">\; c</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; t</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&PartialD;</span>\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; )</span>}_{\mathrm{<span\; class="notranslate">\; lost</span>}<span\; class="notranslate">\; .</span>}$

In the formula: the AI term and the PI term represent the change of conduction band electron density caused by avalanche ionization and photoionization, respectively, and the Los. term represents the loss of carriers.

The avalanche ionization coefficient and photoionization coefficient pretreatment refers to the relationship between light intensity, electric field intensity and energy density, and finally establishes the function between ionization coefficient and energy density and calculates and analyzes it; conduction band free electron density pretreatment It refers to the use of its characterization equations, combined with the spatiotemporal expression of the laser for further derivation, and the simulation calculation is carried out based on the initial conduction band, valence band electron density and other parameters of silicon nitride.

The relation of described light intensity and electric field intensity, energy density is:

$\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; no</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}{\mathrm{<span\; class="notranslate">\; E.</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; R</span>}}_{\mathrm{<span\; class="notranslate">\; 0</span>}}}<span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; f</span>}}{\mathrm{<span\; class="notranslate">\; \&tau;</span>}}$

In the formula: R ₀ =cμ ₀ is the vacuum impedance, c is the speed of light in vacuum, μ ₀ is the vacuum permeability, and n ₀ is the refractive index of the medium.

Compared with the prior art, the present invention has the following beneficial effects:

1. Compared with the schemes reported in the past, the simulation calculation method used in the present invention fully considers the electron ionization mode and electron density evolution law when the ultrashort pulse laser interacts with silicon nitride, and can more truly reflect the combustion process. The mechanism of material removal during etching;

2, the simulated calculation method adopted in the present invention is compared with the scheme reported in the past, and each processing parameter is determined conveniently and is easy to realize in actual production, and corrects after comparing with simulation result, can be included in process parameter database, for follow-up process Optimization has important guiding value;

3. Compared with the schemes reported in the past, the simulation calculation method adopted in the present invention can be flexibly used in the processing of products with different precisions and sizes, which reduces design and development costs while improving efficiency and making it easier to ensure processing quality.

Description of drawings

The present invention has 6 accompanying drawings, wherein:

Fig. 1 is a flow chart of a simulation calculation method for ultrashort pulse laser ablation of silicon nitride;

Fig. 2 is the change law of silicon nitride ablation threshold with pulse width;

Fig. 3 is the change rule of silicon nitride ablation depth with energy density;

Fig. 4 is the change law of silicon nitride ablation volume with energy density;

Figure 5(a) and (b) are schematic diagrams of the influence of scanning speed and line overlap ratio on residual height, respectively;

Fig. 6 is the change law of residual height with scanning speed, energy density and pulse width.

Detailed ways

Below in conjunction with accompanying drawing 1-6 and embodiment, the present invention is further described: as shown in Figure 1, a kind of flow chart of the simulation calculation method of ultrashort pulse laser ablation silicon nitride, the embodiment of the present invention requires: at laser wavelength Under the condition of 780nm, the process parameters used to guide the actual processing are obtained through simulation calculations. It mainly includes damage threshold under the condition of pulse width of 10fs～10ps; energy density of 1J/cm ² ～8J/cm ² , ablation depth and volume and scanning speed of 0～3mm/ s, residual height when the energy density is 4.0J/cm ² -8.0J/cm ² . The specific simulation calculation steps are as follows:

A. Establish the ablation threshold, depth and ablation volume model of ultrashort pulse laser ablation of silicon nitride, and initialize the parameters of the model constants. Through Matlab software, when the electric field strength E=0～350MV/cm, and the forbidden band width Δ=5.0eV, the Townsend coefficient can be calculated; in addition, the free electron saturation drift rate v _s =2×10 ⁷ cm/s, the vacuum impedance R ₀ =376.991Ω, refractive index n0=2.0, the relationship between avalanche ionization rate and energy density can be calculated; for silicon nitride materials, the number of photons to be absorbed by exciting an electron from the valence band to the conduction band is n=4, and then in When the reduced effective mass of electron-hole pairs is m=0.86m _e =7.833×10 ^-31 kg, the photoionization rate is calculated by simulation.

B. Use the conduction band free electron density evolution equation ρ _c (x,t) to describe the electron density during the ablation process of silicon nitride, where: the initial valence band electron density ρ _a0 =1.10×10 ²³ cm ^-3 , the initial Conduction band electron density ρ _c0 =10 ¹⁰ cm ^-3 , critical electron density ρ _cr =1.6×10 ²¹ cm ^-3 . When ρ _c (x, t) > ρ _cr , the ablation threshold condition is satisfied, and then the short pulse ablation threshold of silicon nitride is calculated under the condition of pulse width of 10fs ~ 10ps, the results are shown in Fig. 2, where When τ>10ps, it is the long pulse ablation threshold.

C. Based on the threshold simulation results, the damage thresholds obtained when the pulse width is 12fs, 35fs and 220fs are respectively 2.0J/cm ² , 1.8J/cm ² and 1.4J/cm ² ; the value range of the energy density is set to 1J /cm ² ～8J/cm ² , and calculate the ablation depth, the results are shown in Figure 3. Furthermore, when the beam waist radius ω ₀ =30 μm, the ablation volume at different pulse widths is analyzed through the volume model, as shown in Fig. 4; the ablation depth curve obtained in Fig. 3 can be used to select the ablation depth at different depths The energy density parameter of .

D. Establish the relationship model of scanning velocity v, line overlap rate δ and residual height Δx respectively. The modeling diagrams are shown in Figure 5(a) and (b); at the same time, define the boundary conditions of the ablation residual height, and use Figure 2 The obtained ablation threshold results and the ablation depth results obtained in Fig. 3 were loaded into the residual height model.

E. Set the pulse frequency to 1kHz, the energy density to 4.0J/cm ² to 8.0J/cm ² , and the scanning speed to 0 to 3mm/s, and use the scanning speed model to simulate and analyze the residual height. The results are shown in Figure 6 Show. According to the requirements of residual height, the process parameters such as scanning speed, energy density and pulse width used in processing can be determined. If Δx>Δx _max does not meet the requirements, the variable calculation should be re-entered. If Δx≤Δx _max , the processing accuracy is met. Output guidance parameters to complete the simulation calculation.

The above-mentioned embodiments are only used to illustrate rather than limit the method scheme of the present invention, and any equivalent replacement or change according to the technical scheme of the present invention and its inventive concept shall fall within the protection scope of the present invention.

Claims (1)

1. a simulation method for mechanism of ultrashort-pulse laser ablation silicon nitride, is characterized in that the method comprises the following steps:

A, set up the ablating model of mechanism of ultrashort-pulse laser ablation silicon nitride, and parameter initialization is carried out to model constants, define plasma density boundary condition ρ during silicon nitride generation ablation in institute's established model simultaneously _crbe 1.6 × 10 ²¹cm ^-3;

B, determine optical maser wavelength, setting pulsewidth computer capacity is 10fs ~ 10ps, calculates silicon nitride ablation threshold; When conduction band electron density exceedes critical plasma density, i.e. ρ _c(x, t) > ρ _crtime, silicon nitride produces damage, works as ρ _c(x, t)≤ρ _crtime, assignment must be carried out again to initializaing variable and calculate;

C, setting laser energy density codomain are 1J/cm ²~ 8J/cm ², based on threshold value result, ablation depth is simulated; Under the modeling volume assumed condition that single pulse ablation pattern is cone shape, set up ablation volume-based model, and define laser beam waist radius, and then the ablation volume under utilizing ablation depth analog result to obtain corresponding conditions; Described volume-based model is

$V=\frac{x\&CenterDot;\mathrm{\&pi;}{{\mathrm{\&omega;}}_0}^2}{6}\ln \left(\frac{F}{F_{\mathrm{th}}}\right)$

In formula: V is ablation volume, x is ablation depth, ω ₀for waist radius, F _thfor ablation threshold, F is energy density;

D, set up the relational model of sweep speed v, line overlap rate δ and residual altitude Δ x, definition ablation residual altitude boundary condition, and threshold value and depth model result of calculation are loaded in residual altitude model; The relational model of described sweep speed v and residual altitude Δ x is

$v=\mathrm{\&Delta;L}(f-1)=\mathrm{\&xi;D}(f-1)=\frac{\mathrm{\&Delta;x}}{x}\&times;100\%\&CenterDot;D(f-1)$

In formula: Δ L is two pulse distances on sweep speed direction, and ξ is relative residual altitude, and D is ablation diameter, and f is pulse frequency, and x is ablation depth; Δ x is residual altitude, refers to the distance between ablation profile summit line and profile valley line; The relational model of described line overlap rate δ and residual altitude Δ x is

$\mathrm{\&delta;}=\frac{(d-\mathrm{\&Delta;d})}{d}\&times;100\%=(1-\mathrm{\&xi;})\&times;100\%$

In formula: Δ d is ablation distance between centers of tracks, Δ d=ξ d; D is ablation line width, and ξ is relative residual altitude, d=D;

E, paired pulses frequency carry out assignment, and are 4.0J/cm in energy density ²~ 8.0J/cm ², sweep speed is calculate residual altitude under 0 ~ 3mm/s condition; If Δ x > Δ x _max, do not meet the demands, re-enter variable and calculate; If Δ x≤Δ x _max, meet machining accuracy, export guide parameters, complete analog computation.