Stochastic simulation and calibration of organometallic photoresists for extreme ultraviolet lithography

Zelalem Belete; Peter De Bisschop; Ulrich Welling; Andreas Erdmann

doi:10.1117/1.JMM.20.1.014801

4 March 2021 Stochastic simulation and calibration of organometallic photoresists for extreme ultraviolet lithography

Zelalem Belete, Peter De Bisschop, Ulrich Welling, Andreas Erdmann

Journal of Micro/Nanopatterning, Materials, and Metrology, Vol. 20, Issue 1, 014801 (March 2021). https://doi.org/10.1117/1.JMM.20.1.014801

Abstract

Organometallic photoresists are being pursued as an alternative photoresist material to push the current extreme ultraviolet lithography (EUVL) to the next generation of high-NA EUVL. In order to improve the photoresist performance, an understanding of the photoresist’s response to different process conditions is required. In this endeavor, a stochastic development model is implemented, integrated into full photoresist process steps, and applied for photoresist performance investigations. The model is applied to Inpria-YA photoresist, which works mainly by the process of aggregation. Previously published modeling approaches for metal-organic photoresists assume that the development characteristics of these materials depend only on the size of the created oxo-clusters. In contrast to that, we propose a modeling approach that provides a more detailed description of the interaction among the developer, ligands, and oxo-bonds. Further, the calibration procedures conducted to extract the model parameters to match experimental data are discussed. The model approximated the experimental data with CD RMSE and LWR RMSE of 0.60 and 0.40 nm, respectively. We also investigated the impact of photoresist parameters on the process metrics, line width roughness (LWR), critical dimension (CD), dose-to-size (DtS), and exposure latitude (EL) with the calibrated model. The investigation shows that details of the interaction of photoresist and developer, especially, the so-called development critical value, have a significant impact on the LWR and DtS.

1. Introduction

Extreme ultraviolet lithography (EUVL) increases the resolution using a smaller wavelength (13.5 nm) for exposure compared to DUV lithography (193 nm). However, this improvement using smaller wavelength also led to several problems, such as stochastic effects.¹ These stochastic effects deteriorate the photoresist performance and limit the economically viable scaling. So far, even if the feature size is decreasing, the line width roughness (LWR) remained constant, which makes the LWR as the limiting factor.

Due to the small number of photons (large shot noise) in EUV lithography and the small absorption efficiency of chemically amplified photoresists (CARs), new photoresist materials are essential. As a result, the industry is pursuing different photoresist materials for the next generation of EUV high-volume manufacturing (HVM). Molecular organometallic resists for EUV are one of the alternative photoresist materials. They have a higher absorption ( $\sim 16$ to $20 μ m^{- 1}$ ) compared to CARs ( $\sim 4 μ m^{- 1}$ ).² These photoresists have a core–shell structure where the core contains metal oxide molecules, and the shell contains radiation-sensitive ligands. The absorption efficiency depends mainly on the type of metal used. Inpria-YA photoresists have a tin-oxide core that enables their high photon absorption.²^–⁴ The radiation-sensitive ligands control the reaction behavior of the photoresist, preventing background condensation reaction.⁵ Organometallic photoresists are designed to absorb photons that lead to cleavage of the ligands and subsequent chemical reactions to change the solubility. The cleavage of the ligands site creates an “active” site that can create a bond with another “active” site on an adjacent metal-oxide core. As this condensation reaction continues and active sites on neighbor metal-oxide cores form a bond, it creates a large networked structure (oxo-network) that is resistant to the development.⁶

Previously, several modeling approaches and investigations of organometallic photoresists have been reported.⁶^–¹¹ The presented models provide a good understanding of the organometallic photoresist behavior during exposure. In this paper, a modeling procedure to characterize and quantify the development process is presented together with a calibration of the model with experimental data to contribute further to understanding of the photoresists’ response to different processes. Section 2 explains the modeling procedures applied for exposure, condensation, and development. The calibration of the developed model with experimental Bossung data and verification procedures of the calibration results are discussed in Sec. 3. A simulation study on the impact of photoresist parameters on the lithography metrics, LWR, critical dimension (CD), exposure latitude (EL), and dose-to-size (DtS) is presented in Sec. 4. Finally, the results and findings are summarized, and an outlook for future work is discussed.

2. Modeling Procedure

2.1.

Exposure

The photoresist simulation volume is discretized into lattice cells of size $δ x$ , $δ y$ , and $δ z$ , where $δ x * δ y * δ z$ is the volume of the photoresist molecule (metal-oxide core with ligands). The $x$ and $y$ represent directions perpendicular to the exposure direction, and $- z$ represents the direction of exposure. The intensity distribution inside the photoresist (bulk image) is simulated using the Fraunhofer lithography simulator, Dr. LiTHO.¹² Then the intensity absorbed by the photoresist molecule is computed from the bulk image according to the Beer–Lambert law, as described in the following equation:

Eq. (1)

I_{abs, i} = I_{i} * (1 - e^{- α * δ z}),

where

I_{i}

and

I_{abs, i}

are the bulk image intensity and the absorbed intensity for the

i_{th}

lattice cell, respectively.

α

is the absorption coefficient and

δ z

represents the thickness of a single lattice cell.

From the computed absorbed intensity distribution, the average number of photons absorbed $(⟨ N_{p, i} ⟩)$ in the $i_{th}$ single lattice cell is given by

Eq. (2)

⟨ N_{p, i} ⟩ = \frac{I_{abs, i} * D * δ x * δ y}{h c / λ},

where

D

is the exposure dose and

(δ x * δ y)

is the surface area of a single lattice cell.

h

represents the Plank’s constant (

6.62607004 \times 10^{- 34} Js

),

c

is the speed of light (

3 \times 10^{8} m / s

), and

λ

is the wavelength of the exposure (13.5 nm for EUV).

The stochastic distribution of the photons in the photoresist volume follows a Poisson distribution.¹³^,¹⁴ Therefore, the actual number of photons absorbed by the photoresist molecules is distributed according to a Poisson distribution from the average number of absorbed photons.

The absorbed photons trigger subsequent chemical changes to modify the solubility of the photoresist. For EUVL, the chemical change is mainly initiated by the photoelectric effect where photon absorption leads to the generation of electrons.⁷^,¹⁵ The generated photoelectron, in turn, produces secondary electrons for further ionization of the photoresist molecules at a distance from the point of photon absorption—resulting in a blurring effect on the distribution of the absorbed photons.¹⁶ As tracking each electron and its interaction with the photoresist molecules is computationally intensive, a simplified approach is implemented in our model to simulate these processes.

The average number of generated electrons per absorbed photon is approximated by the quantum yield.

Eq. (3)

N_{e} = N_{p} * Φ,

where

N_{e}

and

N_{p}

represent the number of generated electrons and absorbed photons, respectively.

Φ

is the quantum yield, the actual number of electrons generated per absorbed photon. The actual number of electrons generated per photon is computed by generating the actual quantum yield according to the Poisson probability distribution from the average quantum yield

⟨ Φ ⟩

and Eq. (3).

After the computation of the generated electrons, the electron blur due to the electrons’ random walk is approximated by Gaussian convolution. Such approach enables a significant reduction of computing time. The Gaussian kernel applied for the convolution is defined as in the following equation:

Eq. (4)

G (x, y, z) = \frac{1}{{(2 π σ_{blur}^{2})}^{\frac{3}{2}}} e^{- \frac{x^{2} + y^{2} + z^{2}}{(2 σ_{blur}^{2})}},

where

σ_{blur}

is the electron blur length in nm. Application of this convolution results in a (deterministic) blur of the location of the generated electrons with respect to the point of photon absorption.

Typical photon and electron distributions that are computed from the bulk image (generated by Dr. LiTHO) are shown in Fig. 1.

Fig. 1

EUV lithography exposure simulation of 18 nm vertical lines with 36 nm pitch with dipole illumination and $52 mJ / {cm}^{2}$ dose for a 0.33-NA system. (a) Bulk image computed using Dr. LiTHO. For the computation of the electron distribution, $Φ$ of 8 and blur length of 2.4 nm are used. (b) Stochastic photon distribution and (c) stochastic electron distribution.

2.2.

Condensation Reaction

The photoresist changes its solubility after exposure due to physical or chemical structure changes. Several studies demonstrated that the metal-oxide cores absorb the photons to generate photoelectrons. However, the ligands act only as non-reactive spacers inhibiting the metal-oxide cores from reacting.¹⁷^,¹⁸ During exposure, the protecting ligands are cleaved by electrons. This cleavage of the ligands leads to a creation of what is called an “active site,” i.e., a site where a condensation reaction [i.e., the creation of a bond (M–O–M oxo-bridge) with an active site of another adjacent metal-oxide core] can occur.⁶^,¹⁷

