High-precision lithography thick-mask model based on a decomposition machine learning method

Ziqi Li; Ziqi Li; Lisong Dong; Lisong Dong; Lisong Dong; Lisong Dong; Xuyu Jing; Xuyu Jing; Xu Ma; Yayi Wei; Yayi Wei; Yayi Wei; Yayi Wei; Yayi Wei

doi:10.1364/OE.454513

1. Introduction

Lithography is a key technique to manufacture the very large-scale integration circuit (VLSI), and it is the only process to replicate the pattern of circuit layout from the mask to the wafer [1,2]. Deep ultraviolet (DUV) lithography is a mature technique for the high-volume-manufacturing of advanced integration circuits [3]. The simplified projection structure of DUV lithography system is shown in Fig. 1. The DUV light rays emitted from the illumination system are diffracted by the mask pattern and generate the diffraction near-field (DNF) underneath the mask. Then, the light rays propagate through the projection system and form the interference image on the photoresist, where the chemical reaction is initiated by the photons so that the mask pattern is transferred to the resist stack.

Fig. 1. The simplified projection structure of DUV lithography system

Download Full Size | PDF

As the density of semiconductor devices continuously increases, the wavelength of DUV light becomes larger than the critical dimension (CD) of the integration circuits at advanced technology nodes. In order to improve of the image fidelity on the wafer, resolution enhancement techniques, such as optical proximity correction (OPC) [1,4,5] and source-mask optimization (SMO) [6–9], are utilized in lithography process. Computing the DNF of photomask is a fundamental step in lithography image simulation, and it is indispensable to the OPC and SMO workflows [10]. Thus, it is important to retain the accuracy and speed of the DNF calculation.

However, these two goals are contradictory to each other in general. Kirchhoff approximation model, also called the thin mask assumption, is widely used in the mask DNF simulation due to its low time-consumption property. The method employs scalar diffraction theory and considers the mask an infinitely thin layer. Therefore, the diffraction effect of the three-dimensional (3D) mask structure, also called thick-mask effect, is ignored [11]. The DNF under Kirchhoff approximation is expressed as:

(1)$$\textrm{t}(x,y) = \left\{ {\begin{array}{lc} 1&{x,\textrm{ }y \in \textrm{transparent area}}\\ 0&{x,y \in \textrm{ otherwise}} \end{array}} \right.,$$

where $\textrm{t}(x,y)$ stands for the mask transmissive function, and the amplitude of DNF is assumed to distribute evenly in the transparent area.

At the 45 nm-7 nm technology nodes and beyond, the scattering effect of photomask becomes more and more pronounced. The thin-mask assumption cannot predict accurate results anymore. Considering the polarization of the light wave, the DNF can be represented by four diffraction matrices denoted as ${\mathbf E}(UV)$, where U = X or Y and V = X or Y. The diffraction matrix is defined as the complex amplitude of the DNF polarized in U direction, which is generated by a unit incident electric field polarized in the V direction. The diffraction matrices of thick-mask can be calculated by the rigorous electromagnetic field (EMF) solvers [12–16], including the finite-difference time-domain (FDTD) method [17,18] and rigorous coupled wave analysis (RCWA) method [19]. However, those methods are so computationally expensive that they can be hardly used in the large-scale simulation or inverse mask optimization approaches, for the unacceptable time consumption.

In order to balance the computational cost and simulation accuracy, many approximate thick-mask models were proposed in previous literatures, such as boundary-layer (BL) model [20–22], domain decomposition method (DDM) [23], filter-based methods [24–29], machine learning methods [3,30], deep-learning methods [31–33] and analytical models [34–36]. Among them, the machine learning and deep-learning methods are well studied in recent years. In Ref. [33], authors proposed a deep-learning method for aerial image prediction in EUV lithography. By constructing a neural network with the training data, the mothed can predict the aerial image and reduce the simulation error significantly. On the other hand, the filter-based mask 3D model had been widely used, which can rapidly predict the DNF by using the calibrated convolution kernels. In this kind of methods, the mask pattern is first decomposed into multiple segments, and the DNF of each mask segment is calculated by convoluting with the pre-calculated kernel. Then, the DNF segments are synthesized together to generate the full DNF image [24,27,28]. According to the property of convolution, these methods assume that the diffraction behavior within each segment is shift-invariant, and thus the DNF should have smooth distribution.

However, irradiated by the incident light, the boundaries of mask features generate diffractive waves to the external area. The diffractive wave further interferes with the incident wave and produces fluctuation in the DNF. The phenomenon is named as edge interference effect, which will be further discussed in Sec. 2. As the CD of mask features decreases, the edge interference effect becomes more pronounced. Figures 2(a) and 2(b) respectively show a mask layout and the amplitude of its DNF under the X-polarization illumination, where Fig. 2(b) is obtained by the rigorous EMF simulator. Figure 2(c) provides the amplitude of DNF calculated by the traditional filter-based method. It shows that the filter-based method based on the shift-invariant assumption of diffraction behavior cannot render the accurate prediction result, due to the existence of edge interference effect.

Fig. 2. Simulation of thick-mask DNF under the X-polarization illumination: (a) the mask layout; (b) the amplitude of DNF obtained by the rigorous EMF simulator, where the edge interference effect is observed; (c) the amplitude of DNF calculated by the traditional filter-based method.

Download Full Size | PDF

This paper proposes a high-precision thick-mask model for the polarized DNF simulation. In this model, the diffraction transfer matrix (DTM) is introduced to represent the transformation from the mask pattern to the corresponding DNF. The local DNF is calculated by the matrix multiplication between the mask segment and the associated DTM, which is different form the filter-based method. Since the convolution is replaced by the DTM, the shift-invariant assumption of diffraction behavior is no longer needed, and the fluctuation in DNF can be simulated with high accuracy.

In the proposed method, the training library of thick-mask DNF is first established based on some training mask clips. The DNF images of the training mask clips are pre-calculated by the rigorous EMF simulator. Then, the training layouts and their DNFs are decomposed into multiple groups based on the graphic features, including the concave corners, convex corners and edge segments. The DTM associated to each feature group is calculated from the training dataset based on the least square method. Subsequently, the test layout is also decomposed as mentioned above, and the local DNFs of all feature segments are calculated based on the pre-calibrated DTMs. Finally, the entire DNF is synthesized by stitching up all of the DNF segments. In Sec. 4, we also discuss the influence of the dimension of DTM on the simulation accuracy. Compared with the traditional filter-based method, the proposed method can reduce the root-mean-square (RMS) error of the synthesized DNF by about 95%. The calculation speed of the proposed method is 2.3 times slower than the filter-based method, but it still 300-fold faster than the rigorous EMF method.

The rest of this paper is organized as follows. The mathematical deduction of the edge interference effect is provided in Sec. 2. The traditional filter-based model and the proposed thick-mask model are described in Sec. 3. Simulations and discussions are presented in Sec. 4. Conclusions are given in Sec. 5.

2. Edge interference effect

In this section, we discuss the edge interference effect and its mathematical deduction. The simulation results based on the diffraction model and the rigorous EMF simulator are presented.

Figures 3(b) and 3(c) show the DNF amplitudes of a contact hole under the Y-polarization illumination and X-polarization illumination, respectively. The dimension of the contact hole is 300 nm on the wafer scale, and the illumination wavelength is 193 nm. The thick-mask DNF is calculated using the rigorous EMF simulator. It is observed that, when the CD of mask features is smaller than the wavelength, the amplitude of DNF shows obvious fluctuation. Additionally, the direction of the fluctuation is highly correlated with the polarization state. This phenomenon is named “edge interference effect”, and its generation mechanism is discussed next. Notice that, in this section, the bold symbols denote the vector, and the italic symbol with subscripts x, y, and z denotes the component of the vector along the x, y and z axes, respectively.

Fig. 3. Illustration of the edge interference effect. (a) The diagrammatic sketch of diffraction model, and the DNF images calculated by the rigorous EMF simulator under the (b) Y-polarization illumination and (c) X-polarization illumination, respectively. The fluctuations are observed apparently in (b) and (c).

Download Full Size | PDF

The typical mask structure for DUV lithography is composed of a transparent SiO₂ base layer, and an opaque chromium layer that contains the mask pattern. The diffraction model is illustrated in Fig. 3(a). It is worth noting that, the electric field component of the polarized illumination is parallel to the z-axis under this circumstance. Given the high conductivity of the metallic material (chromium) and the boundary condition of the electromagnetic field, it can be inferenced that, on the surface of metallic layer [37]: (1) the tangential component of the electric field E is 0; (2) the tangential component of the magnetic field H is vertical with the current density J, and the amplitude of H is 4π|J|/c, where c is the speed of light in vacuum; (3) the normal component of H is 0; (4) the normal component of E is proportional to the charge density on the surface.

Irradiated by the incident light field, a scatter field is generated on the surface of metallic layer. Therefore, the total electric field E is composed of the incident electric field and the scatter electric field. Suppose that the polarization direction is paralleled with the boundary of the metallic layer, as shown in Fig. 3(a). In the plane diffraction model, the partial derivatives of Maxwell equations in z-axis are zeros. Therefore, it can be inferred that:

(2)$${E_x} = {E_y} = {H_z} = 0,\,{H_x} = \frac{1}{{ik}}\frac{{\delta {E_z}}}{{\delta y}},\,{H_y} = \frac{1}{{ik}}\frac{{\delta {E_z}}}{{\delta x}}, $$

(3)$$\frac{{{\partial ^2}{E_z}}}{{\partial {x^2}}} + \frac{{{\partial ^2}{E_z}}}{{\partial {y^2}}} + {k^2}{E_z} = 0,$$

where ${E_p}$(p = x, y or z) denotes the component of electric field in the corresponding axis, and ${H_P}$ denotes the component of magnetic field in the corresponding axis. Based on Eq. (2) and Eq. (3), the expression of E and H can be rewritten in the plane wave angular spectra:

(4)$${\mathbf E} = (0,0,1) \cdot \exp ( - ikr\cos (\theta - \alpha )),\,{\mathbf H} = (\sin \alpha , - \cos \alpha \textrm{,}0) \cdot \exp ( - ikr\cos (\theta - \alpha )),$$

where k is the wave number, $\theta $ and r respectively denote the angle and distance between the observation point and the origin of coordinate system, and $\alpha $ is the included angle between the propagation direction of the plane wave and the x axis. According to Eq. (4) and the second boundary condition mentioned above, the z-axis component of current intensity J on the surface can be written as the Fourier integral:

(5)$${J_z}(x) ={-} \frac{c}{{2\pi }}\int_{ - \infty }^\infty {P(\mu )\exp (ikx\mu )d\mu },$$

where $\mu = \cos \alpha $, and formula $P(\mu )$ determines the spectral distribution. Let ${{\mathbf E}^i}$ and ${{\mathbf E}^s}$ denote the incident electric field and the scatter electric field, respectively. The component $E_z^s$ generated by the current distribution in Eq. (5) and the incident component $E_z^i$ can be formulated as:

(6)$$E_z^s = \int {P(\cos \alpha )} \exp (ikr\cos (\theta \mp \alpha ))d\alpha = \int_{ - \infty }^\infty {\frac{{P(\mu )}}{{\sqrt {1 - {\mu ^2}} }}} \exp (ikx\mu )d\mu,$$

(7)$$E_z^i = \exp ( - ikr\cos (\theta - {\alpha _0})),$$

where ${\alpha _0}$ stands for the incident angle.

Considering the first boundary condition mentioned above, we have $E_z^s + E_z^i = 0$. Therefore, based on the Jordan’s lemma and the Cauchy’s residue theorem, it can be further deducted from Eqs. (6) and (7) that:

(8)$$P(\cos \alpha ) = \frac{i}{\pi }\frac{{\sin \frac{{{\alpha _0}}}{2}\sin \frac{\alpha }{2}}}{{\cos \alpha + \cos {\alpha _0}}} , $$

The z-axis component of the total electric field can be expressed as:

(9)$${E_z} = E_z^i + E_z^s = \exp ( - ikr\cos (\theta - {\alpha _0})) - \frac{1}{{i\pi }}\int {\frac{{\sin \frac{{{\alpha _0}}}{2}\sin \frac{\alpha }{2}}}{{\cos \alpha + \cos {\alpha _0}}}\exp (ikr\cos (\theta \mp \alpha ))d\alpha } , $$

Introducing the Fresnel integral, the Eq. (9) can be transformed as follows:

(10)$$\begin{array}{l} {E_\textrm{z}} = \frac{{\exp ( - \frac{{i\pi }}{4})}}{{\sqrt \pi }}\left\{ {\exp ( - ikr\cos (\theta - {\alpha_0}))F\left[ { - \sqrt {2kr} \cos \frac{1}{2}(\theta - {\alpha_0})} \right]} \right.\\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \left. { - \exp ( - ikr\cos (\theta + {\alpha_0}))F\left[ { - \sqrt {2kr} \cos \frac{1}{2}(\theta + {\alpha_0})} \right]} \right\}\\ F(a) = \int_a^\infty {\exp (i{\mu ^2})d\mu } \end{array} , $$

Then, the Eq. (10) can be further simplified by introducing the expressions in Eqs. (11) and (12):

(11)$$u ={-} \sqrt {2kr} \cos \frac{1}{2}(\theta - {\alpha _0}),\,v ={-} \sqrt {2kr} \cos \frac{1}{2}(\theta + {\alpha _0}), $$

(12)$$G(a) = \exp ( - i{a^2})F(a),$$

and Eq. (10) can be rewritten into a more compact form:

(13)$${E_z} = \frac{{\exp ( - \frac{{i\pi }}{4})}}{{\sqrt \pi }}\exp (ikr)\{{G(u) - G(v)} \},$$

The edge interference effect can be described as follows. For a given transparent area, the boundaries of the metallic layer are irradiated by the polarized light beams, thus generating the diffractive waves to the external areas. The diffractive wave further interferes with the incidence wave, and produces the amplitude fluctuation in the thick-mask DNF. It is worth noting that, in the deduction process, the metallic layer is assumed to be infinitely thin, ideally. However, the diffraction model can still come out with the fluctuation in DNF, which is identical to the rigorous EMF simulation. On the other hand, the Kirchhoff approximation, also called the thin-mask approximation, is a rather simple model. It ignores the thickness of the mask, as well as the edge interference effect mentioned above.

Figure 4(a) illustrates the amplitude of electric field outward the mask boundary, and Fig. 4(b) shows the cross-section of Fig. 4(a). The result in Fig. 4 is calculated from Eq. (13). The fluctuation caused by the mask boundary is obvious in the amplitude distribution of electric field.

Fig. 4. The amplitude of electric field outward the mask boundary that is calculated by the diffraction model. (a) The 3D amplitude distribution and (b) the cross-section of the amplitude distribution

Download Full Size | PDF

Figures 5(a) and 5(b) illustrate the 3D distribution and cross-section of the electric field amplitude that is calculated by the rigorous EMF simulator. It is worth noting that, the amplitude distribution in Fig. 4 is deducted based on one-dimensional approximation, which assume the pattern is infinite in y dimension and finite in x dimension. On the other hand, the amplitude distribution in Fig. 5 is calculated with the rigorous EMF simulator. The pattern is in 3-dimensional shape and the physical properties of the metallic layer is identical with the actual material. From the same tendency in Fig. 4 and Fig. 5, it shows the edge interference is an important reason that cause DNF fluctuates.

Fig. 5. The amplitude of electric field outward the mask boundary that is calculated by the rigorous EMF simulator. (a) The 3D amplitude distribution and (b) the cross-section of the amplitude distribution

Download Full Size | PDF

3. Traditional filter-based model and the proposed model

Calculating the DNF of thick mask is a fundamental task in the advanced lithography image simulation. In this section, we describe the traditional filter-based model and the proposed model, which is referred to as the DTM model in the following. We also compare these two models in detail.

3.1 Traditional filter-based model

By applying the fast Fourier transform (FFT), the filter-based model can calculate the thick-mask DNF rapidly. Compared with the Kirchhoff approximation model and the EMF method, the filter-based model strikes a balance between the time consumption and simulation accuracy [23,26,27]. The filter-based model assumes the diffraction behavior is shift-invariant. Therefore, the diffraction effect can be expressed by the convolution between the mask pattern and the pre-calibrated filters:

(14)$${\mathbf A} = {\mathbf M} \otimes f,$$

where ${\mathbf A}$ stands for the synthesized amplitude matrix of DNF, ${\mathbf M}$ represents the mask pattern matrix, ${\otimes} $ is the convolution symbol, and f denotes the convolution filter, respectively.

Moreover, the diffraction behaviors of different geometric patterns are dissimilar. Therefore, Eq. (14) can be further improved. In the algorithm, the mask layout pattern is decomposed into segments based on its geometric features. The convolution filters corresponding to different types of features are calibrated independently. Taking the mask 3D method in [26] as an instance, the mask pattern is decomposed into polygon areas, horizontal edges, vertical edges, and corners. The entire DNF is synthesized by stitching up all of the DNF segments that are calculated by the filter-based method. The calculation process is given as follows:

(15)$$\begin{array}{l} {\mathbf A} = {{\mathbf A}_P} + {{\mathbf A}_{HE}} + {{\mathbf A}_{VE}} + {{\mathbf A}_{Cor}}\\ \left\{ {\begin{array}{{c}} {{{\mathbf A}_P} = {{\mathbf M}_P} \otimes {f_P}}\\ {{{\mathbf A}_{HE}} = {{\mathbf M}_{HE}} \otimes {f_{HE}}}\\ {{{\mathbf A}_{VE}} = {{\mathbf M}_{VE}} \otimes {f_{VE}}}\\ {{{\mathbf A}_{Cor}} = {{\mathbf M}_{Cor}} \otimes {f_{Cor}}} \end{array}} \right. \end{array},$$

where the footmark P, HE, VE, and Cor stand for the polygons, horizontal edges, vertical edges and corners, respectively, A, M and f stand for the calculated DNF amplitudes, mask pattern segments, and convolution filters of different graphic features, respectively.

As the CD of mask feature decreases at advanced technology nodes, the edge interference effect becomes more and more pronounced, which causes the intense inhomogeneity in the DNF amplitude. Figure 6 illustrates the significant difference between the DNF amplitudes of the mask openings with large CD and small CD. The wavelength of the illumination is 193 nm, and the polarization direction of the incident light is displayed in the figure. On the left side, the CD of the mask opening is about 1.6µm on the wafer scale. The influence of the edge interference effect mainly concentrates within the boundary areas. The periodical fluctuations are trivial compared to the size of the pattern. The DNF amplitude can be depicted by the blue dashed line approximately. The image on the right side indicates the result as the CD reduces to 150 nm on the wafer scale. The fluctuation effect has significant impact on the entire distribution of DNF amplitude.

Fig. 6. The cross-sections of the mask openings and DNF amplitude distributions with different CD sizes. The CD sizes of the mask openings are 1.6µm (left) and 150 nm (right), respectively.

Download Full Size | PDF

It is observed that, when the CD of mask feature is large enough, the DNF amplitude distributes inside the pattern evenly. Due to the shift-invariant assumption, the contribution to DNF amplitude of each point at different positions is identical. As the filter shifts to another position inside the pattern, the filter-based model can render a similar DNF amplitude. Therefore, the calculated DNF amplitude is even-distributed inside the pattern, which is identical with the ground truth. However, as the CD of mask feature decreases, the DNF amplitude distributes inside the pattern unevenly. The even-distributed result that rendered by the filter-based model is no longer precise, for the contribution to DNF amplitude of each point at different positions are not identical anymore. Therefore, the shift-invariant assumption is no longer valid, and the traditional filter-based model is inadequate to obtain accurate result.

3.2 Proposed DTM model

To calculate the DNF of thick-mask precisely, the DTM method is proposed. Let matrix ${\mathbf A} \in {{\mathbf {\mathbb R}}^{{N_d},{N_d}}}$ denote the amplitude of DNF segment, where ${\mathbf A}({x_d},{y_d})$ represents the amplitude value in the position $({x_d},{y_d})$, and ${N_d}$ denotes the lateral dimension of the calculated DNF segments. Similarly, let matrix ${\mathbf M} \in {{\mathbf {\mathbb R}}^{{N_m},{N_m}}}$ denote the mask segment, where ${\mathbf M}({x_m},{y_m})$ represents the mask pixel in the position $({x_m},{y_m})$, and ${N_m}$ denotes the lateral dimension of mask segments. For a binary mask, M is consisted of 0 and 1 elements, which stand for the opaque pixels and transparent pixels, respectively. Suppose that each pixel of DNF amplitude can be formulated as the weighted summation of the mask pixels within the mask segment. Therefore, we have:

(16)$${\mathbf A}({x_d},{y_d}) = {T^{{x_d},{y_d}}} \odot {\mathbf M},$$

where ${T^{{x_d},{y_d}}} \in {{\mathbf {\mathbb R}}^{{N_m},{N_m}}}$ is a position-dependent transfer matrix. The elements in ${T^{{x_d},{y_d}}}$ stand for the weight parameters of the mask pixels that contribute to the DNF pixel $({x_d},{y_d})$, and ${\odot} $ represents the element-by-element multiplication. According to Eq. (16), the amplitude of DNF segment can be further expressed as:

(17)$${\mathbf A} = \sum\nolimits_{{x_d} = 1}^{{N_d}} {\sum\nolimits_{{y_d} = 1}^{{N_d}} {{T^{{x_d},{y_d}}} \odot {\mathbf M}} },$$

To further simplify the expression, let $\vec{{\mathbf A}} \in {{\mathbf {\mathbb R}}^{{N_d}^2,1}}$ and $\vec{{\mathbf M}} \in {{\mathbf {\mathbb R}}^{{N_m}^2,1}}$ denote the vectorized representations of the DNF segment A and the mask segment M, respectively. Therefore, Eq. (17) can be rewritten as:

(18)$$\vec{{\mathbf A}} = {\mathbf T\vec{M}},$$

where ${\mathbf T} \in {{\mathbf {\mathbb R}}^{{N_d}^2,{N_m}^2}}$ denotes the DTM, indicating the intensity mapping from the mask pattern to the DNF amplitude. Figure 7 illustrates the transformation from Eq. (17) to Eq. (18). It shows that ${\mathbf T}(i,j)$ stands for the contribution of the j th mask pixel to the i th DNF pixel. It is noted that the row dimension and column dimension of the DTM are equal to ${N_d}^2$ and ${N_m}^2$, respectively. The matrix T is the integration of all the elements ${T^{{x_d},{y_d}}}$ that correspond to the pixels in the DNF segment.

Using the matrix multiplication in Eq. (18), each pixel in DNF segment is separately calculated from the surrounding pixels on the mask, and the shift-invariant assumption of the diffraction behavior is no longer needed. Then, the edge interference effect can be accurately simulated by the proposed model. Simulation and comparison results are provided in Sec. 4.

Fig. 7. The relationship between the DNF segment, the mask segment, and the position-dependent transfer matrices. By vectorized the DNF segment and the mask segment, the expression can be further simplified to the form of matrix multiplication in Eq. (18). The DTM indicating the intensity mapping from the mask pattern to the DNF amplitude, which is the integration of all the position-dependent transfer matrices.

Download Full Size | PDF

The edge interference effect of DNF shows strong directivity. Therefore, the proposed model decomposes the mask pattern based on its geometric features. Specifically, the mask pattern is divided into four categories of non-overlapped segments: convex corners, concave corners, horizontal edges and vertical edges. Each category of features uses a different DTM, and these DTMs are calibrated for different categories independently. The optimal DTMs are calculated using the least square method based on the training data. The similar calibration process is applied to the filter-based model, and the simulation accuracy can be further improved. The proposed method shows promising performance in both of accuracy and speed.

3.3 DNF calculation algorithm based on DTM model

The detailed algorithm of the DTM model is described in this section, and the calculation of DNF amplitude is taken as an instance. Other DNF components (real part and imaginary part) can be calculated with similar algorithm. For the calculation of DNF amplitude, the relationship between the DNF amplitude and the mask pattern can be expressed as following:

(19)$$\begin{array}{l} {\mathbf A} = {{\mathbf A}_{Conv}} + {{\mathbf A}_{Conc}} + {{\mathbf A}_{HE}} + {{\mathbf A}_{VE}}\\ \left\{ {\begin{array}{{c}} {{{\vec{{\mathbf A}}}_{Conv}} = {{\mathbf T}_{Conv}}{{\vec{{\mathbf M}}}_{Conv}}}\\ {{{\vec{{\mathbf A}}}_{Conc}} = {{\mathbf T}_{Conc}}{{\vec{{\mathbf M}}}_{Conc}}}\\ {{{\vec{{\mathbf A}}}_{HE}} = {{\mathbf T}_{HE}}{{\vec{{\mathbf M}}}_{HE}}}\\ {{{\vec{{\mathbf A}}}_{VE}} = {{\mathbf T}_{VE}}{{\vec{{\mathbf M}}}_{VE}}} \end{array}} \right. \end{array},$$

where the subscripts Conv, Conc, VE, HE stand for the convex corner, concave corner, vertical edge and horizontal edge, respectively. In Eq. (19), A stands for the amplitude matrix of DNF, $\vec{{\mathbf A}}$ is the vectorized representation of DNF, $\vec{{\mathbf M}}$ is the vectorized representation of mask, and T is the DTM corresponding to different graphic features.

Figure 8 illustrates the flow chart of the algorithm. First, a training library of thick-mask DNF is established based on some training mask clips. The training DNF data is pre-calculated by the rigorous EMF simulator. Then, the training mask clips and their corresponding DNFs are decomposed into non-overlapped segments in the following procedures. First, a set of sampling points are anchored on the boundaries and corners of the mask clip with a fixed interspace, L₀. All of the sampling points are then extended simultaneously to cover their adjacent regions, which are referred to as mask segments. The extension of a segment terminates as it contacts other segments, or the radius of the segment reaches the preset upper limit. After the extension process, the training mask clips and their DNF images are decomposed into 4 kinds of features: the concave corner segments, convex corner segments, vertical edge segments and horizontal edge segments. In the decomposing process, patterns in different size can be divided into different amounts of segments in similar size, which can be further calculated with the same DTM. Therefore, the pattern-decomposing process increase the generalizability of the algorithm.

Fig. 8. Flow chart of the algorithm in DTM model. (a) The main process of the algorithm, and (b) the DTM calculation process

Download Full Size | PDF

Afterwards, based on the least square method, the DTM associated to each kind of features is calculated from the training dataset, which can be written as:

(20)$${{\mathbf T}_x} = {(\vec{{\mathbf M}}_x^T{\vec{{\mathbf M}}_x})^{ - 1}}\vec{{\mathbf M}}_x^T{\vec{{\mathbf A}}_x},$$

where the subscript x represents the type of segments (concave corner, convex corner, horizontal edge, or vertical edge), and the superscript T is the transpose symbol. It is worth noting that the calibrated DTMs can be saved and repetitively used in the following simulations.

Given a test mask pattern, we first decompose it into four groups of feature segments, similar as the training mask clips. Then, each segment is vectorized and multiplicated with the corresponding DTM in Eq. (20), so that the local DNF of the mask segment is calculated. Finally, the entire DNF image is synthesized by stitching up all of the non-overlapped DNF amplitude segments.

It is noted that, the algorithm mentioned above is not limited to calculate the amplitude of DNF. By changing the library data, the real part, and the imaginary part of the DNF can be calculated using the similar learning-based method. Then, the complex diffraction matrix of DNF can be obtained by composing the real part and imaginary part of DNF together:

(21)$${\mathbf E}({x_d},{y_d}) = {\mathbf A}({x_d},{y_d}) \cdot (\cos {\mathbf P}({x_d},{y_d}) + i\sin {\mathbf P}({x_d},{y_d})) = {\mathbf R}({x_d},{y_d}) + i{\mathbf I}({x_d},{y_d}),$$

where E represents the complex diffraction matrix; A and P respectively stand for the amplitude and phase of the DNF; R and I respectively denote the real part and imaginary part of the DNF; $({x_d},{y_d})$ denotes the spatial coordinate of the DNF image.

4. Result and discussion

In this paper, six mask clips with different geometries are selected as the training samples, and another mask clips is selected as the test sample. The dimension of these mask patterns is 1200 pixels×1200 pixels, where each pixel stands for a 1nm×1 nm square on the wafer scale. The waveguide method, a rigorous EMF simulator, is used to calculate the accurate DNF data of the training masks. The mode coverage factor of the waveguide method is set as 1.5. The simulations use an 193 nm immersion lithography system, where the numerical aperture of the projection optics is 1.35. For each mask pattern, we need to calculate the four diffraction matrices E(XX), E(XY), E(YX) and E(YY), which represent the complex amplitude of DNF. All of the simulation codes are implemented by MATLAB, and the computations are carried out on a personal computer with Intel Core i7-10870H CPU, 2.20 GHz, 32.0 GB of RAM.

In this work, the amplitudes of DNF and the aerial image are used to evaluate the accuracy of different models. The result of waveguide method is used as the benchmark to calculate the simulation error. For the other three methods, the root-mean-square (RMS) error of the simulation result is defined as:

(22)$$\textrm{RMS error} = \sqrt {\frac{1}{{{N_x} \cdot {N_y}}}\sum\nolimits_{x,y} {||{{{\mathbf R}_t}(x,y) - {{\mathbf R}_w}(x,y)} ||_2^2} },$$

where ${{\mathbf R}_t}$ and ${{\mathbf R}_w}$ denote the results obtained by the target method and waveguide method respectively; ${N_x}$ and ${N_\textrm{y}}$ denote the dimensions of mask pattern in two axes. The maximum absolute (MA) error is also used to evaluated the accuracy of methods. The MA error is defined as:

(23)$${\mathbf MA}\textrm{ error} = \max (|{{{\mathbf R}_t}(x,y) - {{\mathbf R}_w}(x,y)} |)$$

which stands for the maximum difference of the results between the target method and the benchmark method.

The test pattern and the training patterns of the library for the proposed model is displayed in Fig. 9(a) and 9(b), respectively. Figure 9(c) illustrates the RMS error of the DTM model using different numbers of training patterns. The RMS error decrease with the increasement of the training patterns, which demonstrate the proposed model is not over-fitting.

Fig. 9. (a) the test pattern and (b) the training patterns for library. (c) the curve of RMS error of the proposed model versus the number of training patterns. The pattern in panel (a) is selected as the test pattern.

Download Full Size | PDF

In the DTM model, the lateral dimensions of the calculated DNF segment (${N_d}$) and the mask segment (${N_m}$) are two key parameters that may affect the simulation accuracy. The size of a pixel is fixed at 1nm×1 nm on the wafer scale. In the proposed model, ${N_m}$ determines the area of the mask segment that taken into account in the DNF calculation, the optimal value of ${N_m}$ is determined by the illumination condition and other environment factors. On the other hand, ${N_d}$ determines the size of calculated DNF segments. In the pattern dividing process, the sampling points are arched at the edges and the corners of the pattern. Therefore, the optimal value of ${N_d}$ should be large enough to cover the main features of the mask pattern. The illustration of the circumstances when ${N_d}$ is too small is provided in the A. Appendix.

Figure 10 provides the simulation results of a test mask layout. From left to right, it shows the results calculated by the waveguide method, the Kirchhoff approximation model, the traditional filter-based model, and the proposed model, respectively. From top to bottom, it shows the amplitudes of the diffraction matrices E(XX), E(XY), E(YY), E(YX), respectively. The Kirchhoff approximation model treats the mask as an infinite thin film, and it ignores the polarization effect that is caused by the mask 3D structure. Therefore, it leads to a significant simulation error.

Fig. 10. The amplitudes of the diffraction matrices obtained by the EMF simulator, the Kirchhoff approximation model, the filter-based model, and the proposed model.

Download Full Size | PDF

As mentioned in Sec. 3.1, the shift-invariant assumption of diffraction behavior in the traditional filter-based model is invalid for small-CD pattern. The result of the filter-based model shows homogeneity in the amplitude distributions, and it fails to simulate the edge interference effect accurately. On the other hand, the proposed model retains more information of the ambient DNF distribution, and the edge interference effect is simulated accurately. Compared with the traditional filter-based model, the proposed model can effectively reduce the RMS error of the synthesized DNF amplitude by about 75%, and the MA error is also reduced significantly. The calculation speed of the proposed model is 2.3 times slower than the filter-based model, but is still about 300 times faster than the rigorous EMF simulator. Table 1 shows the runtimes and RMS errors of the four different models. The simulation results of different patterns are provided in the B. Appendix.

Table 1. The RMS errors, maximum absolute error and runtimes of the four different models.

View Table

Based on the proposed method, the real part and imaginary part of DNF can be calculated.

Figure 11 shows the absolute values of the real part and imaginary part of the four complex diffraction matrices. From top to bottom, it shows the results of E(XX), E(XY), E(YY), E(YX), respectively. From left to right, it shows the real parts and imaginary parts calculated by the EMF method, filter-based method, and proposed DTM method, respectively. As shown in Fig. 11, the homogeneous property of the filter-based model is obvious, which is also observed in the results of DNF amplitude simulation. On the other hand, the proposed model renders the simulation results with higher fidelity. The complex amplitude of each pixel can be further calculated by composing the real part and imagery part together according to Eq. (21). In the Kirchhoff approximation model, the imaginary part of DNF is 0, which is not shown in the Fig. 11.

Based on the Eq. (21), the complex amplitude of four diffraction matrices can be calculated with the real part and imaginary part of DNF, and the aerial images of the lithography system are further calculated to assess the simulation accuracy. Figure 12 illustrates the lithographic aerial images that are simulated based on the EMF method, the Kirchhoff approximation model, the filter-based model and the proposed model, respectively.

Fig. 11. The real part and imaginary part of the diffraction matrices. From left to right, the results are calculated by the EMF method, the filter-based method, and the proposed DTM model, respectively.

Download Full Size | PDF

Fig. 12. The aerial images that are simulated based on (a) the EMF method, (b) the Kirchhoff approximation model, (c) the filter-based model and (d) the proposed model.

Download Full Size | PDF

The aerial image result of the EMF method is used as benchmark to evaluate the accuracy of the other three methods. Figure 13 displays the absolute error maps of the aerial images that are calculated based on the Kirchhoff approximation model, the filter-based model and the proposed model, and the RMS errors are presented underneath the figures. It shows that, compared to the Kirchhoff approximation model and filter-based model, the proposed model can significantly improve the accuracy in aerial image calculation.

5. Conclusion

This paper introduces the edge interference effect that appears in the thick-mask DNF of sub-wavelength lithography. A high-precision thick-mask model is proposed to calculated the DNF taking into account the influence of the edge interference effect. In this method, a training library of the rigorous thick-mask DNFs is established firstly, and the DTM is learned based on the training dataset. Given an arbitrary thick-mask pattern, it is decomposed into a set of non-overlapped segments around the corners and edges, and the local DNF of each segment is quickly calculated using the pre-calibrated DTM. Then, the entire DNF image is synthesized by stitching up all of the local DNFs together. Compared to the traditional filter-based model, the proposed model can significantly improve the simulation accuracy in both of the DNF calculation and the aerial image calculation, meanwhile retaining high computational efficiency compared to the rigorous EMF method. In the future work, we may incorporate the proposed model with the lithography optimization program to improve the computational efficiency.

Fig. 13. The absolute error maps of the aerial images that are calculated by (a) the Kirchhoff approximation model, (b) the filter-based model, and (c) the proposed model.

Download Full Size | PDF

A. Appendix

In the pattern dividing process, the sampling points are arched at the edges and the corners of the pattern. Therefore, if the ${N_d}$ is too small, areas inside the bulk pattern can not be covered and calculated. Figure 14 illustrates the circumstance when the ${N_d}$ is too small.

Fig. 14. The illustration of circumstances when ${N_d}$ is big enough and when ${N_d}$ is too small.

Download Full Size | PDF

B. Appendix

The DNF amplitude simulation of different patterns are conducted, and the results are summarized in Fig. 15 (taking the polarization in E(XX) as an example). In each calculation, the tested pattern is removed from the training library. The performances of results are roughly identical, which confirm the generalizability of the proposed model.

Fig. 15. The simulation of different patterns that calculated with the EMF method and the proposed model. Corresponding RMS errors and MA errors are attached below figures.

Download Full Size | PDF

Funding

National Natural Science Foundation of China (61804174); Youth Innovation Promotion Association of the Chinese Academy of Sciences (No. 2021115).

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

References

1. F. M. Schellenberg, “Resolution Enhancement Techniques in Optical Lithography,” (SPIE, 2001), Chap. 1.

2. A. Erdmann, T. Fühner, P. Evanschitzky, V. Agudelo, C. Freund, P. Michalak, and D. Xu, “Optical and EUV projection lithography: A computational view,” Microelectron. Eng. 132, 21–34 (2015). [CrossRef]

3. X. Ma, X. Zhao, Z. Wang, Y. Li, S. Zhao, and L. Zhang, “Fast lithography aerial image calculation method based on machine learning,” Appl. Opt. 56(23), 6485–6495 (2017). [CrossRef]

4. R. Wu, L. Dong, X. Ma, and Y. Wei, “Compensation of EUV lithography mask blank defect based on an advanced genetic algorithm,” Opt. Express 29(18), 28872 (2021). [CrossRef]

5. Y. Peng, S. X. Shi, and D. Z. Pan, “Process variation aware OPC with variational lithography modeling,” in 2006 43rd ACM/IEEE Design Automation Conference (2006).

6. J. Lin, L. Dong, T. Fan, X. Ma, and T. Ye, “Learning-based compressive sensing method for EUV lithographic source optimization,” Opt. Express 27(16), 22563 (2019). [CrossRef]

7. X. Ma, Z. Wang, X. Chen, Y. Li, and G. R. Arce, “Gradient-Based Source Mask Optimization for Extreme Ultraviolet Lithography,” IEEE Trans. Comput. Imaging 5(1), 120–135 (2019). [CrossRef]

8. L. Pang, G. Xiao, V. Tolani, P. Hu, T. Cecil, T. Dam, K. H. Baik, and B. Gleason, “Considering MEEF in inverse lithography technology (ILT) and source mask optimization (SMO),” in Photomask Technology 7122, 71221W (2008). [CrossRef]

9. X. Ma, C. Han, Y. Li, L. Dong, and G. R. Arce, “Pixelated source and mask optimization for immersion lithography,” J. Opt. Soc. Am. A 30(1), 112–123 (2013). [CrossRef]

10. X. Ma and G. R. Arce, “Computational lithography,” Wiley, 58–82 (2010).

11. K. Lai, A. Erdmann, and J. Finders, “The impact of Mask 3D and Resist 3D effects in optical lithography,” Proc. SPIE 9052, 905205 (2014). [CrossRef]

12. A. Erdmann and P. Evanschitzky, “Rigorous electromagnetic field mask modeling and related lithographic effects in the low k1 and ultrahigh numerical aperture regime,” J. Micro/Nanolith. MEMS MOEMS 6(3), 031002 (2007). [CrossRef]

13. S. Burger, Z. Lin, F. Schmidt, P. Evanschitzky, and A. Erdmann, “Benchmark of Rigorous Methods for Electromagnetic Field Simulations,” Proc. SPIE 7122, 71221S (2008). [CrossRef]

14. A. Erdmann, P. Evanschitzky, G. Citarella, T. Fühner, and P. D. Bisschop, “Rigorous mask modeling using waveguide and FDTD methods: an assessment for typical hyper-NA imaging problems,” Proc. SPIE 6283, 628319 (2006). [CrossRef]

15. S. Burger, P. M. Martin, R. J. Naber, R. Köhle, L. Zschiedrich, H. Nguyen, F. Schmidt, R. März, and C. Nölscher, “Rigorous Simulation of 3D Masks,” Proc. SPIE 6349, 63494Z (2006). [CrossRef]

16. P. Evanschitzky and A. Erdmann, “Fast near field simulation of optical and EUV masks using the waveguide method,” Proc. SPIE 6533, 65330Y (2007). [CrossRef]

17. J. T. Azpiroz, G. W. Burr, A. E. Rosenbluth, and M. Hibbs, “Massively-parallel FDTD simulations to address mask electromagnetic effects in hyper-NA immersion lithography,” Proc. SPIE 6924, 69240Y (2008). [CrossRef]

18. Q. Yan, Z. Deng, and J. Shiely, “Fast synthesis of topographic mask effects based on rigorous solutions,” Proc. SPIE 6730, 67302N (2007). [CrossRef]

19. S. G. Lee, K. I. Lee, J. U. Lee, Y. Choi, and H. S. Kim, “More stable algorithm for rigorous coupled wave analysis applied to topography simulation in optical lithography and its numerical implementation,” Proc. SPIE 2726, 288–298 (1996). [CrossRef]

20. X. Ma and G. R. Arce, “Binary mask optimization for forward lithography based on the boundary layer model in coherent systems: erratum,” J. Opt. Soc. Am. 27(7), 82–84 (2010). [CrossRef]

21. J. Tirapu-Azpiroz, P. Burchard, and E. Yablonovitch, “Boundary layer model to account for thick mask effects in PhotoLithography,” Proc. SPIE 5040, 1611–1619 (2003). [CrossRef]

22. X. Ma, G. R. Arce, and Y. Li, “Optimal 3D phase-shifting masks in partially coherent illumination,” Appl. Opt. 50(28), 5567–5576 (2011). [CrossRef]

23. K. Ad Am and A. R. Neureuther, “Methodology for accurate and rapid simulation of large arbitrary 2D layouts of advanced photomasks,” Proc. SPIE 4562, 1051–1067 (2002). [CrossRef]

24. L. Peng, “Accurate prediction of 3D mask topography induced best focus variation in full-chip photolithography applications,” Proc. SPIE 8166, 816640 (2011). [CrossRef]

25. H. Zhang, Q. Yan, E. Croffie, L. Zhang, Y. Fan, T. B. Faure, and P. W. Ackmann, “An accurate ILT-enabling full- chip mask 3D model for all-angle patterns,” Proc. SPIE 8880, 88800G (2013). [CrossRef]

26. H. Zhang, Q. Yan, Z. Lin, E. Croffie, and Y. Fan, “Eficient full-chip mask 3D model for Off-Axis Illumination,” Proc. SPIE 8880, 888023 (2013). [CrossRef]

27. P. Liu, P. P. Naulleau, X. Xie, W. Liu, and K. Gronlund, “Fast 3D thick mask model for full-chip EUVL simulations,” Proc. SPIE 8679, 86790W (2013). [CrossRef]

28. L. Peng, C. Yu, L. Chen, G. Chen, and S. Lee, “Fast and accurate 3D mask model for full-chip OPC and verification,” Proc. SPIE 6520, 65200R (2007). [CrossRef]

29. P. Liu, R. J. Naber, H. Kawahira, C. Zuniga, Z. Ma, and H. Feng, “Validation of a fast and accurate 3D mask model for SRAF printability analysis at 32 nm node,” Proc. SPIE 6730, 67301R (2007). [CrossRef]

30. X. Ma, J. Gao, X. Chen, L. Dong, and Y. Li, “A fast and flexible library-based thick-mask near-field calculation method,” Proc. SPIE 9426, 94261G (2015).

31. W. Ye, M. B. Alawieh, Y. Watanabe, S. Nojima, Y. Lin, and D. Z. Pan, “TEMPO: Fast mask topography effect modeling with deep learning,” in Proceedings of the 2020 International Symposium on Physical Design (2020), pp. 127–134.

32. X. Zheng, X. Ma, Q. Zhao, Y. Pan, and G. R. Arce, “Model-informed deep learning for computational lithography with partially coherent illumination,” Opt. Express 28(26), 39475–39491 (2020). [CrossRef]

33. A. Awad, P. Brendel, P. Evanschitzky, D. Woldeamanual, A. S. Rosskopf, and A. Erdmann, “Accurate prediction of EUV lithographic images and 3D mask effects using generative networks,” J. Micro/Nanopattern. Mats. Metro. 20(04), 043201 (2021). [CrossRef]

34. Y. Cao, X. Wang, A. Erdmann, B. Peng, and B. Yang, “Analytical model for EUV mask diffraction field calculation,” Proc. SPIE 8171, 81710N (2011). [CrossRef]

35. O. R. Wood, E. M. Panning, X. Liu, X. Wang, S. Li, G. Yan, and A. Erdmann, “Fast rigorous model for mask spectrum simulation and analysis of mask shadowing effects in EUV lithography,” Proc. SPIE 9048, 90483E (2014). [CrossRef]

36. X. Liu, X. Wang, S. Li, G. Yan, and A. Erdmann, “Fast model for mask spectrum simulation and analysis of mask shadowing effects in extreme ultraviolet lithography,” J. Micro/Nanolith. MEMS MOEMS 13(3), 033007 (2014). [CrossRef]

37. M. Born and E. Wolf, “Principles of Optics” (Pergamon Press, 1959), Chap. 11.

method	error (compared with the EMF method)				Time
	XX	XY	YX	YY
EMF method	-	-	-	-	∼4500s
Kirchhoff (RMS error)	0.1356	-	-	0.1603	-
Kirchhoff (MA error)	1.1109	-	-	1.1139
Filter-based (RMS error)	0.1080	0.0035	0.0034	0.1356	6.2s
Filter-based (MA error)	1.5704	0.1297	0.0601	0.9559
Proposed model (RMS error)	0.0279	0.0012	0.0012	0.0071	14.3s
Proposed model (MA error)	0.3658	0.0411	0.0382	0.4154

High-precision lithography thick-mask model based on a decomposition machine learning method

Abstract

1. Introduction

2. Edge interference effect

3. Traditional filter-based model and the proposed model

3.1 Traditional filter-based model

3.2 Proposed DTM model

3.3 DNF calculation algorithm based on DTM model

4. Result and discussion

5. Conclusion

A. Appendix

B. Appendix

Funding

Disclosures

Data availability

References

Data availability

Cited By

Figures (15)

Tables (1)

Equations (23)

Optics Express