The number of active sites generated on each particular metal-oxide core follows directly from the number of electrons that land in the $δ x$ , $δ y$ , $δ z$ volume element of the metal-oxide cluster. Whether or not an active site will lead to the formation of an M–O–M bridge is simulated in a probabilistic approach, using percolation theory.⁶^,⁹^,¹⁹

In the model used in Ref. 6, a metal-oxide core and a neighbor metal-oxide core are randomly selected, and if both cores have active sites, a bond is created. If there are no available active sites, another core and its neighbor are randomly selected. However, this implementation is computationally intensive, and it is not well suited for the large number of computations required for a calibration procedure and for extensive modeling studies.

Therefore, our model uses a different, pseudo-random, approach that is faster without losing the stochastic nature of the process. The operation scheme is summarized in the pseudo-code in Algorithm 1. The metal-oxide cores with active sites are visited sequentially. Then one (and just one) of the six nearest neighbor metal-oxide cores is randomly selected. If this particular neighbor core also has an active site, an oxo-bond between the two cores is created, replacing the two active sites (i.e., the bond consumes the two active sites). Otherwise, the active site remains active. Then the next active site on the given metal-oxide core is selected, and another neighbor metal-oxide core is randomly selected. This process is iterated until all metal-oxide cores, and all of their active sites are visited. The pseudo-code of the implementation is presented in Algorithm 1.

Algorithm 1

Quasi-random implementation of condensation reaction

1: start

2: for $k \leftarrow 1$ to number of cores with active sites do

3: for $l \leftarrow 1$ to number of actives sites do

4: randomly select one single neighbor (from the six nearest neighbors)

5: if randomly selected neighbor has active site then

6: create bond

7: consume active sites from current and neighbor cell

8: else

9: ignore

10: end

A periodic boundary condition is applied for $x$ - and $y$ - directions. If the randomly selected neighbor’s index gets out of the simulation domain, a bond is created with the metal-oxide core on the opposite side.

Post exposure bake (PEB) impacts the condensation reaction and, in turn the sensitivity of the photoresist.²⁰ However, our investigation focuses on the modeling of the stochastic development process for a finished condensation reaction. Therefore, we do not consider explicit modeling of the PEB process.

The oxo-networks computed from the electron distribution, for a contact feature, with exposure doses of 20, 40, and $90 mJ / {cm}^{2}$ are shown in Fig. 2. For a small exposure dose, the created oxo-networks are small in size. As the dose increases, the size of the oxo-networks that are created in the exposed region increase, whereas small oxo-clusters are created in the unexposed regions. For large exposure dose, a large oxo-cluster is created in the exposed region [Figs. 2(c) and 2(d)].

Fig. 2

Size distribution of the oxo-networks (represented by the colors) created for dose values of (a) 20, (b) 40, (c) $90 mJ / {cm}^{2}$ for 18 nm contact with 72 nm pitch. (d) The largest oxo-network created for $90 mJ / {cm}^{2}$ exposure dose [from the figure shown in (c)]. For the simulation, an attenuated phase shift mask with thickness of 31.6 nm and quasar illumination with $σ_{outer} = 0.4$ are used. The parameters are taken from Ref. 21.

2.3.

Development

After the solubility is changed for the exposed regions of the photoresist due to the condensation reaction, the photoresist is developed with a solvent. This development step creates the final photoresist profile.

The development behavior of organometallic photoresist that aggregates during exposure, to create a large oxo-network, is mainly dependent on the size of the created oxo-cluster. Small oxo-clusters can develop quickly, whereas large oxo-clusters are resistant to the developer. There are several alternatives to implement the development process, such as simple threshold,⁸ percolation model,²² or tracking of the oxo-network in contact with the photoresist.⁶ The investigation of stochastic effects during the development requires a more detailed description of the dissolution mechanism. Several experimental investigations²³^–²⁷ demonstrated that the photoresist dissolution behavior depends on the behavior of the ligands on the metal cores.

Our stochastic development model tracks the substitution of the ligand and oxo-bond sites by the solvent molecules based on balance equations. If the number of these solvent substituted sites per metal-oxide core is above a certain number and all the oxo-bonds on the core are broken, the core is washed away. This procedure is iterated until the developer time is reached. For this purpose, we employed the basic approach of the critical ionization model to track the interaction of the ligands with solvent molecules in stochastic development model, based on the implementations described in Refs. 28 and 29. The reasoning and details of the model formulations are explained below.

At the beginning of the development process, the metal-oxide cores in the photoresist have three kinds of sites, as shown in Fig. 3.

1. Ligand site. The ligand is not cleaved, i.e., it is not affected by the exposure.
2. Active site. An active site is created, but it did not form an oxo-bond during condensation.
3. Oxo-bond site. An active site created during exposure is replaced by a bond during condensation.

Fig. 3

(a)–(d) Possible states of the molecules in the photoresist during the development: (a) ${MOL}_{4}$ ; (b) ${MOL}_{2} A B$ ; (c) ${MOL}_{2} AS$ ; and (d) MOLABS. The cube represents a single grid cell that contains the molecule. The circles, the line, and the star represent the components in the molecules. Metal-oxide core (MO), ligands (L), and solvent (S) are represented by green, gray, and purple circles, respectively. Red stars represent active sites (A) and the oxo-bonds (B) are represented by a green line. (e) and (f) Reactions of solvent molecules with ligand and oxo-bond sites during development (active sites are treated as ligands and are therefore not shown).

The illustrations in Fig. 3 place the components on a specific geometrical location in the grid cell. In the simulation, however, no such geometrical-location information is used. In each cell, we only keep track of the number of L, A, B, and S components.

During the development in organic solvents, the ligands and active sites can be dissociated from the metal-oxide core and substituted by a solvent molecule.³⁰^,³¹ Ligands from the developer bulk can also substitute a solvent site on a metal-oxide core again, in a reverse reaction. Oxo-bond sites are normally assumed to be resistant to the developer. With strong developer solutions, however, the bonds can be broken and substituted by solvent molecules. A metal-oxide core is removed from the photoresist into the developer solution, if the number of solvent substituted sites per metal-oxide core is larger than the development critical number of sites ( $L_{th}$ ).

In order to reduce the complexity and the number of parameters in our model, ligand and active sites are assumed to be similar and will be treated identically in the development step. In this case, the active sites will be treated as ligand sites. This limits the types of sites on the metal-oxide cores, at the start of the development process, to only two kinds—ligand and oxo-bond sites. For example, ${MOL}_{2} AB$ is treated as ${MOL}_{3} B$ . Using this assumption, the population balance equations for substitution reactions of the development processes at the interface between the developer and the photoresist can be defined by the following equation:³⁰

Eq. (5a)

{MOL}_{n} B_{6 - n} + S ⇄_{K_{2}}^{K_{1}} {MOL}_{n - 1} B_{6 - n} S + L,

Eq. (5b)

B_{6 - n - 1} L_{n} MO \binom{\underset{⏟}{-}}{oxo-bond} {OML}_{m} B_{6 - m - 1} + 2 S \overset{K_{3}}{\to} {MOL}_{n} B_{6 - n - 1} S + {MOL}_{m} B_{6 - m - 1} S,

where S is the solvent and,

K_{1}

,

K_{2}

, and

K_{3}

are the rate constants for the ligand substitution, the solvent substitution, and the oxo-bond substitution, respectively. M, O, L, and B represent metal, oxygen, ligand, and oxo-bond on the photoresist molecule and

n

and

m

represent the number of ligands on the given metal-oxide cores. The mechanisms of the population balance equations, Eqs. (5a) and (5b) are illustrated in Figs. 3(e) and 3(f), respectively. Note that, as we consider MO-bonds with one of the six direct neighbors of each cell, and our model only counts the number of L, A, B, and S in each cells, we need to distinguish no more than six ligands per core, even if the real resist molecules may well have more than six ligands per core.

Our model describes the process of stochastic events during the development by an implementation of the Gillespie algorithm. This algorithm employs Eq. (5) and given reaction rates to generate stochastically correct, possible trajectories of the developer front. This means that the reactions are simulated explicitly as a stochastic process instead of solving the equation analytically.

In the model’s implementation, the solvent molecules’ reaction with ligand and oxo-bond sites of the cores are simulated by discrete events with probabilities that are governed by the reaction rates. Only the cores that are in contact with the developer are considered. To track the reactions by separate events, a time step, with a probability that a single reaction occurs in this time interval, must be first defined. This time step is computed based on the reaction rates ( $K_{1}$ , $K_{2}$ , and $K_{3}$ ) and the corresponding number of sites on the cores. Then the time is updated. In each time step, substitution reactions are determined with a probability for all the cores on the development front. The reaction can be either ligand or oxo-bond or solvent substitution reaction. The probability of the selected substitution reaction depends on the reaction constants and governs the update of the sites on the cores. If the number of solvent substituted sites on a core is larger than the development critical value ( $L_{th}$ ), and all oxo-bonds are broken, the core is removed, and its neighbor cores are added to the developer front. Then the next time step is computed, and the process is repeated. These processes are iterated until the total development time is reached. The developer front trajectory, tracked by the steps described above, creates the final photoresist profile.

Details on this implementation are given in Appendix A.

The implemented development model transforms an oxo-bonds distribution, as obtained from the condensation reaction, into the final latent profile. The distribution of the oxo-bonds per core is shown in Fig. 4(a) and the developer front at different development time steps are shown in Figs. 4(b)–4(d).

Fig. 4

(a) Distribution of the oxo-bonds per core after the condensation reaction. (b)–(d) Developer front at time = 10, 15, and 30 s from the computed oxo-network (or number of bonds per core) for EUV lithography exposure simulation of 18 nm vertical lines with 36 nm pitch with dipole illumination and $52 mJ / {cm}^{2}$ dose.

For a standard implementation of the model in Python, simulation of a vertical line with 36 nm pitch and 100 nm line length, carried out on only a single-core on a PC with Intel-Core i7-4770 at 3.4 GHz (CentOS Linux 7) and 4 GB RAM, requires $\sim 9.5 h$ . The optimized implementation of the algorithm described in Appendix A reduces the computation time to only $\sim 100 s$ .

2.4.

Illustrations of the Complete Process Flow

A schematics that illustrates the working principle of the model is presented in Fig. 5. The schematics shows the photoresist molecules and the effect of the different processes on the sites of the metal-oxide cores; the activation of the ligands during exposure, random creation of oxo-bonds during condensation, and the evolution of the developer front during development.

Fig. 5

Molecular representation of the modeling process flow (2D scheme).

Representative simulation results for the modeling process flow of the organometallic photoresist are shown in Fig. 6. The process starts with the computation of the bulk image [Fig. 6(a)] based on the exposure conditions, such as the source and mask specifications. The average number of photons absorbed in the photoresist is computed from the bulk image by specifying the exposure dose and and the absorption coefficient ( $α$ ) based on Eq. (2). The average number of photons absorbed in the photoresist is redistributed according to the Poisson probability distribution. This redistribution process generates a stochastic distribution of the photons [Fig. 6(b)].

Fig. 6

Stochastic distributions (cross-section plots) for different process steps for simulation of EUV exposure of 18 nm vertical lines with 36 nm pitch with dipole illumination, $34 mJ / {cm}^{2}$ dose and $0.065 μ m$ focus. Simulation parameters: $α = 20 μ m^{- 1}$ , $Φ = 8.2$ , $σ_{blur} = 3.18 n m$ , $L_{th} = 3$ , $K_{1} = 11.9 s^{- 1}$ , $K_{2} = 5.5 s^{- 1}$ , and $K_{3} = 7.4 s^{- 1}$ . (a) The bulk image; the stochastic distributions of (b) absorbed photons; (c) generated electrons; (d) oxo-bonds per core; and (e) the DArT.

The electrons generated from the absorbed photons are computed from the stochastic distribution of the photons and the quantum yield ( $Φ$ ), based on Eq. (3). To describe the electron blur effect, the distribution of the electrons is convoluted with the Gaussian distribution with a standard deviation equal to the blur length ( $σ_{blur}$ ), as described in Eq. (4). The final electron distribution is shown in Fig. 6(c). Afterward, the electrons’ distribution is related to the activation of ligands on the metal-oxide cores, and oxo-bonds are created between activated sites of neighboring cores by applying the percolation model (see Algorithm 1). This procedure results in the distribution of oxo-bonds, displayed in Fig. 6(d). Finally, the development step is simulated based on the distribution of oxo-bonds, applying the model parameters, the development critical value ( $L_{th}$ ), the rate constants ( $K_{1}$ , $K_{2}$ , and $K_{3}$ ), and the development time (see Algorithm 2 in Appendix A). The final result of the development simulation, the developer arrival time (DArT), is shown in Fig. 6(e).

Table 1 summarizes the model parameters and sources of variability in each process step, together with the implementation procedures that are applied to capture these variabilities. The variability that is expressed by the LWR is a combined result of stochastics in the photon absorption, the electron generation, the bond creation, and the interaction of the developer molecules with the ligands and bonds of the photoresist molecules.

Table 1

Model parameters, sources of variability, and implementation procedure.

Module	Model parameters	Source of variability	Variability implemented by
Exposure	$α$ , $σ_{blur}$ , $Φ$	Probabilistic photon absorption	Poisson probability function
Exposure	$α$ , $σ_{blur}$ , $Φ$	Probabilistic electron generation	Gaussian probability for spatial distribution of electrons
Condensation		Random bond creation	Random probing of neighbor cores (Algorithm 1)
Development	$K_{1}$ , $K_{2}$ , $K_{3}$ , $L_{th}$	Probabilistic substitution of L, B, and S	Gillespie algorithm
			Probabilistic selection of reacting cores
			Probabilistic selection of reaction types
			Poisson probability number of reactions

3. Calibration of Model with Experimental Data

The parameter values for the developed photoresist model, summarized in Table 1, are extracted by a model calibration with experimental data. The calibration procedure is conducted using lines-and-space (L/S) Bossung data (CD and LWR values) with exposure dose values ranging from 41 to $63 mJ / {cm}^{2}$ and focus values from $- 0.05$ to $0.15 μ m$ in steps of $0.02 μ m$ . The Bossung experimental data were obtained from Imec for Inpria-YA photoresist, with pixelated dipole illumination, shown in Fig. 7(b) at $NA = 0.33$ . The mask contains a multi-layer with 40 layers of Mo/Si bilayers with a Ru cap layer. On the top of the multi-layer, a patterned 55 nm Ta absorber is used. The wafer stack is composed of 22 nm Inpria-YA photoresist on 10 nm spin-on-glass underlayer.

The experimental data contain 758 measurement points for vertical L/S patterns with 5 mask CD and pitch variations, as shown in table in Fig. 7(a) for a center-slit position. For our model calibration, only 73 data points are selected to keep the computation time to a minimum but also include enough data points that can define the model behavior, simultaneously. The calibration data points are shown in red X marks in Fig. 10. All measurement data points are used for verification of the calibrated parameter values.

Fig. 7

(a) Feature types and (b) source used for the measurement. The feature types are specified as P “pitch” V “mask CD,” with the specified mask CD and pitch. “V” stands for vertical line.

The experimental data are from CD-SEM measurements conducted for 16 scans (or frames). During multiple scans of the CD-SEM measurement, the electron beam shrinks the photoresist lines. Our calibration is based on CD-SEM data from the single-frame and 16th-frame. These data exhibit a CD shrinkage in the range of 2.5 to 6 nm between single frame and 16th-frame. Single-frame SEM images have a non-negligible but small shrinkage compared to the 16th-frame image. In our calibration and verification, we used the CD data from the single-frame measurement. However, the single-frame measurements are noisy and cannot provide a reliable LWR measurement. Therefore, we used the LWR data from the 16th-frame SEM image for the calibration. The effect of the SEM on the roughness due to the photoresist lines shrinkage is neglected. For further details on SEM measurement impact on the photoresist shrinkage and a correction procedure, refer to the discussion in the Appendix of Ref. 32.

3.1.

Simulation Parameters and Process Conditions

The process conditions and simulation parameters in Table 2 are used to compute the bulk image using the imaging module of Dr. LiTHO. The computed bulk images provide the input for the photoresist model—to compute the final latent image. The molecular volume ( $δ x * δ y * δ z$ ) of the photoresist is assumed to be $2.3 {nm}^{3}$ .⁶

Table 2

Fixed parameters used in the simulation for the calibration of the model.

Input parameters	Symbols	Value	Reference
Molecular volume ( ${nm}^{3}$ )	$V_{cell}$	2.3	6
Refractive index	n	0.998	33
Number of ligands per core	$L_{core}$	12	6
Development time (s)	$t_{dev}$	30
Photoresist thickness (nm)	th	22

Due to computation constraints, our calibration runs are simulated for line lengths of 100 nm to compute the CD and LWR, for this prototype photoresist model implementation, with a loss of small accuracy but with faster computation time. CD and LWR values are extracted from the CD values computed along the line and over a certain thickness range of the final latent image. Maas et al.⁸ investigated the probe depth of the electron beam for 300-eV CD-SEM measurement for organometallic photoresist to be the top 8 nm of the photoresist. For our calibration data, CD-SEM with a 500-eV electron beam was employed, close to double the CD-SEM voltage used for the simulation by Maas et al. Therefore, we assumed that only the top 16 nm of the photoresist, twice the probe depth calculated for 300 eV, could affect the CD-SEM measurement. As a result, the top 16 nm thickness of the photoresist is considered for the computation of CD and LWR, without considering any thickness loss. For a simulation of a line feature with 100 nm line length and 22 nm photoresist thickness, this corresponds to 900 CD values computed from the final latent profile. The corresponding LWR value for the simulated profile is computed from these computed CD values and the selected sampling along the line and versus height, respectively.

Pythmea,³⁴ a python multi-objective evolutionary algorithm from Dr. LiTHO, is used for the calibration of the model to the experimental data. The model parameters are calibrated for CD and the LWR simultaneously using a multi-objective optimizer. The fitness of the calibration was determined, for both LWR and CD, from a root-mean-square error (RMSE) of simulation results from the experimental data. After the model is calibrated, the verification runs were conducted with 300 nm line length for 758 data points for a calibration result selected from the Pareto-optimal solutions.

The photoresist model parameters, absorption coefficient, quantum yield, electron blur length, and development model parameters, were varied during the calibration (Table 4). Because the simulation zero-focus position deviates from the exposure tool zero-focus position with an unknown offset, a focus offset parameter is included in the calibration. Additionally, metrology offset parameters, to compensate for systematic deviations of CD and LWR measurements, are also included in the calibration.

3.2.

Model Options, Calibration Results, and Verification of the Models

Calibration was conducted under three different assumptions.

1. Model 1. In this model, it was assumed that the developer does not affect the oxo-bonds created during the condensation reaction. For this procedure, only Eq. (5a) is considered and the metal-oxide cores are developed if the number of solvent substitute sites is above the critical value, irrespective of the number of oxo-bonds on the core. This assumption is valid for a solvent with a small dielectric constant or low polarity.³⁵ The CD RMSE and LWR RMSE values for the calibration and the verification, for one of the selected solution from the Pareto front, are summarized in Table 3.
2. Model 2. This model includes the substitution reaction for the oxo-bonds. Yeh et al.²³ demonstrated that developers with high dielectric constants, such as ethanol, break oxo-bonds, and damage the oxo-networks created during the condensation reaction. AFM images of the patterned lines, patterned using ethanol as a developer, exhibit a rough surface with apparent damage to the oxo-cluster. In contrast, patterned lines developed in 1-propanol have a smooth surface. Even though the patterned lines developed in 1-propanol have a smooth profile, it does not necessarily mean the breaking of the oxo-bonds does not occur. Therefore, in this model, both balance equations [Eqs. (5a) and (5b)] are considered. Moreover, it is enforced that every oxo-bond on a photoresist core should be broken before the core is considered as developed. This condition has to be satisfied even if the number of solvent substituted sites on the core is above the development critical value. The addition of this assumption increased the sensitivity of the development model with respect to number of oxo-bonds per core. As summarized in Table 3, model 2 fits the experimental data better than model 1. Figure 8 compares the RMSE data that were obtained in the verification of the models with and without breaking of the oxo-bonds.
3. Model 2 with mask bias. An additional calibration parameter, mask bias, is applied. The verification RMSE for CD exhibited an unreasonably high value for P60V27 compared to the other feature types for both models. The occurrence of a large RMSE of CD values at a single pitch is hard to explain by chemical effects in the photoresist. As a result, an additional calibration parameter, a mask bias only for P60V27 feature, is applied and a new calibration run was conducted. The results show a further improvement of the RMSE value for the CD and LWR (Fig. 8). This calibration run results in good accuracy for the CD and LWR (results shown in Table 3).

Table 3

Verification results for the different calibration procedures.

	Calibration RMSE		Verification RMSE
	CD (nm)	LWR (nm)	CD (nm)	LWR (nm)
Model 1	1.17	0.78	1.07	0.75
Model 2	1.15	0.58	1.06	0.47
Model 2 with mask bias	0.64	0.56	0.60	0.40

Fig. 8

RMSE for verification of the model 1 and model 2, in comparison with the calibration procedure with mask bias applied to P60V27 for model 2. (a) CD RMSE (nm) versus feature types and (b) LWR RMSE (nm) versus feature types.

The verification results of the three calibration runs are summarized in Table 3. The model approximates the experimental data with CD RMSE of 0.60 nm and LWR RMSE of 0.40 nm for the verification run. These results are comparable to the observed results of well-established CAR photoresist models.

The relative errors of all data points, shown in Fig. 9, were calculated from the verification errors normalized by the corresponding experimental CD or LWR values. The model can approximate the experimental data with 30% maximum error for LWR and 6% maximum error for the CD. Uncertainties in CD computation for profiles with higher roughness can lead to minor deviations on the final CD and LWR values. The parameter values for the “best” calibration solution selected from the Pareto front are summarized in Table 4.

Fig. 9

Relative RMSE for verification of the calibration parameters for model 2 with mask-bias for P60V27. The percentage value is calculated from the experimental data (a) CD and (b) LWR values. The colors of the dots represent the different feature types of the data points.

Table 4

Model parameters calibrated with experimental data, the ranges used for the calibration, and parameter values of a calibration result selected from the Pareto front.

Calibration parameters	Symbols	Range	Model 1	Model 2	Model 2 with MB
Absorption coefficient ( $μ m^{- 1}$ )	$α$	15.0 to 22.0	18.88	18.2	20.2
Electron blur length (nm)	$σ_{blur}$	1.5 to 4.0	3.41	2.86	3.28
Quantum yield	$Φ$	7.0 to 9.0	8.1	8.55	7.98
Critical number of sites	$L_{thr}$	2 to 5	4	4	3
Ligand substitution rate ( $s^{- 1}$ )	$K_{1}$	1.0 to 50.0	18.13	27.1	26.9
Solvent substitution rate ( $s^{- 1}$ )	$K_{2}$	0.0 to 35.0	19.77	23.99	13.49
Oxo-bond substitution rate ( $s^{- 1}$ )	$K_{3}$	0.0 to 35.0		28.28	13.63
Focus offset ( $μ m$ )		$- 0.08$ to $- 0.06$	$- 0.062$	$- 0.064$	$- 0.062$
CD metrology offset (nm)		$- 2.0$ to 2.0	0.03	0.36	$- 0.94$
LWR metrology offset (nm)		$- 1.0$ to 1.0	$- 0.62$	$- 0.89$	$- 0.97$

The Bossung data for verification of the calibration parameters for the smallest pitch (P36V18) and the largest pitch (P70V27) are shown in Fig. 10. For P70V27, the simulated LWR and CD errors vary from 0 to 1.3 nm and 0 to 1.0 nm, respectively. Similarly, for P36V18, simulated LWR errors vary from 0 to 2.0 nm and simulated CD errors vary from 0 to 1.14 nm. The experimental data for feature type P36V18 contain outlier data that could not be approximated by the model [marked by blue circles in Fig. 10(b)]. Excluding the outlier data, for P36V18, simulated LWR errors vary from 0 to 1.0 nm and simulated CD errors vary from 0 to 1.05 nm.

Fig. 10

Verification results for CD and LWR of the model with the experimental data points (model 2). Bossung data plots are shown for P36V18 [(a) and (b)] and P70V27 [(c) and (d)] feature types. The dots represent the experimental data values and the lines represent the simulated values. The red X represent the data points used for the calibration. The dose values are represented by color; red color represent small doses. The blue circles represent the outlier LWR data points.

4. Relation Between Photoresist Parameters and the Lithography Metrics

This section studies the impact of the parameters of the implemented model on important lithography metrics CD, LWR, dose-to-size, and exposure latitude. To perform a qualitative analysis of observed dependencies and to derive general tendencies, we extend previously published work³⁶ on the application of correlation analysis (CA) for photoresists parameters. The generation of appropriate datasets for CA requires a preprocessing to provide an efficient sampling of the data space and to avoid defects in the simulated photoresist profiles. To focus our investigations on the impact of the photoresist parameters and to include knowledge from previously published work on the scaling of lithography metrics,³⁷^,³⁸ some of the computed lithography metrics are normalized by factors that describe the impact of imaging parameters.

4.1.

Methodology

First, the parameter values of the datasets have to be randomly generated. For this purpose, parameter values are selected from the given bounds with Latin hypercube sampling (LHS) instead of random sampling to cover the entire parameter space with a minimum number of datasets.³⁹ Each dataset contains values of process settings (mask CD, focus, and pitch) and photoresist parameters ( $α$ , $σ_{blur}$ , $Φ$ , $L_{th}$ , $K_{1}$ , $K_{2}$ , $K_{3}$ , and th) summarized in Table 5. 600 samples are generated for the 11 parameters to ensure a good convergence for the CA. The pitch, unlike the other parameters, is not directly generated. Instead, the duty ratio is sampled; and the pitch is computed from the mask CD and the duty ratio.

Table 5

The lower and upper bounds for parameters sampled with LHS.

Parameters	Symbols	Lower bound	Upper bound
Focus (from best focus position) ( $μ m$ )		$- 0.65$	0.45
Mask CD (nm)		18.0	30.0
Duty ratio	$η$	2.0	2.5
Photoresist thickness (nm)	th	18.0	32.0
Absorption coefficient ( $μ m^{- 1}$ )	$α$	17.0	22.0
Quantum yield	$Φ$	7.0	12.0
Electron blur length (nm)	$σ_{blur}$	2.0	6.0
Critical number of sites	$L_{thr}$	2	3
Ligand substitution rate ( $s^{- 1}$ )	$K_{1}$	1.0	50.0
Solvent substitution rate ( $s^{- 1}$ )	$K_{2}$	0.0	35.0
Oxo-bond substitution rate ( $s^{- 1}$ )	$K_{3}$	0.0	35.0

The specifications of the lower and upper bounds used for LHS, summarized in Table 5, are determined by the effect of the parameter values on the final results of the simulation, to avoid defects on simulated profiles. In addition, the bounds are also limited to keep the monotonic relation of the parameter with the metrics because a non-monotonic data can suppress the correct results of the correlation. The relation of the LWR with the blur length is monotonic for a small range of blur length values. Careful selection of the parameter range and the application of the Spearman’s rank correlation (see below) helps to address the characteristics behavior of blur in our CA.

After the generation of parameters, simulations are performed for the computation of the metrics. First, the CD values, corresponding to the randomly generated process settings (mask CD, pitch, and focus), are computed with the calibrated parameters of model 2 from Table 4 with the nominal dose of $52 mJ / {cm}^{2}$ . In the following simulation with variable photoresist parameters, the CD values, obtained from the calibrated model 2 parameters, are used as target-CD for the computation of DtS and EL. LWR was extracted at the computed DtS. Finally, the CD was computed with the nominal dose of $52 mJ / {cm}^{2}$ .

The final step in the preprocessing of the data before CA addresses the simultaneous variation of process settings (mask CD, pitch, and focus) and photoresist parameters and their impact on the obtained results. In theory, it would be possible to vary only the photoresist parameters and apply the CA. However, the study of the photoresist response for fixed exposure conditions will not provide representative data for a wide range of applications. In order to include a wide range of process variations, the process settings are varied as well, but their dominating impact on certain lithography metrics has to be considered. This is done by appropriate normalization techniques and application of scaling rules for lithography metrics that have been described by several authors.³⁷^,³⁸^,⁴⁰ For example, it is known that the LWR decreases with increasing normalized image log-slope (NILS) and dose. Therefore, we normalize the variation of exposure conditions by application of appropriate scaling rules.

These impacts of the process settings are normalized in two steps.

1. The impact of the investigated process settings on the image quality and lithography metrics is typically quantified by image log-slope (ILS) and NILS. Here we employ analytical dependencies from Refs. 37, 38, and 40, to account for the impact of image quality and exposure dose on LWR and EL. Specifically, we normalize LWR according to
Eq. (6)
${LWR}_{n} = LWR * \sqrt{\frac{DtS}{{DtS}_{ref}}} * ILS,$
where ${DtS}_{ref}$ is the nominal dose ( $52 mJ / {cm}^{2}$ ) applied for computation of the target-CD.
The values of the EL are normalized to
Eq. (7)
${EL}_{n} = \frac{EL}{NILS} .$
In the case of CD, we use the change in CD ( $Δ C D$ ) instead of the absolute value, deviation of the simulated CD from the target-CD. For DtS, we use the data without any normalization.
2. The impact of some photoresist parameters on the metrics is dependent on the process settings, and these indirect contributions of the process settings should be normalized. For example, the effect of blur length on the metrics depends on the pitch. In order to remove the pitch dependency and treat only the blur length impact, the modulation transfer function (MTF), defined in Refs 38 and 40, is fed to the CA instead of the blur length values. MTF, defined as $MTF = e^{- 2 {(\frac{π σ_{blur}}{P})}^{2}}$ , is the Fourier transform of the Gaussian kernel (Eq. 4) that we applied for the approximation of the electron blur effect.⁴¹ However, the inverse relation $1 / MTF$ is applied because the increase in blur length values corresponds to a decrease in MTF values. Otherwise, the sign of the correlation coefficients for the blur length will be inverted.

Figure 11 demonstrates the impact of these postprocessing procedures. It can be seen that the application of the normalization increases the correlation between the photoresist parameters and the lithography metrics.

Fig. 11

Comparison for the impact of blur length ( $σ_{blur}$ ) variation on the LWR and quantum yield ( $Φ$ ) variation on the CD, before [(a) LWR versus $σ_{blur}$ and (b) CD versus $Φ$ ] and after [(c) LWR versus 1/MTF and (d) $Δ CD$ versus $Φ$ ] postprocessing. The red lines exhibit linear fits of the data before and after postprocessing. $r$ and $p$ represent correlation coefficient and $p$ -value, respectively.

Finally, the combined impact of the parameters on the lithography metrics has to be considered in the CA of the lithography metrics. To decorrelate the combined impacts of the parameters on the lithography metrics, we compute the semipartial correlation coefficients.³⁶ To consider the monotonic, but non-linear relation between some of the parameters and the lithography metrics, Spearman’s rank correlation coefficient is computed instead of Pearson’s linear correlation coefficients. To address both non-linear dependencies and decorrelate the impact of the photoresist parameters, we calculate semipartial rank correlation coefficients (SRCCs).³⁹^,⁴²

SRCC values are in the range from $- 1$ to 1. Negative SRCC values mean the increase in the photoresist parameter value leads to a reduction of the lithography metrics, whereas positive SRCC values mean the increase in the photoresist parameter value leads to an increase in the metrics. On the other hand, an absolute value of SRCC above 0.2 indicates that the photoresist parameter’s variation significantly impacts the metric for the defined parameter space. Otherwise, the impact is insignificant.

4.2.

Results of the Analysis

The results of our CA, shown in Fig. 12, provide several expected findings that are consistent with the literature. As shown in Fig. 12(a), the quantum yield ( $Φ$ ) has a significant effect on the reduction of the LWR. The increase in the blur length ( $σ_{blur}$ ) and the photoresist thickness (th) leads to an increase in LWR. The CD of lines increases with increasing absorption coefficient ( $α$ ) and quantum yield [Fig. 12(b)]. In turn, the increase in the CD for the increase in these parameters means a decrease in the DtS. EL is mainly affected by the blur length; the increase in the blur length leads to a reduction of the EL.

Fig. 12

The SRCCs for (a) LWR; (b) CD; (c) DtS; and (d) EL. The green color bars represent the development model parameters and the blue color bars represent the other remaining photoresist parameters.

The findings on the development parameters are less obvious. The development critical value ( $L_{th}$ ) exhibits a significant correlation with LWR and DtS and increases the CD of the lines. The correlation of the other development parameters with the process metrics is insignificant. The development critical value has a linear effect on the lithography metrics—its value decides the number of oxo-bonds per core that create the edges of the profiles. Its variation leads to a variation in the shift of the edge position. The other development process parameters impact the path of the development,¹⁴ but their impact on the CD and LWR is less seen in the final results of the lithography metrics variation. Nevertheless, the inclusion of these parameters and the corresponding path is important to obtain a fitness of the model calibration. The full understanding of the impact of these parameters on the process metrics needs further investigation.

Finally, we have to emphasize that the CA results are dependent on the parameter bounds chosen for the investigations. Consideration of different ranges of values can lead to deviations from the presented results.³⁶ Nevertheless, the observed tendencies especially for the significant impact of the development critical value on LWR are observed for all reasonable choices of parameters ranges.

5. Conclusion and Summary

A development model for organometallic photoresists that tracks the developer’s effect on the ligand and oxo-bond sites of the metal-oxide cores for negative-tone development was proposed and implemented in a complete imaging and resist process simulation flow. The model approximates the experimental Bossung data, of lines and spaces with different features and pitches, with CD RMSE of 0.60 nm and LWR RMSE of 0.4 nm. Notably, the application of multi-objective optimization and the detailed description of the dissolution process by the model resulted in good approximations of both CD and LWR simultaneously. The model also demonstrated that the development behavior of the investigated organometallic photoresist can be approximated based on the number of created oxo-bonds on the cores. The calibration results show that, for Inpria-YA photoresists, three ligands on the photoresist core should be substituted by solvent molecules before the photoresist core is assumed to be developed.

The correlation analysis confirms that both exposure and development photoresist parameters have a significant impact on the lithography metric. The most important impact of the development parameters is seen for the development critical value that impacts the LWR, CD, and DtS. The increase in the development critical value leads to a reduction in LWR and DtS, and an increase in the CD of the lines. The increase in quantum yield has a significant influence on the reduction of the LWR and the DtS. However, the quantum yield has no significant effect on the EL. The increase in the quantum yield also leads to an increase in the CD.

Application of the model to other feature types including contact holes and line ends could provide further insights on its predictivity. In addition, calibration of the models including the LER in addition to the LWR and CD data can improve the model and give better understanding of the development process. Further investigation and extensions of the photoresist parameter CA could help to separate and optimize the impact of exposure and development effects on the trade-off among resolution, LWR, and sensitivity.

6. Appendix A: Implementation of the Development Model

A discrete stochastic development model is implemented by simulating each reaction explicitly using the Gillespie algorithm.²⁸^,⁴³ The evolution of the development process through time can be tracked by computing the reaction of the solvent molecules with the sites of the metal-oxide cores that are in contact with the developer, based on Eq. (5).

The probability that a reaction occurs in time interval $[t, t + δ t]$ is given by $R_{i} * δ t$ , where $i$ is a reaction type and $R_{i}$ is a reaction rate. At the start of the simulation, the developer front contains only the cores at the top of the photoresist as the cores in contact with the developer constitute the developer front.²⁸ The reaction rates of the population balance equations in Eq. (5) are defined based on the number of ligand, oxo-bond, and solvent sites on the metal-oxide cores. These reactions rates determine not only the “stochastic state” of the developer front but also the next reaction. For our case, the state of the developer front is computed by the following equation:

Eq. (8)

R_{tot} = R_{L} + R_{S} + R_{B} R_{L} = K_{1} N_{L} R_{S} = K_{2} N_{S} R_{B} = K_{3} N_{B},

where

N_{L}

,

N_{B}

, and

N_{S}

are the total number of ligand sites, oxo-bond sites, and solvent sites on the developer front, respectively. These total number of sites are the sums of the respective number of sites on metal-oxide cores that are on the developer front. The cell for the next reaction and the reaction type are randomly selected based on reactivity rates,

R_{L}

,

R_{B}

, and

R_{S}

.

First, a core $i$ is randomly selected from the cores on the developer front. Then the procedure based on Von Neumann rejection from Ref. 28, instead of the original Gillespie algorithm procedure, is applied for acceptance of the randomly selected core. The Von Neumann rejection is adapted for our application as shown in the following equation:

Eq. (9)

if r_{0} \leq \frac{R_{L, i} + R_{S, i}}{\max_{i} {R_{L, i}, S_{L, i}} \forall i} accept if r_{0} > \frac{R_{L, i} + R_{S, i}}{\max_{i} {R_{L, i}, S_{L, i}} \forall i} reject,

where

r_{0}

is generated from a uniformly distributed random numbers.

For the acceptance of the randomly selected core, the sum of the reaction rates of the core, excluding the oxo-bond substitution reaction rate, is compared with the maximum of the ligand or the solvent substitution rate. This is because a metal-oxide core without an oxo-bond is expected to be developed—probability of reaction is 1—while a core that created oxo-bonds on all of the possible sites is resistant to the developer—probability of reaction is 0. As a result, the probability that the core reacts is computed from the solvent and ligand substitution rates compared to the maximum value of the ligand substitution rate.

If the selected core is accepted, the reaction type is randomly selected and the number of sites on the core ( $i$ ) is updated, as described in Eq. (10):

Eq. (10a)

if r_{1} < R_{L} / R_{tot} ligand substitution N_{L, i} (t + δ t) = N_{L, i} (t) - 1 N_{S, i} (t + δ t) = N_{S, i} (t) + 1,

Eq. (10b)

if r_{1} > R_{L} / R_{tot} and r_{1} < (R_{L} + R_{B}) / R_{tot} oxo-bond substitution N_{S, i} (t + δ t) = N_{S, i} (t) + 1 N_{L, i} (t + δ t) = N_{L, i} (t) - 1,

Eq. (10c)

if r_{1} > (R_{L} + R_{B}) / R_{tot} solvent substitution N_{B, i} (t + δ t) = N_{B, i} (t) - 1 N_{S, i} (t + δ t) = N_{S, i} (t) + 1,

where

i

represents the selected core on the developer front.

r_{1}

is generated from the uniformly distributed random numbers and

N_{j, i}

is the number of sites of type

j

(ligand or oxo-bond or solvent sites) on core

i

selected from the developer front.

The time step for the reaction, in the Gillespie algorithm, is calculated from the state of the system using Eq. (11). Then the time will advance from $t$ to $t + δ t$ :²⁸

Eq. (11)

δ t = - \ln (r_{2}) / R_{tot},

where

r_{2}

is generated from a uniformly distributed random numbers.

If the number of solvent substituted sites per core is larger than the critical number of sites ( $L_{th}$ ), and all the oxo-bonds on the core are broken, the core is developed. If the core meets the development criteria, the developer front is updated by removing the core and inserting the neighbor cores to the developer front.

These steps—random selection of the core and the reaction type, update of the sites on the core, and update of the developer front—are executed for each time step. These processes are iterated until the development time is reached.

Similar to the previous process steps of exposure and condensation, a periodic boundary condition is applied for $x$ and $y$ directions.

6.1.

Computational Optimization Approach

The exact stochastic simulation approach for the Gillespie algorithm (the naïve implementation), discussed above, is computationally intensive.⁴⁴ Randomly selecting a single core and the reaction type for the selected core requires several iterations until all the cores on the developer front are visited. Also due to the large $R_{tot}$ , because the developer front contains several cores, the update time ( $δ t$ ) is small. Until the development time is reached, the simulation requires several iterations. Furthermore, after each iteration, the reaction rate is recalculated, and the developer front is updated. The developer front update is also a computationally intensive process. As a result, especially for simulation of profiles with long lines, as it is required for proper computation of the LWR, the development model becomes inefficient. Therefore, to reduce the computation time, optimization of the processes is required.

In order to reduce the number of updates, the tau-leaping approximation procedure⁴⁴^,⁴⁵ is included the implementation. In the tau-leap method, the number of sites updated at a time is randomly selected based on Poisson probability, instead of executing a single reaction for each iteration. The generated random number of reactions is used to update the sites on the core based on the selected reaction type. Due to the Poisson distribution’s unbounded nature, the maximum numbers of reactions are limited to the number of available sites on the core corresponding to the reaction type, to avoid negative values for the sites on the core.⁴⁶ Note that, in our implementation, a single reaction is executed if the Poisson generated number of reactions is zero.

However, the tau-leaping method-based optimization is still insufficient for the model to be fast enough to be applied for calibration. For a single metal-oxide core, six reaction sites are available with three possible reactions types. Each core requires several reaction steps—at a minimum, reactions equal to the critical value—before it is developed. In addition, due to the large size of the developer front, several iterations are still required to change the state of the developer front. As a result, an aggressive optimization approach is required.

The development process is a weakly coupled reaction network.⁴⁷^,⁴⁸ This means that the development behavior of a core on the developer front is affected only by the immediate neighbor cores, irrespective of the developer front’s size. In addition, a reaction on a metal-oxide core affects the neighbor cores or changes the state of the developer front significantly if the reaction leads to the development of the core. Due to these reasons, the metal-oxide cores on the developer front can be treated to be spatially independent. In other words, the developer front can be subdivided into smaller systems where the reactions types are probabilistically selected for each subsystem.⁴⁹ As a result, the developer front is partitioned into subsystems of single cores. In this assumption, each core can be treated as a separate system, and the reaction type is sampled independently for each core. Therefore, all the cores can simultaneously undergo substitution by a solvent or by a ligand.

In the current implementation, the cores on the developer front are visited sequentially, and depending on the reaction rates of the corresponding core, the reaction type is randomly selected. The next reaction’s time step is computed separately for each core from its total reaction rate based on Eq. (11). The computed update time is multiplied by the Poisson generated number of reactions, to account for the number of reactions generated based on Poisson probability. Then based on the selected reaction type, the sites on the current core are updated. After all the cores are visited, the reaction rates are recalculated, and the developer front is updated. These processes are iterated until the development time is reached. The pseudo-code for the optimized procedure is presented in Algorithm 2.

Algorithm 2

Development algorithm

1: start

2: initialize the developer front

3: while $t < t_{dev}$ do

4: calculate reaction rates for each core on the front

5: for $k \leftarrow 1$ to number of cores on the developer front do

6: generate random numbers,

r_{0}

,

r_{1}

,

r_{2}

7: compute the update time according to Eq. (11)

8: accept or reject the reaction on the current core according to Eq. (9)

9: choose the reaction type and update sites on the core according to Eq. (10)

10: end for

11: Update the developer front

12: end while

13: end

Acknowledgments

The authors would like to express their gratitude to Inpria for providing the experimental data for the calibration of our model.

References

1.

J. J. Biafore et al., “Resist pattern prediction at EUV,” Proc. SPIE, 7636 76360R (2010). https://doi.org/10.1117/12.846535 PSISDG 0277-786X Google Scholar

2.

R. Fallica et al., “Dynamic absorption coefficients of chemically amplified resists and nonchemically amplified resists at extreme ultraviolet,” J. Micro/Nanolithogr. MEMS MOEMS, 15 (3), 033506 (2016). https://doi.org/10.1117/1.JMM.15.3.033506 Google Scholar

3.

A. Grenville et al., “Integrated FAB process for metal oxide EUV photoresist,” Proc. SPIE, 9425 94250S (2015). https://doi.org/10.1117/12.2086006 PSISDG 0277-786X Google Scholar

4.

S. T. Meyers et al., “Organotin oxide hydroxide patterning compositions, precursors, and patterning,” (2016). Google Scholar

5.

A. Robinson and R. Lawson, Materials and Processes for Next Generation Lithography, Elsevier Science(2016). Google Scholar

6.

W. D. Hinsberg and S. Meyers, “A numeric model for the imaging mechanism of metal oxide EUV resists,” Proc. SPIE, 10146 1014604 (2017). https://doi.org/10.1117/12.2260265 PSISDG 0277-786X Google Scholar

7.

A. V. Pret et al., “Modeling and simulation of low-energy electron scattering in organic and inorganic EUV photoresists,” Proc. SPIE, 10146 1014609 (2017). https://doi.org/10.1117/12.2261434 PSISDG 0277-786X Google Scholar

8.

R. Maas et al., “Stochastics in extreme ultraviolet lithography: investigating the role of microscopic resist properties for metal-oxide-based resists,” J. Micro/Nanolithogr. MEMS MOEMS, 17 (4), (2018). https://doi.org/10.1117/1.JMM.17.4.041003 Google Scholar

9.

A. Sasaki, M. Ishino and M. Nishikino, “An estimation of line width roughness of photoresists due to photon shot noise for extreme ultraviolet lithography using the percolation model,” Jpn. J. Appl. Phys., 58 (5), 055002 (2019). https://doi.org/10.7567/1347-4065/ab06bc Google Scholar

10.

U. Welling et al., “Modeling of CAR alternatives for EUV lithography,” in 35th Eur. Mask and Lithogr. Conf. (EML19), (2019). Google Scholar

11.

C. D. Needham et al., “Calibration of a mox-specific EUV photoresist lithography model,” Proc. SPIE, 11323 113230G (2020). https://doi.org/10.1117/12.2552151 PSISDG 0277-786X Google Scholar

12.

T. Fuehner et al., “Dr. LiTHO: a development and research lithography simulator,” Proc. SPIE, 6520 65203F (2007). https://doi.org/10.1117/12.709535 PSISDG 0277-786X Google Scholar

13.

J. J. Biafore et al., “Statistical simulation of photoresists at EUV and ArF,” Proc. SPIE, 7273 727343 (2009). https://doi.org/10.1117/12.813551 PSISDG 0277-786X Google Scholar

14.

C. Mack, “Stochastic approach to modeling photoresist development,” J. Vac. Sci. Technol. B, 27 (3), (2009). https://doi.org/10.1116/1.3117346 JVTBD9 1071-1023 Google Scholar

15.

T. Kozawa, J. J. Santillan and T. Itani, “Analysis of stochastic effect in line-and-space resist patterns fabricated by extreme ultraviolet lithography,” Appl. Phys. Express, 6 (2), 026502 (2013). https://doi.org/10.7567/APEX.6.026502 APEPC4 1882-0778 Google Scholar

16.

C. A. Mack et al., “Stochastic exposure kinetics of extreme ultraviolet photoresists: simulation study,” J. Micro/Nanolithogr. MEMS MOEMS, 10 (3), 033019 (2011). https://doi.org/10.1117/1.3631753 Google Scholar

17.

B. Cardineau et al., “EUV resists based on tin-oxo clusters,” Proc. SPIE, 9051 90511B (2014). https://doi.org/10.1117/12.2046536 PSISDG 0277-786X Google Scholar

18.

J. Passarelli et al., “Organometallic carboxylate resists for extreme ultraviolet with high sensitivity,” J. Micro/Nanolithogr. MEMS MOEMS, 14 (4), 043503 (2015). https://doi.org/10.1117/1.JMM.14.4.043503 Google Scholar

19.

A. Aharony and D. Stauffer, Introduction To Percolation Theory, Taylor & Francis(2003). Google Scholar

20.

O. Yildirim et al., “Improvements in resist performance towards EUV HVM,” Proc. SPIE, 101430Q (2017). https://doi.org/10.1117/12.2257415 PSISDG 0277-786X Google Scholar

21.

A. Erdmann et al., “Attenuated phase shift mask for extreme ultraviolet: can they mitigate three-dimensional mask effects?,” J. Micro/Nanolithogr. MEMS MOEMS, 18 (1), 011005 (2019). https://doi.org/10.1117/1.JMM.18.1.011005 Google Scholar

22.

Y. Ma, J. Shin and F. Cerrina, “Line edge roughness and photoresist percolation development model,” J. Vac. Sci. Technol. B, 21 (1), 112 (2002). https://doi.org/10.1116/1.1534572 JVTBD9 1071-1023 Google Scholar

23.

C.-C. Yeh et al., “Chemical and structural investigation of zinc-oxo cluster photoresists for DUV lithography,” J. Mat. Chem. C, 5 (10), 2611 –2619 (2017). https://doi.org/10.1039/C6TC05201K Google Scholar

24.

S. Chakrabarty et al., “Oxide nanoparticle EUV resists: toward understanding the mechanism of positive and negative tone patterning,” Proc. SPIE, 8679 867906 (2013). https://doi.org/10.1117/12.2011490 PSISDG 0277-786X Google Scholar

25.

S. Chakrabarty et al., “Increasing sensitivity of oxide nanoparticle photoresists,” Proc. SPIE, 9048 90481C (2014). https://doi.org/10.1117/12.2046555 PSISDG 0277-786X Google Scholar

26.

C. Ober et al., “New developments in ligand-stabilized metal oxide nanoparticle photoresists for EUV lithography,” Proc. SPIE, 9422 942207 (2015). https://doi.org/10.1117/12.2086488 PSISDG 0277-786X Google Scholar

27.

J. Jiang et al., “Systematic study of ligand structures of metal oxide EUV nanoparticle photoresists,” Proc. SPIE, 9422 942222 (2015). https://doi.org/10.1117/12.2084896 PSISDG 0277-786X Google Scholar

28.

A. Philippou, “Mesoscopic and continous simulation of photoresist development simulation,” Technical University of Berlin, Berlin, (2008). Google Scholar

29.

L. W. Flanagin, V. K. Singh and C. G. Willson, “Surface roughness development during photoresist dissolution,” J. Vac. Sci. Technol. B, 17 (4), 1371 (1999). https://doi.org/10.1116/1.590763 JVTBD9 1071-1023 Google Scholar

30.

E. F. Caldin and H. P. Bennetto, “The role of solvent structure in ligand substitution and solvent exchange at some divalent transition-metal cations,” J. Sol. Chem., 2 (2), 217 –238 (1973). https://doi.org/10.1007/BF00651974 Google Scholar

31.

J. M. Santos et al., “NMR kinetic study of ligand exchange on a (μ-oxo)bis(μ-acetato)diruthenium(iii) complex. A special kinetic case of introduction of a more labile ligand,” Can. J. Chem., 75 (6), 890 –898 (1997). https://doi.org/10.1139/v97-107 CJCHAG 0008-4042 Google Scholar

32.

P. De Bisschop, “Stochastic effects in EUV lithography: random, local CD variability, and printing failures,” J. Micro/Nanolithogr. MEMS MOEMS, 16 (4), 041013 (2017). https://doi.org/10.1117/1.JMM.16.4.041013 Google Scholar

33.

A. V. Pret et al., “Characterizing and modeling electrical response to light for metal-based EUV photoresists,” Proc. SPIE, 9779 977906 (2016). https://doi.org/10.1117/12.2219736 PSISDG 0277-786X Google Scholar

34.

T. Fuehner et al., Genetic Algorithms to Improve Mask and Illumination Geometries in Lithographic Imaging Systems, 208 –218 Springer Berlin Heidelberg, Berlin, Heidelberg (2004). Google Scholar

35.

N. Thakur et al., “Stability studies on a sensitive EUV photoresist based on zinc metal oxoclusters,” J. Micro/Nanolithogr. MEMS MOEMS, 18 (4), 043504 (2019). https://doi.org/10.1117/1.JMM.18.4.043504 Google Scholar

36.

B. Tollkuhn et al., “Correlation analysis: a fast and reliable method for a better understanding of simulation models in optical lithography,” Proc. SPIE, 5755 37 –47 (2005). https://doi.org/10.1117/12.599390 PSISDG 0277-786X Google Scholar

37.

D. van Steenwinckel et al., “Novel method for characterizing resist performance,” J. Micro/Nanolithogr. MEMS MOEMS, 7 (2), 023002 (2008). https://doi.org/10.1117/1.2909204 Google Scholar

38.

J. G. Santaclara et al., “One metric to rule them all: new k4 definition for photoresist characterization,” Proc. SPIE, 11323 321 –330 (2020). https://doi.org/10.1117/12.2554493 PSISDG 0277-786X Google Scholar

39.

B. Gomero, “Latin hypercube sampling and partial rank correlation coefficient analysis applied to an optimal control problem,” University of Tennessee, (2012). Google Scholar

40.

C. A. Mack, “Reducing roughness in extreme ultraviolet lithography,” J. Micro/Nanolithogr. MEMS MOEMS, 17 (4), 041006 (2018). https://doi.org/10.1117/1.JMM.17.4.041006 Google Scholar

41.

T. A. Brunner et al., “Impact of resist blur on MEF, OPC, and CD control,” Proc. SPIE, 5377 141 –149 (2004). https://doi.org/10.1117/12.537472 PSISDG 0277-786X Google Scholar

42.

P. Schober, C. Boer and L. A. Schwarte, “Correlation coefficients: appropriate use and interpretation,” Anesth. Analg., 126 (5), 1763 –1768 (2018). https://doi.org/10.1213/ANE.0000000000002864 Google Scholar

43.

D. T. Gillespie, “Exact stochastic simulation of coupled chemical reactions,” J. Phys. Chem., 81 (25), 2340 –2361 (1977). https://doi.org/10.1021/j100540a008 JPCHAX 0022-3654 Google Scholar

44.

D. T. Gillespie, “Approximate accelerated stochastic simulation of chemically reacting systems,” J. Chem. Phys., 115 (4), 1716 –1733 (2001). https://doi.org/10.1063/1.1378322 JCPSA6 0021-9606 Google Scholar

45.

Y. Cao, D. T. Gillespie and L. R. Petzold, “Efficient step size selection for the tau-leaping simulation method,” J. Chem. Phys., 124 (4), 044109 (2006). https://doi.org/10.1063/1.2159468 JCPSA6 0021-9606 Google Scholar

46.

Y. Cao and D. C. Samuels, “Discrete stochastic simulation methods for chemically reacting systems,” Methods Enzymol., 454 115 –140 (2009). https://doi.org/10.1016/S0076-6879(08)03805-6 MENZAU 0076-6879 Google Scholar

47.

R. Ramaswamy, N. González-Segredo and I. F. Sbalzarini, “A new class of highly efficient exact stochastic simulation algorithms for chemical reaction networks,” J. Chem. Phys., 130 (24), 244104 (2009). https://doi.org/10.1063/1.3154624 JCPSA6 0021-9606 Google Scholar

48.

R. Ramaswamy and I. F. Sbalzarini, “A partial-propensity variant of the composition-rejection stochastic simulation algorithm for chemical reaction networks,” J. Chem. Phys., 132 (4), 044102 (2010). https://doi.org/10.1063/1.3297948 JCPSA6 0021-9606 Google Scholar

49.

J. Elf and M. Ehrenberg, “Spontaneous separation of bi-stable biochemical systems into spatial domains of opposite phases,” Syst. Biol., 1 (2), 230 –236 (2004). https://doi.org/10.1049/sb:20045021 SYBIER 1076-836X Google Scholar

Biography

Zelalem Belete is a doctoral candidate at the University of Erlangen-Nuremberg. He received his MSc degree in microelectronics and his MSc degree in advanced optical technologies from Addis Ababa University in 2012 and the University of Erlangen-Nuremberg in 2015, respectively. His current research interests include photoresist modeling, optics, lithography, physical modeling, stochastic simulations, and topography simulations.

Andreas Erdmann is the head of the Fraunhofer IISB Computational Lithography and Optics Group and teaches as “Privatdozent” at the University of Erlangen. He has more than 20 years of experience in optical and EUV lithography. He chaired SPIE conferences on “Optical Microlithography” and “Optical Design” and is an organizer of the International Fraunhofer Lithography Simulation Workshop. He contributed to the development of several advanced lithography simulators including Dr. LiTHO. He is a fellow of SPIE.

Biographies of the other authors are not available.

CC BY: © The Authors. Published by SPIE under a Creative Commons Attribution 4.0 Unported License. Distribution or reproduction of this work in whole or in part requires full attribution of the original publication, including its DOI.

Citation Download Citation

Zelalem Belete, Peter De Bisschop, Ulrich Welling, and Andreas Erdmann "Stochastic simulation and calibration of organometallic photoresists for extreme ultraviolet lithography," Journal of Micro/Nanopatterning, Materials, and Metrology 20(1), 014801 (4 March 2021). https://doi.org/10.1117/1.JMM.20.1.014801

Received: 16 November 2020; Accepted: 2 February 2021; Published: 4 March 2021

Access the abstract

JOURNAL ARTICLE
22 PAGES

DOWNLOAD PAPER SAVE TO MY LIBRARY

GET CITATION

CITATIONS

Cited by 6 scholarly publications.

Explore citations on Lens.org

KEYWORDS

Photoresist materials

Photoresist developing

Calibration

Line width roughness

Stochastic processes

Data modeling

Photons

1.

Introduction

2.

Modeling Procedure

2.1.

Exposure

Eq. (1)

Eq. (2)

Eq. (3)

Eq. (4)

Fig. 1

2.2.

Condensation Reaction

Algorithm 1

Fig. 2

2.3.

Development

Fig. 3

Eq. (5a)

Eq. (5b)

Fig. 4

2.4.

Illustrations of the Complete Process Flow

Fig. 5

Fig. 6

Table 1

3.

Calibration of Model with Experimental Data

Fig. 7

3.1.

Simulation Parameters and Process Conditions

Table 2

3.2.

Model Options, Calibration Results, and Verification of the Models

Table 3

Fig. 8

Fig. 9

Table 4

Fig. 10

4.

Relation Between Photoresist Parameters and the Lithography Metrics

4.1.

Methodology

Table 5

Eq. (6)

Eq. (7)

Fig. 11

4.2.

Results of the Analysis

Fig. 12

5.

Conclusion and Summary

6.

Appendix A: Implementation of the Development Model

Eq. (8)

Eq. (9)

Eq. (10a)

Eq. (10b)

Eq. (10c)

Eq. (11)

6.1.

Computational Optimization Approach

Algorithm 2

Acknowledgments

References

Biography

Show All Keywords

Keywords/Phrases

Search In:

Publication Years