Modeling of Reactive Sputtering—History and Development

This work critically reviews the evolution of reactive sputtering modeling that has taken place over the last 50 years. The review summarizes the main features of the deposition of simple metal compound films (nitrides, oxides, oxynitrides, carbides, etc.) that were experimentally found by different researchers. The above features include significant non-linearity and hysteresis. At the beginning of the 1970s, specific chemisorption models were proposed. These models were based on the assumption that a compound film was formed on the target due to chemisorption. Their development led to the appearance of the general isothermal chemisorption model, which was supplemented by the processes on the surfaces of the vacuum chamber wall and the substrate. The model has undergone numerous transformations for application to various problems of reactive sputtering. At the next step in the development of modeling, the reactive sputtering deposition (RSD) model was proposed, which was based on the implantation of reactive gas molecules into the target, bulk chemical reaction, chemisorption, and the “knock-on effect”. Another direction of the modeling development is represented by the nonisothermal physicochemical model, in which the Langmuir isotherm and the law of mass action are used. Various modifications of this model allowed describing reactive sputtering processes in more complex cases when the sputtering unit included a hot target or a sandwich one.


Introduction
Films of simple compounds of metals (oxides, nitrides, carbides, oxynitrides, etc.) and their solid solutions are of great interest in many fields of technology. This is due to the fact that they exhibit many properties (semiconductor, ferromagnetic, ferroelectric, electrochromic, photochromic, etc.) that open up new opportunities for innovation. The attention of scientists and engineers to such films is also associated with new areas of their application, which include ecology, medicine, and alternative energy.
For the deposition of these films, there are widely used methods summarized by the term "reactive sputtering" [1,2], which, according to the author [2], first appeared in [3]. Let us denote here the metal, the reactive gas, and their compound with stoichiometric coefficients m and n by M, X 2, and M m X n , respectively. Regardless of the sputtering system type, several methods are applied to deposit M m X n films [4]:

•
Reactive sputtering of the M-target in the Ar + X 2 environment. The method produces a film of stoichiometric composition; • Sputtering of the M m X n target in the Ar + X 2 environment. The method also leads to a stoichiometric film; • Sputtering of the M m X n target in the Ar environment. Due to the preferential sputtering of the light component X 2 , the surface layer of the target is reduced from M m X n to MX x , which allows the depositing of films of intermediate composition.
In early experimental works, the influence of the reactive gas concentration in the gas mixture or its partial pressure on the film growth rate, discharge voltage, and film

•
The target is free from reaction products; • Metal is sputtered from the target surface. The sputtering rate is high, and the ioninduced electron emission yield is not high; • The reactive gas adsorbs on the target, the inner surfaces of the vacuum chamber, on which atoms of the sputtered metal are deposited, and is pumped out; • Current-voltage (I-V) characteristic of the discharge corresponds to the discharge in pure argon.
In the reactive mode: • The target surface is completely covered by the M m X n compound; • M m X n molecules are sputtered from the target surface. The sputtering rate is low, and the ion-induced electron emission yield is high; • The target consumes an insignificant part of the reactive gas to maintain the steady state, the rest of the gas is pumped out by a vacuum pump; • The I-V characteristics of the discharge differ from that of the metallic mode-the direction and magnitude of the changes depend on the material of a target and the type of reactive gas.
Getting to the point of instability, the process of reactive sputtering spontaneously passes into a new steady state. In a large number of experiments, when changing Q 0 or the discharge current (power), the hysteresis effect was found. It was observed during the deposition of oxide films TiO 2 [25,26], Al 2 O 3 [27,28], Cd 2 SnO 4 [29,30], ZnO [31], Ta 2 O 5 [32], Ta 2 O 5 -TiO 2 [33] and nitrides AlN [34] and TiN [22,24,35]. The essence of the hysteresis resides in the following: there are two points of instability in the dependence of the X 2 reactive gas partial pressure on Q 0 at I = const (Figure 1a, curve 1). With an increase in Q 0 from zero, the target remains in the metallic mode up to point A. With this value of Q 0 , the target avalanches into the reactive mode. A return to the metallic mode with a decrease in Q 0 occurs at another point C. Curve 2 in Figure 1a reflects the corresponding dependence when the magnetron is turned off.
A similar nonlinearity was found in processes in which the discharge current is controlled at Q 0 = const (Figure 1b). In this case, as I increases from zero, the target remains in the reactive mode up to point A. At this value of the current, the target avalanches into the metallic mode. A return to the reactive mode with a decrease in I occurs at another point, C. deposition of oxide films TiO2 [25,26], Al2O3 [27,28], Cd2SnO4 [29,30], ZnO [31], Ta2O5 [32], Ta2O5-TiO2 [33] and nitrides AlN [34] and TiN [22,24,35]. The essence of the hysteresis resides in the following: there are two points of instability in the dependence of the X2 reactive gas partial pressure on Q0 at I = const (Figure 1a, curve 1). With an increase in Q0 from zero, the target remains in the metallic mode up to point A. With this value of Q0, the target avalanches into the reactive mode. A return to the metallic mode with a decrease in Q0 occurs at another point C. Curve 2 in Figure 1a reflects the corresponding dependence when the magnetron is turned off. Thus, numerous experiments have shown, firstly, that the main independent variables of the reactive sputtering process are the reactive gas flow Q 0 and the discharge current I (the power released on the target). Secondly, the curves of the dependence of the partial pressure of the reactive gas p on Q 0 ( Figure 1a) at I = const contain several sections:

•
In the initial section 0A on curve 1, the value of p is close to zero (the total pressure of the gas mixture is equal to the partial pressure of argon); • As Q 0 increases to point A, a noticeable change in p is not observed; • When point A is reached, the process avalanches into a new steady state (point B); • With a further increase in Q 0 , the pressure p grows in proportion to Q 0 , but for any Q 0 it is always less than the value measured in the absence of a discharge (curve 2 in Figure 1a); • As Q 0 decreases to point C, a proportional decrease in p is observed; • When point C is reached, the process returns like an avalanche to the initial steady state (point D).
When controlling the process of reactive sputtering by incoming flow Q 0 or discharge current I, it is possible to deposit the M m X n film only in the reactive mode. The points in Figure 1a, corresponding to this mode, are located on the line connecting points B and C. In Figure 1b, similar points are on the line connecting points A and D. Points A and C in Figure 1 are connected by a dashed line, which can only be reached by the sputtering process under specially created conditions. These conditions are described in [20,36]. The hysteresis in them was overcome by changing the controlled variable and introducing negative feedback. When sputtering a titanium target in the Ar + N 2 environment, instead of controlling the N 2 incoming flow, an adjustable partial pressure system was used. As a result, continuous dependencies were obtained that relate the partial pressure of N 2 to its incoming flow. The curves were N-shaped with areas of negative slope. This area was not available when controlling the incoming flow of nitrogen. Similar dependences were obtained by sputtering an aluminum target in the Ar + N 2 environment with discharge voltage control [34,37]. The presented results show that in systems with negative feedback on reactive gas pressure or discharge voltage, the sputtering process can be maintained at any point of the sections AC in Figure 1.
A brief review of the experimental results, indicating the features of reactive sputtering, showed the complexity of these processes. Therefore, it was difficult to choose the film deposition mode, especially at the initial stage of technology development. For a more detailed study of the processes of reactive sputtering, many researchers worked on the creation of their models.
Today it is obvious that the key processes in reactive sputtering occur on the target. On its surface, the processes of formation and sputtering of the M m X n film compete, and depending on the ratio of their rates, the target, as indicated earlier, can be in one of two Materials 2023, 16, 3258 4 of 50 steady states. The state of the target is quantitatively described using the relative fraction of its surface θ t covered by the M m X n film.
In the most general form, the kinetic equation for the target can be written as: With the help of the first term on the right side of the written equation, a set of processes that increase the value of θ t is taken into account. The second term describes the decrease in θ t .
A description of the history and development of these works is the purpose of this article. Let us first consider the early models, which we will call "specific" (Section 2). Next, let us turn our attention to publications presenting the Berg isothermal chemisorption models (Section 3), the similar Depla models (Section 4), and the Barybin nonisothermal physicochemical models (Section 5).

Specific Isothermal Models
All the varieties of known models of the reactive sputtering process are actually based on two assumptions: • Two processes compete on the excited target surface: the formation of a thin layer of metal compound with reactive gas and sputtering of this layer by accelerated argon ions; • On the substrate and walls of the vacuum chamber, the sputtering target material is deposited, and the reactive gas molecules are chemisorbed.
The physical model of processes on the target surface was first discussed in 1973 [15]. In 1975 [21,38], kinetic equations for the target based on chemisorption appeared.
In early works [19,21,22,26,29,37,[39][40][41], when constructing a reactive sputtering model, the partial pressure was taken as an independent variable p of the reactive gas and took into account only the processes occurring on the target without taking into account its temperature. The difference between the models of different authors was not fundamental. The main unifying element of all the noted publications was the kinetic equation for the target surface, which in the steady mode took the form: where N is the density of adsorption centers occupied by the reactive gas; α is the coefficient of adsorption of reactive gas molecules to the metal surface of the target; p eff is the effective pressure during sputtering; m 0 is the mass of a gas molecule; k is the Boltzmann constant; T is the absolute temperature; θ t = N/N t is the fraction of the target surface covered by the compound film; N t is the density of adsorption centers on the target; S C is the M m X n film sputtering yield; j is the discharge current density on the target. The effective pressure in (1) was defined as: where a is coefficient converting the number of sputtered atoms into the rate of gas adsorption; S p is the pumping speed of the chamber; A t is the target area; S M is the metal sputtering yield. The first term in (1) determines the growth rate of the M m X n compound, the factor p eff / √ 2πm 0 kT, according to the kinetic theory of gases, sets the density of the reactive gas molecules flowing to the target. The second term determines the sputtering rate of the target surface covered by the M m X n film by argon ions. The second term in expression (2) takes into account the decrease in pressure due to the adsorption of the reactive gas on the chamber walls. Solution (1) with respect to (1 − θ t ) has the form: where A, B, and C are the coefficients depending on the sputtering process parameters.
In [29], Equation (1) was written by its authors using the M m X n film thickness on the x target surface. The unbalance in the sputtering system, at which the proportion of the target surface occupied by the compound begins to increase, is described by the equation where (dx/dt) ch is the compound film formation rate; (dx/dt) sp is the sputtering rate of the M m X n layer. According to [15], the film formation rate is expressed in exponential form: where A(p) is an increasing function of the reactive gas pressure; ρ is the compound density; x 0 is a temperature-dependent parameter. The change in the partial pressure of oxygen upon unbalance is defined as: where A and B are the constants, depending on system parameters. In [37], for the development of the model described by Equation (1), two mechanisms were proposed that form a compound of metal and reactive gas atoms on the target surface: • Chemisorption of neutral reactive gas molecules (which can occur without a glow discharge); • Covering the target with ions and atoms of reactive gas, which are activated by a glow discharge. This process is called by the authors "ion plating".
Kinetic Equation (1), in this case, took the form where F is the density of the flow of reactive gas molecules to the target, which is defined by the Hertz-Knudsen equation; f (p/p t ) is the proportion of reactive gas positive ions in the discharge current; p t is the total pressure. The rest of the parameters have different values on the clean and coated parts of the M m X n target. In this regard, their values, when averaged over the entire surface of the target, will depend on θ t ; α i (θ t ) is the adsorption coefficient for reactive gas ions incident on the target surface; γ(θ t ) is the average ioninduced secondary electron emission yield for the total target surface; S C (θ t ) is the average sputtering yield of reactive gas molecules on the target. The first term on the right side of Equation (3) takes into account the chemisorption of neutral reactive gas molecules on the target surface, the second term-ion deposition associated with the discharge current, and the third term-the sputtering of reactive gas molecules from the target surface. The considered model, according to the authors, allows for determining the composition of the compound film from the intensity of the metal line in the discharge emission spectrum and other process parameters.
Summing up, we note that all models of reactive sputtering, which are described by Equation (1), combine several assumptions:

•
The partial pressure of the reactive gas and the film deposition rate are determined by the degree of coverage of the target surface with the M m X n film; • Chemisorption was adopted as the mechanism of the M m X n film formation on the target surface; • The process is considered isothermal in the sense that the gas environment and all internal surfaces of the vacuum chamber have the same temperature; •  The model equations describe only the process on the target; there is no separate  equation for the chamber wall and the balance equation for gas particles. In [42], when modeling reactive sputtering, adsorption of the reactive gas on the wall of the vacuum chamber was physically correctly taken into account. Remaining within the framework of the isothermal model, the author expressed the formation of the compound on the target through adsorption, writing the balance equation for gas flows as follows: In Equation (4), the reactive gas flows on the target surface due to adsorption and sputtering are denoted by Q t + and Q t − , respectively; on the chamber wall-Q w ; due to the operation of the vacuum pump-Q p . Next, the values of the effective pumping speed for each sorbing surface were introduced. The resulting equations allowed the author to construct kinetic curves and describe the hysteresis effect during the deposition of a TiO 2 film. However, the kinetic dependences presented in [42] have inflection points; this does not correspond to the experimental results presented, for example, in [43]. Another feature of the result obtained in [42] is the possibility of a steady state with an intermediate value 0 < θ t < 1. This contradicts all known experimental results.
In [44,45], in order to develop the model proposed in [21], the change in the partial pressure of the reactive gas was expressed in terms of gas flows: where V is the vacuum chamber volume; F i is the density of the reactive gas flow on the i-th surface (target t, substrate s, chamber wall w); A i is the area of the i-th surface (i = t, s, w) consuming the reactive gas, the flow of which is equal to: where θ i is the degree of coverage of the i-th surface by the compound. As a result, the balance equation for gas flows for the steady state appeared: where the term pS p /kT defines the flow that passes through the pump. On the basis of expressions (5)-(7) and the low-pressure gas discharge model, the kinetics of Al sputtering in the Ar + O 2 environment was studied in [44,45]. Thus, in [37,42,44,45], the reactive sputtering model was further developed. In its analytical description, equations for gas flows to all surfaces of the chamber and the balance equation for gas particles appeared.
Further development of reactive sputtering modeling was carried out in a series of papers [35,36,[46][47][48][49][50][51][52]. The authors most consistently developed a model that we will call "general", in contrast to the specific models of the previous section.

Initial Model
The model that the authors described in [46] was later called the Berg model. It should be noted that the group of Professor S. Berg has been working hard on the development of the reactive sputtering model over the past thirty-five years and has published more than a dozen articles on this topic. Further, when discussing different models, we use notations that partially differ from those adopted by the authors of the referred works and used in Section 1 of this article. This was done in order to unify the mathematical descriptions of different models.
In the first work of 1987 [46], the ideas expressed in the mentioned works were developed for the reactive sputtering of a one-component metal target M in the environment of one reactive gas X 2 . Describing the physical model, the authors proposed a number of constraints: 1.
The model should take into account the processes occurring on the surfaces of the target and the wall of the vacuum chamber with areas A t and A w , respectively (the area of the substrate A s was included in A w due to its insignificance); 2.
The sputtering process occurs under isothermal conditions. This means that the target and the wall have the same temperature, T, which is equal to the temperature of the gas environment; 3.
Two processes compete on the target surface: • Formation of the M m X n compound due to chemisorption with the kinetics (dθ t /dt) chem (see Figure 2). Quantitatively, it is characterized by the sticking coefficient α, which specifies the fraction of the adsorbed gas flow with density J X2 incident on the surface (in Formulas (3) and (5)-(7), this value is denoted by F): where p X2 and m X2 are the partial pressure and mass of a reactive gas molecule, respectively; • Sputtering of the M m X n compound by argon ions with kinetics (dθ t /dt) sp .

4.
Chemisorption on surfaces covered by the M m X n compound is negligible; 5.
Sputtering of the M m X n compound from the target surface occurs in the form of molecules that are deposited on the walls of the chamber; 6.
The process has two independent variables: the discharge current density of argon ions on the target j and the Q 0 reactive gas incoming flow (the flow introduced into the vacuum chamber).

Initial Model
The model that the authors described in [46] was later called the Berg model. It should be noted that the group of Professor S. Berg has been working hard on the development of the reactive sputtering model over the past thirty-five years and has published more than a dozen articles on this topic.
Further, when discussing different models, we use notations that partially differ from those adopted by the authors of the referred works and used in Section 1 of this article. This was done in order to unify the mathematical descriptions of different models.
In the first work of 1987 [46], the ideas expressed in the mentioned works were developed for the reactive sputtering of a one-component metal target M in the environment of one reactive gas X2. Describing the physical model, the authors proposed a number of constraints: 1. The model should take into account the processes occurring on the surfaces of the target and the wall of the vacuum chamber with areas At and Aw, respectively (the area of the substrate As was included in Aw due to its insignificance); 2. The sputtering process occurs under isothermal conditions. This means that the target and the wall have the same temperature, T, which is equal to the temperature of the gas environment; 3. Two processes compete on the target surface: • Formation of the MmXn compound due to chemisorption with the kinetics (dθt/dt)chem (see Figure 2). Quantitatively, it is characterized by the sticking coefficient α, which specifies the fraction of the adsorbed gas flow with density JX 2 incident on the surface (in Formulas (3) and (5)-(7), this value is denoted by F): where pX 2 and mX 2 are the partial pressure and mass of a reactive gas molecule, respectively; JtM is the flux density of sputtered metal atoms; JtC is the flux density of sputtered MmXn molecules; j is the discharge current density; θt is the fraction of the target surface covered by the compound film.

•
Sputtering of the MmXn compound by argon ions with kinetics (dθt/dt)sp.
4. Chemisorption on surfaces covered by the MmXn compound is negligible; 5. Sputtering of the MmXn compound from the target surface occurs in the form of molecules that are deposited on the walls of the chamber; 6. The process has two independent variables: the discharge current density of argon ions on the target j and the Q0 reactive gas incoming flow (the flow introduced into the vacuum chamber).

Figure 2.
The flux densities on the target surface. Notations: J X2 is the flux density of X 2 molecules; J tM is the flux density of sputtered metal atoms; J tC is the flux density of sputtered M m X n molecules; j is the discharge current density; θ t is the fraction of the target surface covered by the compound film.
The physical model presented above in items 1-6 was described by a series of algebraic equations. Two of them set the kinetics of processes on selected surfaces. For the target surface, the kinetic equation similar to Equation (1) is set: where For the steady state, Equation (9), taking into account (10), has the form: During target sputtering, two fluxes arise (see Figure 2): M atoms with density J tM and M m X n molecules with density J tC , which are deposited on the wall surface. Sputtering is absent there, and the kinetic equation for it was written as: where the term (dθ w /dt) chem describes the increase in the fraction of the wall surface θ w covered by the compound film, resulting from chemisorption (see Figure 3); (dθ w /dt) spC describes the increase in θ w due to the flux of M m X n molecules sputtered from the target surface and falling onto a fraction of the wall covered by metal; (dθ w /dt) spM describes the decrease in θ w due to the flux of M atoms sputtered from the target surface and falling onto a fraction of the wall covered by the M m X n compound.
where the term (dθw/dt)chem describes the increase in the fraction of the wall surface θw covered by the compound film, resulting from chemisorption (see Figure 3); (dθw/dt)spC describes the increase in θw due to the flux of MmXn molecules sputtered from the target surface and falling onto a fraction of the wall covered by metal; (dθw/dt)spM describes the decrease in θw due to the flux of M atoms sputtered from the target surface and falling onto a fraction of the wall covered by the MmXn compound. For the stationary state of the sputtering process, Equation (12), written taking into account all the details, takes the form In addition to (11) and (13), three equations define the reactive gas flows in the sputtering system. Two of them are similar to (6) at F = JX 2 with subscripts for the target surfaces i = t and for the wall i = w. In addition, a part of the gas Qp is pumped out by a vacuum pump at a speed Sp: For the stationary state of the sputtering process, Equation (12), written taking into account all the details, takes the form In addition to (11) and (13), three equations define the reactive gas flows in the sputtering system. Two of them are similar to (6) at F = J X2 with subscripts for the target surfaces i = t and for the wall i = w. In addition, a part of the gas Q p is pumped out by a vacuum pump at a speed S p : The sixth equation of the analytical description of the reactive sputtering model, which is illustrated in Figure 4, sets the balance of gas flows: The numerical solution of the system of Equation (6) for i = t and w, (11), (13)-(15) allowed the authors of [46] to determine the dependences p = f (Q 0 ) at j = const. The authors qualitatively demonstrated the adequacy of the model for sputtering a titanium target in a nitrogen-containing environment.
The sixth equation of the analytical description of the reactive sputtering model, which is illustrated in Figure 4, sets the balance of gas flows: The numerical solution of the system of Equation (6) for i = t and w, (11), (13)-(15) allowed the authors of [46] to determine the dependences p = f(Q0) at j = const. The authors qualitatively demonstrated the adequacy of the model for sputtering a titanium target in a nitrogen-containing environment.

Development of the Model
Developing the initial model, the authors in [47] proposed to take into account the reactive gas consumption on the full surfaces of the target and the chamber wall. This meant that gas adsorption was also allowed on parts of the surfaces covered by the compound. In the model with this assumption, the authors did not use the wall, replacing it with a substrate with area As. In essence, under isothermal conditions, the name of the surface on which the particles of the substance sputtered from the target surface are deposited does not matter. It is important that its area has a significant effect on the sputtering process. Therefore, if in [46] the authors assumed Aw >> As, then in [47], their obvious assumption was Aw << As. In this regard, further, when writing the corresponding equations, the subscript symbol w will be replaced by the symbol s. Taking into account the introduced assumption, the equations for reactive gas flows on the surface of the target and substrate of type (6) take the form where αtM, αtC, αsM, and αsC-adsorption coefficients of reactive gas atoms to metal M and compound MmXn on the target and substrate surfaces, respectively. Obviously, at αtC = αsC = 0.01, the reactive gas consumption will increase by no more than 2%. However, the equations describing the steady states of both surfaces of the type (11) and (13) cannot change since the values of θt and θs are influenced by the sorption of X2 gas only on the metal (fractions of the target surface 1-θt and 1-θs). In other words, in the equations of the steady state of both surfaces given in [47], the terms written in the accepted notation

Development of the Model
Developing the initial model, the authors in [47] proposed to take into account the reactive gas consumption on the full surfaces of the target and the chamber wall. This meant that gas adsorption was also allowed on parts of the surfaces covered by the compound. In the model with this assumption, the authors did not use the wall, replacing it with a substrate with area A s . In essence, under isothermal conditions, the name of the surface on which the particles of the substance sputtered from the target surface are deposited does not matter. It is important that its area has a significant effect on the sputtering process. Therefore, if in [46] the authors assumed A w >> A s , then in [47], their obvious assumption was A w << A s . In this regard, further, when writing the corresponding equations, the subscript symbol w will be replaced by the symbol s. Taking into account the introduced assumption, the equations for reactive gas flows on the surface of the target and substrate of type (6) take the form where α tM , α tC , α sM, and α sC -adsorption coefficients of reactive gas atoms to metal M and compound M m X n on the target and substrate surfaces, respectively. Obviously, at α tC = α sC = 0.01, the reactive gas consumption will increase by no more than 2%. However, the equations describing the steady states of both surfaces of the type (11) and (13) cannot change since the values of θ t and θ s are influenced by the sorption of X 2 gas only on the metal (fractions of the target surface 1 − θ t and 1 − θ s ). In other words, in the equations of the steady state of both surfaces given in [47], the terms written in the accepted notation 2α tC J X 2 (1 − θ t ) and 2α sC J X 2 (1 − θ s ), can be discarded. In [47,48], the authors made a very useful suggestion. They extended the analytical description of the model with an expression for the target sputtering rate: The family of curves R = f (Q 0 ) j = const is qualitatively shown in Figure 5. The authors assumed that the points corresponding to the deposition modes of a constant composition film lie in the R − Q 0 plane on a straight line. Therefore, if we add straight lines coming from the origin of coordinates to the family R = f (Q 0 ) j = const , then each of them will represent the deposition process of the constant composition film. Figure 5 shows one of these straight lines corresponding to a film of stoichiometric composition (for θ t = 1). A very important conclusion follows from this hypothesis: the growth rate of a film of a given composition can be changed by varying the discharge current, but it is also necessary to change the value of Q 0 in such a way that the operating point remains on the corresponding straight line in the R − Q 0 plane. the deposition process of the constant composition film. Figure 5 shows one of these straight lines corresponding to a film of stoichiometric composition (for θt = 1). A very important conclusion follows from this hypothesis: the growth rate of a film of a given composition can be changed by varying the discharge current, but it is also necessary to change the value of Q0 in such a way that the operating point remains on the corresponding straight line in the R − Q0 plane. The work [49] was the next step in the development of the Berg model. In this work, the authors proposed models of two reactive sputtering processes. In the first one, the target of the two-element metal alloy yM1 + (1 -y)M2 was sputtered in the mixture of Ar + X2. Figures 6 and 7 illustrate the model of this process proposed by the authors. Processes of forming and sputtering M1X and M2X films compete on the target containing M1 and M2 metals ( Figure 6). Fractions of the targe surface covered by the above films are θtC 1 and θtC 2 , respectively. In this problem, the equation of the target steady state can be written in the form: where JX 2 is the reactive gas flow density, which is defined by an expression like (8); αM 1 and αM 2 are the adhesion coefficients of X2 molecules to fractions of the target surface where there are no M1X and M2X films, respectively; θtM 1 and θtM 1 are the fractions of the The work [49] was the next step in the development of the Berg model. In this work, the authors proposed models of two reactive sputtering processes. In the first one, the target of the two-element metal alloy yM 1 + (1 − y)M 2 was sputtered in the mixture of Ar + X 2 . Figures 6 and 7 illustrate the model of this process proposed by the authors. Processes of forming and sputtering M 1 X and M 2 X films compete on the target containing M 1 and M 2 metals ( Figure 6). Fractions of the targe surface covered by the above films are θ tC 1 and θ tC 2 , respectively. In this problem, the equation of the target steady state can be written in the form: where J X2 is the reactive gas flow density, which is defined by an expression like (8); α M 1 and α M 2 are the adhesion coefficients of X 2 molecules to fractions of the target surface where there are no M 1 X and M 2 X films, respectively; θ tM 1 and θ tM 1 are the fractions of the target surface where there are no M 1 X and M 2 X films, respectively; S C 1 and S C 2 are the sputtering yields of M 1 X and M 2 X films, respectively. The flux sputtered from the target surface consists of fluxes of atoms M 1 , M 2 and molecules M 1 X, M 2 X with the corresponding densities, which are denoted by J tM 1 , J tM 2 , J tC 1 , and J tC 2 in Figure 6. The same fluxes fall onto the substrate; their densities are denoted by J sM 1 , J sM 2 , J sC 1 , and J sC2 in Figure 7.
Other equations of the analytical description of the model are compiled by analogy with Equations (13)- (17).     In the second process, the co-sputtering of two targets made from different metals was studied. The physical model of this process is illustrated in Figure 8. The targets of metals M1, M2 are sputtered in the mixture of Ar + X2. As in all previous cases, the processes of formation and sputtering of compound films compete on their surfaces. The processes occurring on the targets are not related; therefore, in the steady state, two equations arise with two independent variables: J sM 2 is the flux density of metal atoms M 2 ; J sC 1 is the flux density of M 1 X molecules; J sC2 is the flux density of M 2 X molecules; θ sM 1 is the fraction of the substrate surface covered by the metal film M 1 ; θ sM 2 is the fraction of the substrate surface covered by the metal film M 2 ; θ sC 1 is the fraction of the substrate surface covered by the compound film M 1 X; θ sC2 is the fraction of the substrate surface covered by the compound film M 2 X.
In the second process, the co-sputtering of two targets made from different metals was studied. The physical model of this process is illustrated in Figure 8. The targets of metals M 1 , M 2 are sputtered in the mixture of Ar + X 2 . As in all previous cases, the processes of formation and sputtering of compound films compete on their surfaces. The processes occurring on the targets are not related; therefore, in the steady state, two equations arise with two independent variables: Materials 2023, 16, x FOR PEER REVIEW 13 of 56 ( ) The rules for compiling other equations of a system describing reactive sputtering of this type should not cause difficulties. The derivation of these equations is similar to the procedures used by Professor Berg and co-workers in the publications we reviewed earlier.
In this work, the authors admitted that they could not reliably demonstrate the adequacy of the proposed model for the case in which Ti, Al, and the Ar + O2 gas mixture were used. In addition, they admitted that neither theoretically nor experimentally it was possible to detect any individual influence of the used metals on the dependencies p = f(Q0). Let us note that the assumption of the authors that the composition of the metal alloy film deposited on the substrate will always be identical to the bulk composition of the target is not always valid. In some cases, the effect of preferential sputtering occurs, in which the target surface is continuously depleted of a component with a higher sputtering yield [63].
The solution to an interesting problem in the development of the chemisorption model is given in [52]. In contrast to [46], the authors complicated the process of reactive sputtering by placing one metal target in an Ar + O2 + N2 gas mixture. In the model of this process, as in [47], there is also no wall, and the following main assumptions were made: The rules for compiling other equations of a system describing reactive sputtering of this type should not cause difficulties. The derivation of these equations is similar to the procedures used by Professor Berg and co-workers in the publications we reviewed earlier.
In this work, the authors admitted that they could not reliably demonstrate the adequacy of the proposed model for the case in which Ti, Al, and the Ar + O 2 gas mixture were used. In addition, they admitted that neither theoretically nor experimentally it was possible to detect any individual influence of the used metals on the dependencies p = f (Q 0 ). Let us note that the assumption of the authors that the composition of the metal alloy film deposited on the substrate will always be identical to the bulk composition of the target is not always valid. In some cases, the effect of preferential sputtering occurs, in which the target surface is continuously depleted of a component with a higher sputtering yield [63].
The solution to an interesting problem in the development of the chemisorption model is given in [52]. In contrast to [46], the authors complicated the process of reactive sputtering by placing one metal target in an Ar + O 2 + N 2 gas mixture. In the model of this process, as in [47], there is also no wall, and the following main assumptions were made:

1.
Two processes compete on the target surface: • Formation of the MO + MN compound due to chemisorption. The process was quantitatively characterized by the sticking coefficients of gases to the metal α O and α N , the molecules of which fall on the target surface with flux densities J O 2 and J N 2 , each of which is given by Formula (8); • Sputtering of the MO + MN compound with argon ions.

2.
Sputtering from the target surface of the MO + MN compound occurs in the form of MO and MN molecules, which are uniformly deposited on the substrate (the chamber wall is excluded from the model); 3.
The process has three independent variables: the discharge current density of the argon ions on the target j and incoming flows of reactive gases Q O 2 and Q N2 .
To facilitate understanding of the analytical description of the model, we present it in the notation already used above. They differ from the author's notations in the referred article. The steady state of the target surface in the new problem is now described by two equations. The steady state of the fraction of the target surface θ tC covered by oxide MO gives the equation: where θ tC 2 is the fraction of the target surface covered by nitride; α 12 is the coefficient describing the probability of replacing the nitrogen with the oxygen in MN, as the metaloxygen bond is energetically more favorable than the metal-nitrogen bond; S C 1 is the oxide sputtering yield. For a fraction of the target surface θ tC 2 covered by MN nitride, the following holds true in the steady state: where S C 2 is the sputtering yield of nitride. The steady state of the substrate surface in this problem is described by two equations similar to (13), where the subscript w is replaced by s, separately for fractions covered by oxide θ sC 1 and nitride θ sC2 .
The flows of reactive gases in the sputtering system set: • Oxygen flows on the surfaces of the target i = t and wall i = w: • Nitrogen flows on the surfaces of the target i = t and wall i = w: • Gas flows, which are pumped out: where S p is a pumping speed of the vacuum chamber (does not depend on the type of gas); • Balance equation for flows: Thus, the analytical description of the physical model of the reactive sputtering of a single metal target in an Ar + O 2 + N 2 mixture contains 10 equations. Calculations performed by the authors using expressions (19)-(26) made it possible to establish that the addition of nitrogen to the gas environment contributes to a decrease in the hysteresis effect up to its complete disappearance. In addition, the rate of film deposition increases in comparison with the sputtering of oxides without the addition of nitrogen.
An experimental confirmation of the effects predicted by the model in [52] was performed in [53]. The authors showed that the calculated dependences of the partial pressures of both reactive gases on the incoming nitrogen flow are close to those found experimentally.
In [60], an even more complicated problem of reactive sputtering of two targets (aluminum and zirconium) in an Ar + O 2 + N 2 mixture was studied. Detailed coverage of the evolution of the initial model allows for refraining from presenting the details of [60] here. It is obvious that a second target has appeared in the physical model of this process. Two flows of reactive gases fall on each target in the same way as for a single target in the previous problem. Therefore, in the analytical description of the model, equations of the steady state of each target arose. Fluxes of particles sputtered from their surfaces are deposited on the substrate, which leads to a corresponding change in the equation of its steady state. According to the authors, the model is in qualitative agreement with the experimental results.
In [64], an improved model for reactive sputtering in a vacuum chamber with an Ar + X 2 mixture is discussed. The magnetron, equipped with an aluminum target, is fixed at the top. The Berg simple isothermal chemisorption model was refined by expanding the list of assumptions: 1.
The target is sputtered only in an annular area A r , called the racetrack. Its shape and position on the target are determined by the magnetic field; 2.
Particles sputtered from the target surface are deposited on: • A substrate with area A s located in the center part of the bottom of the chamber; • The peripheral part of the bottom with area A pb ; • The chamber walls with area A w ; • The top of the chamber with area A top , including the non-sputtered part of the target.

3.
There is no redeposition of sputtered particles in the upper part of the chamber.
To describe the steady state of the target in this problem, Equation (11) can be used, replacing the variable θ t with θ r . For a similar state of the i-th (i = s, pb, w, and top) surface on which the sputtered particles are deposited, the following equation should be used: which is constructed by analogy with Equation (13). The main part of the notation in (27) was used by us earlier. In addition, in (27), it is assumed: θ i is the fraction of the area of the i-th (i = s, pb, w, and top) surface covered by the M m X n compound; θ r is the fraction of the area of the sputtered annular target area covered by the M m X n compound; f i is the fraction of the material sputtered from the target that is deposited on the i-th (i = s, pb, w, and top) surface.
Reactive gas flows to selected surfaces are given by equations of the form (6): when the gas balance is equal to In all previous modeling works, the authors used the incoming flow of the reactive gas and the discharge current as independent variables. In [64], it was proposed to supplement where f 0 -ratio of the sputtered flux in the target plane to the scattered flux at a distance d 0 from the target. In conclusion of the discussion of work [64], in which Equations (27)- (30) are the main ones, we note that it presents the next step in the development of the reactive sputtering model. For the first time, the surface of the substrate appeared in it as an independent object. If only this innovation is accepted, excluding other surfaces, then two more equations should be added to the system of equations of the initial model containing six equations. The first of them should describe the steady state of the substrate surface. It can be obtained from Equation (13) by replacing the subscript w with the index s. The second equation for the gas flow consumed by the substrate can be obtained from (6) at i = s.

Adequacy of the Initial Model
Until now, the Berg initial isothermal chemisorption model has been a popular tool for studying the processes of film deposition of simple compounds with the reactive sputtering method [65][66][67]. Undoubtedly, it was an advance in the modeling of the reactive sputtering process. Its analytical description allows us to evaluate the relationship between the process parameters and shows that a hysteresis can exist.
However, as noted in [65], this model has very rough assumptions. In agreement with this statement, let us pay attention to only two circumstances. First, the deposited film, in terms of its physical properties, is a chemical compound that is formed on the substrate mainly due to the plasma-chemical reaction occurring on the target. The replacement of the reaction by the mechanism of chemisorption should be recognized as a very strong simplification. The second strong simplification is the assumption that the temperatures of all surfaces are equal. In real processes, the surfaces of the target, substrate, and chamber walls have different temperatures: the target surface can be heated up to 700-900 K, the substrate is heated above 600 K, and the wall temperature can be maintained at 300 K. In addition, the substrate temperature can steadily increase during deposition due to the release of condensation energy, the kinetic energy of particles, discharge radiation, and also the energy released by chemical reactions [16,[68][69][70].
Let us show by examples the possibilities of the initial Berg model. For this, we use the experimental results presented in [25,32]. In [32], the authors studied the deposition of Ta 2 O 5 films at sputtering system parameters equal to S p = 8.6 L/s; A t = 100 cm 2 ; A w = 300 cm 2 ; j = 0.05 A/cm 2 . The result of solving the system of Equation (6) at i = t and w, (12), (14)- (16) for this experiment at α = 1, S C = 0.024, S M = 0.6 is shown by the solid line in Figure 9a.
In [25], the authors studied the deposition of TiO 2 films at the sputtering system parameters of S = 190 L/s; A t = 190 cm 2 ; A w = 9000 cm 2 ; j = 0.011 A/cm 2 . The result of solving the system of Equation (6) at i = t and w, (12), (14)- (16) for this experiment at α = 1, S c = 0.016, S m = 0.32 is shown by the solid line in Figure 9b. It can be seen from both figures that the Berg model does not adequately predict the points at which the operating modes of the target change.
The obtained results may indicate both the inadequacy of the initial Berg model and the inaccuracy of the values of the model parameters that were used in the calculation. Let us carry out a more detailed analysis of the example of an experiment on the deposition of the Ta 2 O 5 film in [32]. Figure 9a shows that the analytical dependence gives overestimated values of the flow Q 0 at both points of instability A and C. In addition, the linear section of the curve in the region of the reactive mode is significantly shifted relative to the experimental points and has a slight difference in the slope angle. The first difference indicates an overestimated gas flow to the target surface, which maintains the steady oxide mode of its operation. It significantly depends on the target area A t . The angle of inclination of this section reflects the operation of the vacuum pump only. The tangent of this angle is equal to 1/S p .
simplification. The second strong simplification is the assumption that the temperatures of all surfaces are equal. In real processes, the surfaces of the target, substrate, and chamber walls have different temperatures: the target surface can be heated up to 700-900 K, the substrate is heated above 600 K, and the wall temperature can be maintained at 300 K. In addition, the substrate temperature can steadily increase during deposition due to the release of condensation energy, the kinetic energy of particles, discharge radiation, and also the energy released by chemical reactions [16,[68][69][70].
Let us show by examples the possibilities of the initial Berg model. For this, we use the experimental results presented in [25,32]. In [32], the authors studied the deposition of Ta2O5 films at sputtering system parameters equal to Sp = 8.6 L/s; At = 100 cm 2 ; Aw = 300 cm 2 ; j = 0.05 A/cm 2 . The result of solving the system of Equation (6) at i = t and w, (12), (14)- (16) for this experiment at α = 1, SC = 0.024, SM = 0.6 is shown by the solid line in Figure  9a. In [25], the authors studied the deposition of TiO2 films at the sputtering system parameters of S = 190 L/s; At = 190 cm 2 ; Aw = 9000 cm 2 ; j = 0.011 A/cm 2 . The result of solving the system of Equation (6) at i = t and w, (12), (14)- (16) for this experiment at α = 1, Sc = 0.016, Sm = 0.32 is shown by the solid line in Figure 9b. It can be seen from both figures that the Berg model does not adequately predict the points at which the operating modes of the target change. These differences between the experiment and the simulation result were eliminated by selecting adequate values of A t and S p . Figure 10 shows the dependencies built on the basis of the changed data, the values of which are indicated in the text below the figure. Figure 10 demonstrates that the target area indicated in [32] is overestimated by about 2.9 times, and the pumping speed of the vacuum chamber is overestimated by 1.08 times. Here, it can be assumed that the authors of [32] took into account the total target area, which was used in calculating the current density under the assumption of its uniform distribution along the target surface. While the effective area involved in the process of reactive sputtering, which sets the configuration of the magnetic field, is significantly smaller. The obtained results may indicate both the inadequacy of the initial Berg model and the inaccuracy of the values of the model parameters that were used in the calculation. Let us carry out a more detailed analysis of the example of an experiment on the deposition of the Ta2O5 film in [32]. Figure 9a shows that the analytical dependence gives overestimated values of the flow Q0 at both points of instability A and C. In addition, the linear section of the curve in the region of the reactive mode is significantly shifted relative to the experimental points and has a slight difference in the slope angle. The first difference indicates an overestimated gas flow to the target surface, which maintains the steady oxide mode of its operation. It significantly depends on the target area At. The angle of inclination of this section reflects the operation of the vacuum pump only. The tangent of this angle is equal to 1/Sp.
These differences between the experiment and the simulation result were eliminated by selecting adequate values of At and Sp. Figure 10 shows the dependencies built on the basis of the changed data, the values of which are indicated in the text below the figure. Figure 10 demonstrates that the target area indicated in [32] is overestimated by about 2.9 times, and the pumping speed of the vacuum chamber is overestimated by 1.08 times. Here, it can be assumed that the authors of [32] took into account the total target area, which was used in calculating the current density under the assumption of its uniform distribution along the target surface. While the effective area involved in the process of reactive sputtering, which sets the configuration of the magnetic field, is significantly smaller. We retain all the approximations adopted in the model under consideration, except for the second one, in which the temperature of the target and the walls of the vacuum chamber are assumed to be equal. We denote them using Tt and T, respectively, assuming that the temperatures of the walls and the gas are the same. The substrate in this model is not separated into an independent object; its area is included in the wall area. Therefore, the substrate temperature is also considered to be equal to T, although in real processes, this parameter is set at a level of 600 K and higher, which significantly affects the film deposition process.
The separation of temperatures Tt and T is, in fact, a partial transition to a non-iso- We retain all the approximations adopted in the model under consideration, except for the second one, in which the temperature of the target and the walls of the vacuum chamber are assumed to be equal. We denote them using T t and T, respectively, assuming that the temperatures of the walls and the gas are the same. The substrate in this model is not separated into an independent object; its area is included in the wall area. Therefore, the substrate temperature is also considered to be equal to T, although in real processes, this parameter is set at a level of 600 K and higher, which significantly affects the film deposition process.
The separation of temperatures T t and T is, in fact, a partial transition to a nonisothermal model. Characterizing chemisorption on different surfaces, we introduce the sticking coefficients for the target α t and for the wall α w , assuming α t = α w . Let us assume that as the temperature increases, the sticking coefficient of reactive gases decreases. Therefore, when choosing the values of α t and α w , two conditions must be met: α t < 1, α w < 1 and α t < α w , as T t > T.
Let us rewrite Equation (6) at i = t, w, (11) and (13), taking into account accepted assumptions and notations: Based on the solutions of these equations, we can assume that the elimination of the isothermal condition from the initial Berg model makes it possible to increase its adequacy. However, another improvement of the reactive sputtering model is also acceptable. As Figure 10 shows, it is important to change the physical model in order to shift the transition point of the target from the metallic to the reactive mode to the right. To do this, it is necessary to increase the growth rate of the compound film by adding an additional mechanism of the film formation to the model.

New Mechanism
For the first time Professor Berg and co-workers in [71,72] with reference to [73], proposed to introduce an assumption about a significant influence of the implantation of reactive gas ions into the near-surface layer of the target on reactive sputtering. In [74], equations were given that describe a new model of reactive sputtering of a metal target in an Ar + X 2 mixture. In [75], the updated Berg model is described in detail. The authors supplemented the initial model with new assumptions:

•
The M m X n compound is sputtered from the target in the form of atoms; • The formation of the M m X n compound on the target occurs not only on the surface due to the chemisorption of X 2 molecules but also in the volume under the surface due to the interaction of free metal atoms with X atoms that penetrated into it; • X atoms penetrated into the target volume due to the implantation of accelerated X + 2 ions and adsorbed X atoms, which received momentum due to the impact with the Ar + ion. The latter mechanism was called the "knock-on effect"; • Sputtering of the target by X + 2 ions is negligible. The target in the updated model only serves as a source of the metal flux with density J tM (Figure 11). At the same time, the kinetic Equation (9) acquired new components: where (dθ t /dt) impl is the growth rate of the M m X n compound film due to ion implantation; (dθ t /dt) kn is the rate of removal of reactive gas atoms due to the knock-on effect.  For the steady state, Equation (34) takes the form where p is the partial pressure of the reactive gas; pA is the partial pressure of Ar; Sk is a coefficient similar to the sputtering yield, which determines the number of X atoms implanted into the target volume due to impact with one Ar + ion; αi is the probability of implantation of a reactive gas ion into the target volume. The rest of the designations in (34) were used earlier. The arrangement of terms in (35) corresponds to formula (34). This makes it possible to understand the physical meaning of each of them. For example, the first term describes the formation of a film on the target surface due to chemisorption, the second term-due to implantation, and so on. In this problem, only fluxes of reactive gas molecules and metal atoms with densities specified below fall on the substrate. In this case, the respective densities are JX 2 (see (8)) and ( ) where SMC is a coefficient determining the partial sputtering of metal atoms from the MmXn layer. Taking into account (36), the stationary equation for the substrate becomes simpler than (13): As in all the models considered here, the system of equations describing the physical model should include equations for reactive gas flows to the target and substrate, which can be used as (31), the pumping Equation (14), and the gas balance equation of the type (15).
In conclusion, the authors noted that the inclusion of atomic MmXn sputtering and implantation in the initial Berg model was intended to expand the number of modeling parameters. However, from a practical point of view, these interesting proposals complicate the modeling since, in each case, they require a search for the values of these parameters.
In addition, in [75], the authors described a two-layer model based on a more complex dynamical model, which they described earlier in [74]. In this publication, the authors develop the Berg model by assuming that the compound layer on the target surface is several tens of angstroms thick. In the new model, this continuous layer is represented as N discrete layers. For the steady state, Equation (34) takes the form where p is the partial pressure of the reactive gas; p A is the partial pressure of Ar; S k is a coefficient similar to the sputtering yield, which determines the number of X atoms implanted into the target volume due to impact with one Ar + ion; α i is the probability of implantation of a reactive gas ion into the target volume. The rest of the designations in (34) were used earlier. The arrangement of terms in (35) corresponds to formula (34). This makes it possible to understand the physical meaning of each of them. For example, the first term describes the formation of a film on the target surface due to chemisorption, the second term-due to implantation, and so on.
In this problem, only fluxes of reactive gas molecules and metal atoms with densities specified below fall on the substrate. In this case, the respective densities are J X2 (see (8)) and where S MC is a coefficient determining the partial sputtering of metal atoms from the M m X n layer. Taking into account (36), the stationary equation for the substrate becomes simpler than (13): As in all the models considered here, the system of equations describing the physical model should include equations for reactive gas flows to the target and substrate, which can be used as (31), the pumping Equation (14), and the gas balance equation of the type (15).
In conclusion, the authors noted that the inclusion of atomic M m X n sputtering and implantation in the initial Berg model was intended to expand the number of modeling parameters. However, from a practical point of view, these interesting proposals complicate the modeling since, in each case, they require a search for the values of these parameters.
In addition, in [75], the authors described a two-layer model based on a more complex dynamical model, which they described earlier in [74]. In this publication, the authors develop the Berg model by assuming that the compound layer on the target surface is several tens of angstroms thick. In the new model, this continuous layer is represented as N discrete layers.
In [76], the authors applied a new model to study the hysteresis. Studies of the influence on the dependence p = f (Q 0 ) of various model parameters allowed them to establish the possibilities for reducing or completely eliminating hysteresis.
Let us pay attention to one more work by Berg et al. [54]. When sputtering a titanium target in an Ar + CH 4 mixture, the authors found that depending on Q 0 , a film with a carbon content from 0 to 100% can be formed on the substrate. The authors expected similar behavior when carrying out the processes in Ar + B 2 H 6 and Ar + Si n H 2n+2 . For such processes, which the authors called "non-saturated", the initial Berg model turned out to be unsuitable.
The new model proposed by the authors will be explained with the help of Figures 12 and 13. On the target (Figure 12), in contrast to the initial model, the processes of formation and sputtering of TiC and C films compete. In accordance with this, the steady state equation of the target takes the form (38) where J CH4 is the methane flux density, which is defined by an expression like (8); α Ti , α TiC , and α C are the adhesion coefficients of methane molecules to parts of the target surface with open titanium and covered by titanium carbide TiC and carbon C, respectively; θ Ti , θ TiC , and θ C are the fractions of the target surface with open titanium and covered by carbide TiC and carbon C, respectively; S TiC and S C are the sputtering yields of titanium and carbon carbide films, respectively. The flux sputtered from the target surface consists of fluxes of Ti and C atoms and TiC molecules with the corresponding densities, which are denoted by J tTi , J tC , and J tTiC in Figure 12. The same fluxes fall onto the substrate with densities, which are denoted by J sTi , J sC , and J sTiC in Figure 13. The equation of the steady state of the substrate in this problem, as well as other equations of the system, are compiled by analogy with Equations (13)- (17). The model developed in this way, according to the authors, is adequate for the experimental results.
In [76], the authors applied a new model to study the hysteresis. Studies of the influence on the dependence p = f(Q0) of various model parameters allowed them to establish the possibilities for reducing or completely eliminating hysteresis.
Let us pay attention to one more work by Berg et al. [54]. When sputtering a titanium target in an Ar + CH4 mixture, the authors found that depending on Q0, a film with a carbon content from 0 to 100% can be formed on the substrate. The authors expected similar behavior when carrying out the processes in Ar + B2H6 and Ar + SinH2n+2. For such processes, which the authors called "non-saturated", the initial Berg model turned out to be unsuitable.
The new model proposed by the authors will be explained with the help of Figures  12 and 13. On the target (Figure 12), in contrast to the initial model, the processes of formation and sputtering of TiC and C films compete. In accordance with this, the steady state equation of the target takes the form Ti CH tTiC tC TiC CH tTiC C CH tC TiC tTiC C tC where JCH 4 is the methane flux density, which is defined by an expression like (8)   In a recent paper [77], Berg et al. proposed a model of reactive high-power pulsed sputtering, demonstrating it using the example of Ti sputtering in an Ar + O2 mixture. The system of equations describing the proposed model is similar to the system for the initial model and contains equations of the type (6), (11), and (13)- (15).
Ease of learning and use provided the Berg model with multilateral interest in the last twenty years. The isothermal chemisorption model was used by experts to describe practical problems [78][79][80][81][82][83][84][85]. Some of them suggested options for its development [86][87][88][89][90]. In a recent paper [77], Berg et al. proposed a model of reactive high-power pulsed sputtering, demonstrating it using the example of Ti sputtering in an Ar + O 2 mixture. The system of equations describing the proposed model is similar to the system for the initial model and contains equations of the type (6), (11), and (13)- (15).

Isothermal Model of Depla (the RSD Model)
Let us draw the reader's attention to the work of Professor Depla and his group. This team has been actively involved in studies of reactive sputtering and its modeling for the last twenty years. In terms of the intensity of publications in this area, the group from Ghent University can confidently be given second place after the group of Professor Berg. The results of experimental studies of reactive sputtering and the properties of the obtained films of oxides, nitrides, etc., were published in articles [73,[97][98][99][100][101][102][103][104][105][106][107][108][109][110][111][112][113][114][115]. Various aspects of reactive sputtering modeling are described by the researchers of the team in publications [65,. Let us take a look at the main ones.

Initial RSD Model
The work [116] provides data on the influence of the incoming oxygen flow on the change in voltage on an aluminum target sputtered in an Ar + O 2 mixture. In this article, the authors, for the first time, put forward the hypothesis that chemisorption alone cannot explain the modification of the target surface. The observed changes in the discharge voltage, in their opinion, can be additionally initiated by the chemical reaction of target oxidation in the near-surface layers by implanted oxygen atoms. Implantation and subsurface bulk response can be seen as the most important extension of the Berg original model. Its inclusion in the reactive sputtering model led to the appearance of the original model, which in [127] was called the RSD (Reactive Sputtering Deposition) model.
The same hypothesis was put forward in [73], devoted to the study of reactive sputtering of a copper target in Ar + N 2 . A similar conclusion was made in [97] based on the results of a study of silicon target sputtering in the same mixture. Detailing the process, the authors in [97] suggested that molecular nitrogen ions N + 2 are implanted into the target, which, losing their charge, dissociate into individual atoms that react with the target material to form the M m X n compound. Based on these studies, the physical model of reactive sputtering was supplemented by implantation.
In order to form the rationale for their hypothesis about implantation, the authors referred to the results of a study from [140]. In this article, using the method of nuclear reactions, the reactive sputtering of a titanium target in a nitrogen-containing environment was studied.
It was found that the maximum number of retained nitrogen atoms after turning off the magnetron significantly exceeds the adsorbed monolayer. The result obtained, according to the authors, is associated with the implantation of nitrogen molecular ions with subsequent decay into atoms. In addition, the authors assumed that there was an implantation of atoms from adsorbed nitrogen molecules, which received a direct impact from the accelerated argon ion. Such particles that have arisen during ion bombardment of a target are called recoil atoms. Having received an impulse, they can leave the target or create a cascade of collisions in its surface layers. Subsequently, the mechanism of implantation of recoil atoms was called the "knock-on effect" in the model of Depla. The results of the experiment were confirmed by computer simulation.
In [117,118], the authors took the first steps in developing a new model. In [117], the influence of only implanted atoms, which enter into chemical interaction with target atoms, on the process was shown.
In essence, the RSD model can be expressed by the kinetic Equation (34), although it is not written explicitly in Depla's works. The first attempt to write such an equation was made in [117] in the form where n 0 is the target material density; n is the number of gas atoms in the compound; α tC is the probability of the chemical reaction between implanted atoms X and target M; θ tC the degree of surface coverage by the M m X n compound; n(0, t) is the surface concentration of implanted atoms. Equation (39) describes a model of the sputtering process based only on the implantation of ions X + 2 into the target. In (39) there is no need to take into account sputtering of the compound by argon ions. Indeed, dθ tC /dt is the effective compound formation rate since the sputtering effect was taken into account in the calculation of dn/dt. In (39), the chemical reaction is given in the simplest form, which does not take into account the law of mass action and the Arrhenius equation. It is set using the coefficient α tC . In subsequent works, the authors complicated this part of the model.
In [118], the chemisorption mechanism from the Berg initial model was added to the model that takes into account only ion implantation. Thus, volumetric and surface mechanisms of M m X n formation appeared in the model, which can be expressed by the equation: If in Equation (40), as in (39), we assume that sputtering is already taken into account in the second term on the right side, then the third term should be considered as sputtering of a surface film formed due to chemisorption.
This new model with the addition of the knock-on effect is described in more detail in [119]. It takes into account the processes occurring on the surfaces of the target and substrate. One of the system equations describing sputtering establishes, as in all previous cases, the balance of gas flows of the type (15). Three more equations define the gas flows to both surfaces of the type (6) and the pumped-out flow (14). The equations of the steady state of the target and the substrate are the root equations.
Let us consider the work [119] in more detail since it is the key publication of Professor Depla and his team on modeling. The authors detailed the main assumptions of their model in this article. Subsequently, some variations of the RSD model were described in numerous articles. These include the monograph [122] and the tutorial [65], which is a review. In these works, one can get acquainted in detail with both the process of reactive magnetron sputtering and its modeling.
Initially, we highlight the main assumptions of the model in a way that can be understood from studying Depla's team publications. Using the assumptions described earlier in this work, the authors formed a physical model of the process, which was further described analytically using known physical or chemical laws, equations, and formulas. We believe that this is very important in order to see the essence of a particular model, as most often, its understanding is hindered for the reader by the dry lines of formulas.
Therefore, from the publication [119] and other articles of these authors, the RSD physical model can be expressed using a number of assumptions: (1) The process is isothermal, i.e., temperatures of all surfaces and gas are equal; (2) The process is accompanied by: (2.1) Direct implantation of X atoms into a target with a dose D ≈ 2f j + t/e = 2j + tp/e(p + p A ), determined by the ion current density j and the molar fraction of the reactive gas in the vacuum chamber f (e is the electron charge); (2.2) Implantation of X atoms due to the destruction of M m X n molecules on the target surface by a direct impact of an Ar + ion (knock-on effect with coefficient β representing the number of implanted atoms per ion in a flow with density j/e); (2.3) Chemisorption of X atoms resulting from the dissociation of adsorbed X 2 molecules on the target surface determined by the sticking coefficient α.
(3) The normalized distribution of projected ranges c(x) of ions X + 2 along the normal to the target surface is a Gaussian: where R p is the average ion range; ∆R p is the straggle; (4) The distribution of implanted atoms n(x, t) taking into account sputtering at a speed v 0 , is determined by the integral: where c(x − v 0 τ) is the normalized distribution taking into account sputtering; (5) The M m X n compound is formed: (5.1) On the surface due to chemisorption; (5.2) In the near-surface layers of the target due to the bulk chemical reaction mM + nX ↔ M m X n between a part of implanted X atoms with metal atoms. The reaction rate is determined by the expression where k is the rate constant of the bulk reaction; n X (x, t) is the concentration of unreacted implanted X atoms, taking into account sputtering (see expression (42)); n M (x, t) is the concentration of free target atoms.
(6) The formation of the M m X n layer on the surface is accompanied by the process: (6.1) Sputtering of M m X n molecules with argon ions; (6.2) Removal of X atoms from M m X n molecules due to the knock-on effect.
Let us clarify that the M m X n compound, which occupies a fraction of the target surface θ tC , was formed in the subsurface layer due to a bulk chemical reaction with the participation of implanted X atoms. It appeared on the surface as a result of sputtering of the surface and subsequent layers. In fact, due to this process, the surface layer moves inside the target, constantly increasing the value of θ tC . (7) In the area of the target surface, there are ( Figure 14): (7.1) The layer from which the sputtering takes place directly (Surface). It has a thickness on the order of a monomolecular layer s. In the general case, this layer consists of sections: θ tM , which is free of the M m X n compound film; θ tC and θ tCch , which are occupied by the M m X n compound formed due to direct implantation of X atoms and their chemisorption, respectively; (7.2) A subsurface region with a thickness L ≈ R p ± 3∆R p , into which reactive gas atoms are implanted, chemically interacting with the target material. In the model, this area is represented as a thin layer called Subsurface, as shown in Figure 14. It consists of the θ tCs region occupied by the M m X n compound formed due to the bulk chemical reaction, and 1 − θ tCs region occupied by free atoms of metal M.
of X atoms and their chemisorption, respectively; (7.2) A subsurface region with a thickness L ≈ Rp ± 3ΔRp, into which reactive gas atoms are implanted, chemically interacting with the target material. In the model, this area is represented as a thin layer called Subsurface, as shown in Figure 14. It consists of the θtC s region occupied by the MmXn compound formed due to the bulk chemical reaction, and 1 -θtC s region occupied by free atoms of metal M. It was noted in [119] that chemisorption and implantation are identical in their effect on the formation of a compound film on a target. The second process increases the thickness of this film, which appears not only on the surface but also in near-surface layers with a thickness of a few nanometers. The rate of film growth on the surface also increases. Figure 14. Model of target. Notations: J X2 is the flux density of X 2 molecules; J tM is the flux density of sputtered metal atoms; J tC is the flux density of sputtered M m X n molecules; θ tM is the fraction of the surface free from the film of the M m X n compound; θ tC is the fraction of the surface occupied by the M m X n compound formed due to the bulk chemical reaction with the participation of implanted X atoms; θ tCch is the fraction of the surface occupied by the M m X n compound formed due to chemisorption of X atoms; θ tCs is the fraction of the subsurface region in the s plane occupied by the M m X n compound formed due to the bulk chemical reaction mM + nX ↔ M m X n ; θ tMs = 1 − θ tCs is the fraction of the subsurface region occupied by free M metal atoms; n X (x), n M (x) is the volume concentration of free X and M atoms, respectively, which did not take part in the chemical reaction.
It was noted in [119] that chemisorption and implantation are identical in their effect on the formation of a compound film on a target. The second process increases the thickness of this film, which appears not only on the surface but also in near-surface layers with a thickness of a few nanometers. The rate of film growth on the surface also increases.
A new reactive sputtering model was presented in [119] in stages. In the first step, the authors took into account only the process of implantation of reactive gas ions, which is accompanied by a chemical reaction and sputtering. Under these conditions, θ tCch = 0 (see Figure 14), and from the assumptions of the model, items 2.2, 2.3, 5.1, and 6.2. In this problem, as in all previous models, the main one is the kinetic equation describing the change in the state of the target surface. It is easy to write down simplifying (40): However, in [119], the target was described taking into account expressions (40), (42), and (43) by the kinetics of bulk processes: In this case the integral sputtering yield of the target was defined as: where θ tC and θ tM are the surface fractions covered by film and without it, respectively (see Figure 14); S C and S M are the partial sputtering yields of the M m X n compound molecules and M target atoms, respectively. It is obvious that θ tC + θ tM = 1. The correctness of the first part of the model was verified by applying it to the results of studying the angular dependence of the oxidation degree of surface θ tC and subsurface θ tC layers of a silicon target under bombardment with only oxygen ions with an energy of 5 kev. By solving the system of Equation (45), it was shown that the model based on ion implantation is able to adequately describe experiments of this type. This result was achieved by adjusting the reaction rate constant k = 1 × 10 −22 cm 3 ×c −1 .
Further, in [119], using the model, the dependences of θ tC and θ tCs on f = p/(p + p A ) were estimated for sputtering of a silicon target in an Ar + O 2 mixture. It was found that as f increases, the dependences θ tC = f (f ) and θ tCs = f (f ) increase nonlinearly with saturation up to θ tC = 1 and θ tCs = 1. Similar results were also obtained for reactive sputtering of a silicon target in Ar + N 2 .
In the next step in [119], the RSD model was completely formed by including chemisorption and knock-on effect. In accordance with assumption 6 of the model in the steady state, three processes are involved in the formation of the target surface fraction θ tCch covered by the M m X n compound:

•
Formation of M m X n molecules on the fraction of the target surface θ tM by chemisorption; • Removal of M m X n molecules by sputtering; • Removal of X atoms from M m X n molecules due to the knock-on effect. In fact, this process is also understood by the authors of the model as the removal of M m X n molecules.
The addition of the model with chemisorption and the knock-on effect required changing the first equation in system (45).
where p c (x) is the distribution of atoms implanted due to the knock-on effect. In the first approximation, it is assumed that this distribution is also a Gaussian. After formulation of the equations describing the target state in the new model, the authors studied the influence of chemisorption and the knock-on effect on the degree of target oxidation during silicon sputtering in Ar + O 2 . It appeared that without the knock-on effect (β = 0), as the adhesion coefficient increases, the S-shaped dependence θ rb on the reactive gas mole fraction f turns into an exponential with saturation and an inflection point. On the other hand, at a constant value of α = 0.01, an increase in β from 0 to 1 has the opposite influence on the degree of oxidation.
A detailed analysis of the target state in the new model was supplemented by the equation of the steady state of the substrate (wall), which differed from the analogous equation in the Berg chemisorption model (33) in one detail: This detail is visible in the second term on the left side of Equation (48). In square brackets, the fraction θ tC with the M m X n film formed due to implantation is added to the fraction of the i-th surface θ tCch covered with the M m X n compound film formed due to chemisorption. Formula (48) is written under the assumption that both films have the same sputtering coefficient, S C .
Further, the new model was supplemented with equations for reactive gas flows to all surfaces inside the vacuum chamber. The equations for the flow incident on the substrate (chamber wall), pumping out by a vacuum pump, and the balance of gas flows were written in the form (6), (14), and (15), respectively. The flow incident on the target surface included both chemisorbed and implanted components: The first term in (49) describes the amount of gas consumed during chemisorption. The second term specifies the part of the flow from which atoms are implanted into the target as a result of the knock-on effect. Their contribution is subtracted because the model already includes these atoms in the first term. The last term in (49) takes into account the compound molecules formed in the bulk due to implantation.
The main result established in [119] and repeated in [122] is that the dependences p = f (Q 0 ) obtained using the Berg chemisorption model and the RSD model are comparable. However, a difference was noted in the kinetics of gas flow during the first 10 s after turning on the magnetron. The Berg model predicts the transition from the metallic target operation mode to the oxide one with a shorter time constant. Calculations were performed for the sputtering of an aluminum target in Ar + O 2 .

Development of the RSD Model
In subsequent publications [65,[120][121][122][123][124][125][126][127][128][129][130][131][132][133][134][135][136][137][138][139], Depla et al. reported on the development of the original RSD model. For example, in [65], kinetic equations appeared for all parts of the target surface ∂θ tM /∂t, ∂θ tC ch /∂t and ∂θ tC /∂t (see Figure 14), in which all considered physical processes are taken into account. Writing the kinetic equations in a different form did not change anything in the RSD model. This form was more common for writing equations describing the kinetics of processes on the target surface, which was adopted in describing the Berg model. However, a significant amendment has been made to them. If, in all previous cases, one of the main independent variables was the discharge current, then in [65], it was the ion current with the density: where γ i , as in (3), is the ion-induced electron emission yield for the i-th part (i = M, C, C ch ) of the target surface (as in (3)). As it follows from (50), the discharge current, including the electronic component, can be 5-10% higher than the ion current. The authors of the RSD model paid quite a lot of attention to the simulation of the magnetron discharge [116,127,128,137]. To understand reactive sputtering, these problems are of independent importance. A change in the operating mode of the target, in particular, the transition from the metallic mode to the reactive mode and back, was reflected in U = f (j) and U = f (Q 0 ). Jump-like changes appeared in them, corresponding to a change in the operating mode of the target, in which the ability of the surface to emit electrons under the influence of accelerated ions changed.
However, the simulation of the discharge, contributing to a better understanding of the processes of reactive sputtering, was not added to the previously considered models of new physical processes. Therefore, we leave aside the details of these publications and focus only on models that reflect the physicochemical aspects of reactive sputtering. We only note that the authors of these publications, when analyzing the current-voltage characteristics of magnetrons, consider the discharge voltage to be an independent variable. In the physics of a gas discharge, the discharge current is taken as this variable, which corresponds to the physical essence of the process of current flow through a gas.
The development of the RSD model is the subject of article [129], which presents a model of reactive sputtering of a cylindrical magnetron. In it, a cylindrical cathode rotated in a constant magnetic field. The authors of [129] noted that in the steady mode, the cylindrical magnetron behaved similarly to a planar magnetron. However, the change points of the target operation mode depended on the target rotation speed. It was found that the accepted form of the RSD model is insufficient to describe the reactive mode of operation of a rotating target. A more correct description of such a process was achieved by including the re-deposition of sputtered material on the target in the model. A similar correction for the target was proposed in [132,133]. A new modification of the RSD model described in [132,133] is called RSD2013. An interesting work is devoted to reactive highpower pulsed sputtering (R-HiPIMS), [135]. The difference between this method and DC sputtering or RF sputtering was the high proportion of ionization of the sputtered target material. As a result, the target sputtering rate increased. It turned out that in the R-HiPIMS process, the probability of observing hysteresis is much lower compared to DC magnetron sputtering. The reason for this could be related to the implantation of sputtered atoms into the target. However, the authors are right that for any reliable judgment about the features of R-HiPIMS, experimental data, which has been accumulated so far, is not enough.
In [129] and subsequent publications [130,132], the efforts of the authors were aimed at improving the computational procedures for solving problems of reactive sputtering modeling. For example, in [129], the surfaces of the target, substrate, and vacuum chamber are singled out to describe the process of reactive sputtering. Each of them is represented in the form of one or more cells. In this case, the model, as before, took into account the chemisorption of reactive gas molecules on the target, the implantation of chemisorbed particles due to the knock-on effect, and the direct implantation of reactive gas ions. It was pointed out that reactive gas atoms are implanted into the target volume, usually several nanometers below its surface. During direct ion implantation, molecular ions are neutralized immediately in front of the target surface, and upon deceleration in the target, they are split into two atoms. To describe the depth distribution of the concentration of compound molecules, the target was divided into 50 cells. In addition, to take into account the uneven flow of ions incident on the target surface, the sputtering area was divided into 10 × 25 cells.
Finally, let us turn our attention to publications devoted to the observation of double hysteresis in reactive sputtering [136,139]. This effect was observed in [139] when studying the dependence of the partial pressure of a reactive gas or discharge voltage on its flow. Under normal conditions, as follows from all the preceding, the resulting curve is S-shaped, indicating a possible hysteresis in the experiments. The authors of the mentioned works, under certain conditions, observed two S-shaped curves. One of them arose with an increase in the pressure of the reactive gas, the other-when returning to its original state due to a decrease in the pressure of the reactive gas. This behavior has been described as double hysteresis behavior. The origin of the double hysteresis behavior was studied by high-performance computing using a previously developed model. The influence of various process and material parameters has been evaluated on the basis of recently developed measures to characterize process design curves. This high throughput analysis showed that the double hysteresis behavior is related to the difference in the rate of removal of unreacted implanted ions as the reactive gas pressure is increased and decreased. In the parameter space, a region can be defined for which the double hysteresis behavior is strong. The latter can not only help in further experiments to study this behavior but also determine the conditions that limit its influence. It was found that for aluminum, a discharge current density of about 0.025 A/cm 2 provides the maximum double hysteresis.

Implantation
The involvement of new physical processes in the model may become necessary if its initial version is inadequate. In 3.3, the inadequacy of the initial Berg model applied to the experimental results in [25,32] is shown. The reasons for the obtained result, as mentioned earlier, could be both the inaccuracy of the values of the model parameters that were used in the calculation and the incomplete understanding of the physical processes occurring during reactive sputtering. The selection of model parameters is an optimization problem in which the minimum of a function is found in a space of several degrees of freedom. This function, in this case, serves as a criterion of adequacy when it is given by the quadratic form of the proximity between the experiment and the result of prediction by the model. An increase in the number of degrees of freedom, for example, in [130] there were five, significantly complicates the problem. However, if it is not possible to obtain an adequate model in this way, then one should look for additional physical processes or effects that affect the formation of a compound film on the target.
This was done in [75] for the initial Berg model when implantation and knock-on effect were added to it. However, in [119], the authors of the RSD model acted differently. Their initial model contained, as a mechanism for the formation of the M m X n compound on the target surface, the implantation of X atoms accompanied by a bulk chemical reaction and sputtering of the surface layers. The next steps in the development of the initial model were the inclusion of chemisorption and the knock-on effect in it. However, for the reader, the result was unexpected. Recall that the main result established in [119] was that the dependences p = f (Q 0 ) obtained using the Berg chemisorption model and the RSD model were comparable. This is very strange since chemisorption and implantation should have an additive effect on the reactive sputtering process. In the RSD model, each process generates its own area on the target surface covered with the M m X n compound (see Figure 14).
In this regard, let us perform a more detailed analysis of the implantation of X atoms into the target. At the same time, let us pay attention to the result of the study given in [140], which the authors of the RSD model accepted as a reliable confirmation of the correctness of using implantation in the reactive sputtering model. This study was carried out by the method of nuclear reactions after turning off the magnetron. It was found that the maximum number of retained reactive gas atoms significantly exceeded the adsorbed monolayer. From our point of view, the physical reason for such a result could be different. This could be due, firstly, to polymolecular adsorption [141]. Second, the cause could be the diffusion of reactive gas atoms into the target volume. The significance of the second process can be expected to be quite high since, with a target thickness of 5-6 mm, the temperature in the region exposed to a powerful ion flux can reach 700-900 • C.
Let us discuss the primary RSD model based only on implantation. Initially, assuming that the formation of the M m X n film and its sputtering traditionally compete on the target surface, we write the kinetic equation for it in the form where θ tC is the fraction of the target surface covered by the compound (see Figure 14); (dθ tC /dt) impl is the compound film growth rate due to ion implantation; (dθ tC /dt) sput is the sputtering rate of the compound film sputtered by argon ions.
To write the first term in (51), we use the assumption of the authors of the RSD model that only a part of the reactive gas atoms α i implanted in the volume takes part in the chemical reaction with the target metal atoms. Let us write this part of Equation (51) in the form adopted in [75]: The notation adopted in (52) was used in the previous sections of the paper. Applying expression (52) in the model, the user encounters the first hard-to-resolve question: "What value of α i should be taken in the calculations?". The simplest α i = 1 is unsuitable. The reason for this statement can be understood by performing a simple quantitative analysis.
Let us assume that the ion flux X + 2 has a uniform spatial distribution, bombarding the entire surface of the target. Then the result of implantation can be presented in Figure 15, where the distribution of implanted atoms is shown with a gradient color, which has a maximum in the plane denoted by x max .
unsuitable. The reason for this statement can be understood by performing a simple quantitative analysis.
Let us assume that the ion flux 2 X + has a uniform spatial distribution, bombarding the entire surface of the target. Then the result of implantation can be presented in Figure  15, where the distribution of implanted atoms is shown with a gradient color, which has a maximum in the plane denoted by xmax. Figure 15. Distribution of implanted atoms in the target.
If αi = 0, then in Equation (51) θtC ≡ 0, since the bulk chemical reaction does not occur in the target. Let us now assume that αi = 1 and analyze the process of reactive sputtering of a silicon target in Ar + N2 at an ion energy of 0.45 keV, which was considered in [97]. In this case, participation of all implanted nitrogen atoms N in the chemical reaction is allowed. Then Figure 15 will reflect the change in the chemical composition of the nitrogen nitride film in the surface layer of the target. In dark areas, it can be close to the stoichiometric composition of Si3N4; in light areas, it will have nitrogen defects. Assuming the distribution of the concentration of nitrogen atoms N in silicon in the form of the Gaussian distribution (41), we determine the change in the relative concentration of implanted N atoms with time, taking sputtering into account, in the form of an integral with a variable upper limit: If α i = 0, then in Equation (51) θ tC ≡ 0, since the bulk chemical reaction does not occur in the target. Let us now assume that α i = 1 and analyze the process of reactive sputtering of a silicon target in Ar + N 2 at an ion energy of 0.45 keV, which was considered in [97]. In this case, participation of all implanted nitrogen atoms N in the chemical reaction is allowed. Then Figure 15 will reflect the change in the chemical composition of the nitrogen nitride film in the surface layer of the target. In dark areas, it can be close to the stoichiometric composition of Si 3 N 4 ; in light areas, it will have nitrogen defects. Assuming the distribution of the concentration of nitrogen atoms N in silicon in the form of the Gaussian distribution (41), we determine the change in the relative concentration of implanted N atoms with time, taking sputtering into account, in the form of an integral with a variable upper limit: where a is the dimensional constant sputtering rate of the Si 3 N 4 film. The result of integration (53) has the form The parameters of the Gaussian distribution in (53) were determined using the SRIM program: R p = 2.2 nm and ∆R p = 1.1 nm. Figure 16 shows distribution (54) constructed with allowance for an arbitrarily chosen nitrogen nitride film sputtering rate (dx/dt) sput = a = 0.2 nm/s.
where a is the dimensional constant sputtering rate of the Si3N4 film. The result of integration (53) has the form The parameters of the Gaussian distribution in (53) were determined using the SRIM program: Rp = 2.2 nm and ΔRp = 1.1 nm. Figure 16 shows distribution (54) constructed with allowance for an arbitrarily chosen nitrogen nitride film sputtering rate (dx/dt)sput = a = 0.2 nm/s. Despite the fact that the ion with the energy of 0.45 keV was represented in the molecular form 2 N + , in SRIM calculations, Figure 16 describes the change in the concentration of nitrogen atoms N implanted in the target. This remark is based on the assumptions of the RSD model. Recall that it assumes deceleration in the target of an ion that dissociates into two atoms. If an atomic nitrogen ion is set as an ionized particle, then its projected range Rp will be equal to 28 nm at ΔRp = 1.6 nm. The surface concentration of atoms implanted into the target, taking into account sputtering, as Figure 16 demonstrates, will increase. Its change is shown in Figure 17; after 40 s, the process will enter the steady mode when the value c(0, t)/D reaches its maximum value. This can only mean one thing, the target cannot have a metallic mode. After turning on the magnetron, a compound film will be formed on the entire surface of the target, the chemical composition of which will continuously change, tending towards the stoichiometric one. In this case, the film thickness will also remain unchanged in Equation (51) θtC ≡ 1 since all implanted nitrogen atoms participate in the bulk chemical reaction.  Despite the fact that the ion with the energy of 0.45 keV was represented in the molecular form N + RSD model. Recall that it assumes deceleration in the target of an ion that dissociates into two atoms. If an atomic nitrogen ion is set as an ionized particle, then its projected range R p will be equal to 28 nm at ∆R p = 1.6 nm.
The surface concentration of atoms implanted into the target, taking into account sputtering, as Figure 16 demonstrates, will increase. Its change is shown in Figure 17; after 40 s, the process will enter the steady mode when the value c(0, t)/D reaches its maximum value. This can only mean one thing, the target cannot have a metallic mode. After turning on the magnetron, a compound film will be formed on the entire surface of the target, the chemical composition of which will continuously change, tending towards the stoichiometric one. In this case, the film thickness will also remain unchanged in Equation (51) θ tC ≡ 1 since all implanted nitrogen atoms participate in the bulk chemical reaction. From all that has been said, it follows that the inequality αi < 1 is valid in (52), and the final result largely depends on the choice of αi.
The results presented by the authors of the RSD model are interesting, but they do not make it possible to assess the state of the target. For example, it does not matter at all why there was no jump in the dependence θtCs = f (f) and why the model result differed from the experimental one in [97]. The main thing, at first glance, is the ambiguity of how to use the obtained result for modeling since it does not provide information about the main function of the reactive sputtering process p = f (Q0), which was previously accepted in [44,45] and many other publications (see Section 3.1). In [119], it was not possible to find an answer to the question of how the functions θtC = f (f) and θtCs = f (f) can be used to determine the boundary between the metallic and reactive modes of target operation.

Nonisothermal Physicochemical Model (the Barybin Model)
Let us return to the initial isothermal Berg model in [46] based on chemisorption. Chemisorption under conditions of constant temperature is used in all models considered in this work. Although in [65], the authors noted that the Berg model uses "…crude ap-proximations…". However, they did not abandon this model. Section 3.3 of this paper notes two very strong simplifications of the Berg model. They are the replacement of the chemical reaction of MmOn film formation by the chemisorption mechanism and the assumption that the temperatures of all surfaces are equal.
For the last twenty years at Saint Petersburg Electrotechnical University (LETI, Russia), our group has been studying transition metal oxide, nitride, and oxynitride films deposited by reactive magnetron sputtering [142][143][144][145][146][147][148][149]. Initial interest in a traditional magnetron with a well-cooled metal target (cold target) made it possible to create a sputter assembly with a hot target [150][151][152][153] and then with a sandwich target [154][155][156][157][158][159]. The study of the processes that occur during film deposition using these tools inevitably led to modeling [160][161][162][163][164][165]. The starting point in this work was the article [46]. Developing the main assumptions of the Berg model outlined earlier, we proposed to bring the model closer to reality. Firstly, chemisorption was replaced by chemical reactions, and secondly, the limitation on surface temperatures inside the vacuum chamber was removed. The model was named "nonisothermal physicochemical" or the Barybin model after the name of the group leader in the first decade of the 21st century. From all that has been said, it follows that the inequality α i < 1 is valid in (52), and the final result largely depends on the choice of α i .
The results presented by the authors of the RSD model are interesting, but they do not make it possible to assess the state of the target. For example, it does not matter at all why there was no jump in the dependence θ tCs = f (f ) and why the model result differed from the experimental one in [97]. The main thing, at first glance, is the ambiguity of how to use the obtained result for modeling since it does not provide information about the main function of the reactive sputtering process p = f (Q 0 ), which was previously accepted in [44,45] and many other publications (see Section 3.1). In [119], it was not possible to find an answer to the question of how the functions θ tC = f (f ) and θ tCs = f (f ) can be used to determine the boundary between the metallic and reactive modes of target operation.

Nonisothermal Physicochemical Model (the Barybin Model)
Let us return to the initial isothermal Berg model in [46] based on chemisorption. Chemisorption under conditions of constant temperature is used in all models considered in this work. Although in [65], the authors noted that the Berg model uses " . . . crude approximations . . . ". However, they did not abandon this model. Section 3.3 of this paper notes two very strong simplifications of the Berg model. They are the replacement of the chemical reaction of M m O n film formation by the chemisorption mechanism and the assumption that the temperatures of all surfaces are equal.
For the last twenty years at Saint Petersburg Electrotechnical University (LETI, Russia), our group has been studying transition metal oxide, nitride, and oxynitride films deposited by reactive magnetron sputtering [142][143][144][145][146][147][148][149]. Initial interest in a traditional magnetron with a well-cooled metal target (cold target) made it possible to create a sputter assembly with a hot target [150][151][152][153] and then with a sandwich target [154][155][156][157][158][159]. The study of the processes that occur during film deposition using these tools inevitably led to modeling [160][161][162][163][164][165]. The starting point in this work was the article [46]. Developing the main assumptions of the Berg model outlined earlier, we proposed to bring the model closer to reality. Firstly, chemisorption was replaced by chemical reactions, and secondly, the limitation on surface temperatures inside the vacuum chamber was removed. The model was named "nonisothermal physicochemical" or the Barybin model after the name of the group leader in the first decade of the 21st century.

Single Cold Target in Ar + X 2
The first model was developed to study the process of reactive sputtering of a single cold target in an Ar + X 2 mixture [160]. It adopted the following assumptions: 1.
The vacuum chamber contains a target, a substrate, and a wall. The areas of these elements are equal to A t , A s , and A w (or A i , i = t, s, w), respectively; 2.
During film synthesis, the surface temperatures T t , T s , and T w (or T i , i = t, s, w) are different, and the temperature of the gas environment T 0 is equal to the wall temperature T w ; 3.
Every surface can undergo a chemical reaction: Surface reaction (55) proceeds between the adsorbed X 2 reactive gas molecules and the metal on the target surface fractions free from M m X n . In (55), the quantity k(T i ) is the Arrhenius reaction rate constant, which has the dimension of the flux density: where k 0 and E a are the constant and activation energy of the reaction, respectively.

4.
At any time (see Figure 2): • Relative fraction θ t of the sputtered target surface is covered with M m X n . The rest (1 − θ t ) is pure metal M. This state can arise due to the competition of two processes, including the formation of the M m X n film according to surface reaction (55) and its sputtering with argon ions in the form of molecules. The flux sputtered from the target surface includes J tM (pure metal atoms M) and J tC (M m X n molecules); • Due to the surface reaction (55) and fluxes generated by the target, a solid solution M + M m X n is formed on the i-th surface (i = s, w). Let us represent it on each surface as two regions with relative areas θ i and 1 − θ i , containing the compound M m X n and the metal M, respectively.

5.
Each surface consumes reactive gas to maintain surface reaction (55). Let Q i denote the flux incident on the i-th surface (i = t, s, w) and participating in the formation of the M m X n compound on it.
The independent variables in this problem are the discharge current density j and the reactive gas flow Q 0 introduced into the chamber. The main dependent variable of the process is the partial pressure p of the gas X 2 .
Following [46], we will compose a system of algebraic equations describing the kinetics of processes occurring on all surfaces and gas flows inside the vacuum chamber.
The kinetic equation for the target surface, by analogy with (9), takes the form However, unlike (10), the first term in (57) describes the formation of the M m X n film due to the surface chemical reaction (55). In accordance with the law of mass action [166], we set the reaction rate with the expression here N ch is the concentration of chemical reaction centers on the target surface; θ 0t is the fraction of the target surface covered with adsorbed gas molecules. Physical adsorption is described by the Langmuir isotherm, which for nonisotemic conditions can be represented as: In (59), parameter b(T 0 , T t ) is the function of gas temperature T 0 , and target temperature T t : where α 0 is the condensation coefficient of reactive gas molecules; N ph is the concentration of physical adsorption centers τ a ≈ 10 −13 s is the average lifetime of a molecule in the adsorbed state; Q ph is the heat of physical adsorption (cal/mol); R = 1.99 cal/(mol· K) is the universal gas constant. With the second term in (57) (dθ t /dt) sp = −jS C θ t /eN ch , the steady state equation for the target surface takes the form Using (13), the kinetic equation of the steady state of the wall and substrate surfaces can be written as: In (62), θ 0i is given by Equations (59) and (60), in which T t is replaced by T i , i = w, s. The gas flows supporting the chemical reaction on all surfaces are described by the expression: with the dimensional factor c 0 = 2.436 × 10 −18 , (cm 3 ×s)/min, which helps to convert the gas flow from the number of particles per second to cubic centimeters per minute (sccm), standard for technological tasks. When calculating c 0 , as usual, T 0 = 298 K and p 0 =10 5 Pa were chosen as standard conditions. Denoting the volume of one sputtered particle per one minute through V 1 , we obtain: The balance equations for gas flows and pumping out close the system: where S p is the chamber pumping speed (m 3 /s); c = 600, (s×cm 3 )/(Pa×m 3 ×min) is the dimensional factor that converts pascal per cubic meter per second to cubic centimeters per minute under standard conditions. The system of 8 algebraic Equations (61)-(65), taking into account (59) and (60) with respect to the function p = f (Q 0 , j), should be solved numerically, determining the onedimensional dependences p = f (Q 0 ) for j = const or p = f (j) for Q 0 = const. The initial data for calculations are the values of the following parameters: j, α 0 , τ 0 , S M , S C , S p , A t , A s , A w , N ph , Q ph , T t , T s , T w , E a , k 0 . The model makes it possible to estimate the influence of independent variables on the values θ t , θ w , θ s , Q t , Q w , Q s , and Q p . Direct measurement of these variables is not possible.
The proposed model was used to simulate the reactive sputtering of a tantalum target. The experimental results of the process are taken from [32]. The result of the numerical solution of Equations (61)-(65) with allowance for (59) and (60) at a current density j = 0.05 A/cm 2 and parameters from Table 1 is shown in Figure 4 of [160]. This result shows at a qualitative level that the Barybin nonisothermal physicochemical model of reactive sputtering can adequately describe the experimental results. Note. The value of T t = 700 K is an estimate. The values of T s and T w are typical for real processes.

Single Cold Target in Ar + O 2 + N 2
Films of transition metal oxynitrides with composition expressed by the chemical formula M m O n N p in many cases form a continuous series of compounds with different concentrations of nitrogen and oxygen (from M m O n oxide to M m N p nitride) [167].
For the deposition of oxynitride films, the method of DC reactive magnetron sputtering is largely used. This method is most often based on the sputtering of a metal target in a mixture of argon, oxygen, and nitrogen [168][169][170][171][172][173][174]. Let us apply the nonisothermal physicochemical model for such processes. The main assumptions 1 and 2 of this model, set out in Section 5.1, are also valid in this case. Assumptions 3, 4, and 5 should be changed. In this regard, the physical model of this process will be expressed by a number of the following assumptions.

1.
The vacuum chamber contains a target, a substrate, and a wall with surface areas A i , i = t, s, w.

2.
The temperatures of the surfaces T i , i = t, s, w are different, and the temperature of the gas environment T 0 is equal to the wall temperature T w .

3.
Metal oxynitride appears as a result of a surface chemical reaction on A t , A s , and A w surfaces. We assume that this compound arises in the form of a solid solution consisting of M m O n oxide and MN nitride. In this case, the components of the solution are formed due to the occurrence of two independent reactions: were k 1 (T) and k 2 (T) are the reaction rate constants. Here and below, index 1 will correspond to oxide or oxygen, and index 2 to nitride or nitrogen. The accepted chemical formula for nitride in the form of MN corresponds to many compounds of this type (TiN, AlN, TaN, etc.).

4.
At any time: • Relative parts θ t 1 and θ t 2 of the sputtered target surface are covered by M m O n oxide and MN nitride, respectively ( Figure 18). The rest (1 − θ t 1 − θ t 2 ) is pure metal M. This state can arise due to the competition of two processes, including the formation of the M m O n + MN solid solution film according to reactions (66) and (67) were k1(T) and k2(T) are the reaction rate constants. Here and below, index 1 will correspond to oxide or oxygen, and index 2 to nitride or nitrogen. The accepted chemical formula for nitride in the form of MN corresponds to many compounds of this type (TiN,  AlN, TaN, etc.).

At any time:
• Relative parts θt 1 and θt 2 of the sputtered target surface are covered by MmOn oxide and MN nitride, respectively ( Figure 18). The rest (1 − θt 1 − θt 2 ) is pure metal M. This state can arise due to the competition of two processes, including the formation of the MmOn + MN solid solution film according to reactions (66) and (67)  The analytical description of the physical model (items 1-5) is expressed, as in Section 5.1, by a system of equations. The solution of this system will determine the relationship between the dependent and independent variables of the reactive sputtering process. In this case, there are three independent variables: the current density of argon ions on the target j and the incoming gas flows Q0 1 and Q0 2 . The main dependent variables of the process are the partial pressures of oxygen p1 and nitrogen p2. From the further presentation, it will become clear that the process is characterized by several more dependent variables, namely: Qt 1 , Qs 1 , Qw 1 , Qt 2 , Qs 2, and Qw 2 .
The state of the target surface for this process, in contrast to (56), must be described by two kinetic equations:

5.
Each surface consumes reactive gas to maintain surface reactions (66) and (67). Let Q i j denote the fluxes incident on the i-th surface (i = t, s, w) and participating in the formation of oxide M m O n (j = 1) and nitride MN (j = 2).
The analytical description of the physical model (items 1-5) is expressed, as in Section 5.1, by a system of equations. The solution of this system will determine the relationship between the dependent and independent variables of the reactive sputtering process. In this case, there are three independent variables: the current density of argon ions on the target j and the incoming gas flows Q 0 1 and Q 0 2 . The main dependent variables of the process are the partial pressures of oxygen p 1 and nitrogen p 2 . From the further presentation, it will become clear that the process is characterized by several more dependent variables, namely: Q t 1 , Q s 1 , Q w 1 , Q t 2 , Q s 2 , and Q w 2 .
Let us consider the process in more detail and derive equations for finding the dependencies p 1 = f (j, Q 0 1 , Q 0 2 ) and p 2 = f (j, Q 0 1 , Q 0 2 ).
The state of the target surface for this process, in contrast to (56), must be described by two kinetic equations: The kinetics of the compound formation on the target surface in (68) is given by two equations: where N ch j (j = 1, 2) are the concentrations of the chemical reaction centers on the target surface; θ 0t j (j = 1, 2) are the fractions of the target surface covered by the adsorbed molecules of oxygen and nitrogen. Let us describe physical adsorption in the form of the Langmuir isotherm by analogy with (59): In (71), parameter b(T 0 , T t ) is the function of gas temperature T 0 , and target temperature T t : where Q ph j is heat of physical adsorption of the j-th gas. The remaining designations fully correspond to Equation (60). The second terms in (68) are obvious; therefore, following (21) and (22), we write the equations for the steady state of the target surface at the sputtering yields of oxide S C 1 and nitride S C 2 in the form Further, using (62), the kinetic equations for the steady state of the wall and substrate surfaces can be expressed in the form In (75) and (76) θ 0i j (j = 1, 2) are given by Equations (71) and (72), where T t is replaced by T i , i = w, s.
The incoming reactive gas flows Q 0 1 and Q 0 2 have a drain on all surfaces: • Pumped gas flows Q p 1 and Q p 2 , which for j =1, 2, are written by analogy with Equation (64): • Flows to the target Q t 1 and Q t 2 ; • Flows to the substrate, where the compound film is formed from the target material, Q s 1 and Q s 2 ; • Flows to the chamber internal walls, where the target material is deposited, Q w 1 and Q w 2 .
By analogy with (62), for the flows of oxygen and nitrogen to the i-th surface (i = t, s, w), we can write: Now the balance equations of the gas flows (j = 1, 2) take the form The balance equations for gas flows (80) at j = 1, 2 close the system of equations describing the reactive sputtering process under study. The system, taking into account (71) and (72), includes 16 Equations (73)- (80) with the values of subscripts i = t, s, w and j = 1, 2.
In particular cases, at Q 0 1 = 0 or Q 0 2 = 0, the resulting system of equations describes the process of sputtering a metal target in an environment of one reactive gas (nitrogen or oxygen). It is easy to show that, in these cases, the system (73)-(80) reduces to the system of Equations (61)- (65). To do this, it is necessary to substitute O 2 or N 2 as the reactive gas X 2 into the equation of the surface chemical reaction (55).
The system of Equations (73)- (80) was applied to simulate the reactive sputtering of a titanium target in an Ar + O 2 + N 2 gas mixture. The process parameters are given in Table 2. The results of numerical solution of Equations (73)-(80), taking into account (71) and (72) at the current density of j = 200 A/m 2 , are shown in Figure 19.   Figure 19 shows that the proposed model makes it possible to predict the points of change in the operating mode of the target, indicating the hysteresis effect. With an increase in the incoming nitrogen flow, both points of change in the operating mode of the target shift to the region of lower values of the incoming oxygen flow, which is quite expected. In this case, the form of dependences p1 = f(Q0 1 ) changes with increasing Q0 2 .

Single Hot Target in Ar + X2
In recent years, the hot-target magnetron has gained increased attention among specialists [150][151][152][155][156][157][177][178][179][180][181][182]. A feature of this device is that its target can be heated to melting. In a heated state, the target becomes a source of thermoelectrons and evaporated particles. These two fluxes change the conditions of reactive sputtering, so they must be taken into account in the physical model of the process [162,164,165].
Let us compose a model of reactive sputtering of a single hot target. As its basis, by analogy with Section 5.1, we will take several assumptions.
1. The vacuum chamber contains a target, a substrate, and a wall with surface areas Ai, i = t, s, w.  Figure 19 shows that the proposed model makes it possible to predict the points of change in the operating mode of the target, indicating the hysteresis effect. With an increase in the incoming nitrogen flow, both points of change in the operating mode of the target shift to the region of lower values of the incoming oxygen flow, which is quite expected. In this case, the form of dependences p 1 = f (Q 0 1 ) changes with increasing Q 0 2 .

Single Hot Target in Ar + X 2
In recent years, the hot-target magnetron has gained increased attention among specialists [150][151][152][155][156][157][177][178][179][180][181][182]. A feature of this device is that its target can be heated to melting. In a heated state, the target becomes a source of thermoelectrons and evaporated particles. These two fluxes change the conditions of reactive sputtering, so they must be taken into account in the physical model of the process [162,164,165].
Let us compose a model of reactive sputtering of a single hot target. As its basis, by analogy with Section 5.1, we will take several assumptions.

1.
The vacuum chamber contains a target, a substrate, and a wall with surface areas A i , i = t, s, w.

2.
The temperatures of the surfaces T i , i = t, s, w are different, and the temperature of the gas environment is T 0 = T w . The target surface temperature T t depends on the current density j.

3.
On a part of each surface free from the M m X n film, chemical reaction (55) can occur, which is preceded by the physical adsorption of molecules of reactive gas X 2 .

4.
At any moment in time (see Figure 20): • The relative fraction θ t of the sputtered target surface, as in any similar cases, is covered by M m X n . The rest (1 − θ t ) is a pure metal M. The flux sputtered from the target surface includes fluxes J tM (atoms of pure metal M) and J tC (molecules M m X n ) The fluxes of sputtered metal atoms and compound molecules with densities J sp tM (j + ) and J sp tC (j + ) are proportional to the ion current density j + . The corresponding fluxes of evaporated particles with densities J ev tM [T t (j + )] and J ev tC [T t (j + )] depend indirectly on j: where A M , B M , A C , and B C are constants that determine the saturated vapor pressure of the metal M and the compound M m X n in pascals, respectively; m M and m C are the masses of the metal atom M and the M m X n molecule, respectively; • On the i-th surface (i = s, w), due to the surface reaction (55) and fluxes generated by the target, a solid solution M + M m X n is formed. Let us represent it on each surface as two regions with relative areas θ i and 1 − θ i containing the compound M m X n and the metal M, respectively (see Figure 20). ; (81) [ ] The fluxes of sputtered metal atoms and compound molecules with densities

5.
A flow of reactive gas Q i falls on each surface, supporting the formation of the compound M m X n according to reaction (55). 6.
The discharge current density is where γ is the potential ion-induced secondary electron emission yield; j + is the ion current density; j − is the current density of thermal electron emission (described by the Richardson-Dushman equation): where A ≈ 120 A·cm −2 ·K −2 ; ϕ t is the work function of electrons for the target material; k is the Boltzmann constant. 7.
The influence of the temperature on the parameters of the sputtering process (sputtering yield and ion-induced electron emission yield, etc.) is insignificant.
The independent variables in this problem, as before, are the discharge current density j and the reactive gas flow Q 0 introduced into the chamber. The main dependent variable of the process is the partial pressure p of the gas X 2 .
Let us describe the physical model of reactive sputtering, expressed by assumptions 1-7, by a system of equations. It, as in all cases, should include two groups of equations that describe the steady state of three surfaces and gas flows in the vacuum chamber, respectively.
The kinetic equation for the target surface has the form: The first term in (87) describes the formation of the M m X n film according to reaction (55) and has the form of (58). Let us repeat here (58) for the convenience of readers: In Equation (88), θ 0t is given by Equations (59) and (60). The second term in (87) equals where for the compound M m X n : γ C is the potential ion-induced secondary electron emission yield; S C is the sputtering yield; j − C (T t ) is the current density of thermal electron emission defined by Equation (86). The factor [j − j − C (T t )]/[1 + γ C ] in the right part of (89) sets the density of ion current on the target. The third term in (87) is defined by the Hertz-Knudsen equation in the form [183] Substituting (88), (89), and (90) to Equation (87), we can write the kinetic equation for the target steady state (dθ t /dt = 0): Let us express the compound flux density J tC in (91) in more detail using (82), (89), and (90): The steady state equation for the surface of the substrate and the wall during sputtering of a hot target has a form that is similar to Equation (62). Figure 20 shows its difference; that is, when describing the processes on these surfaces, it is necessary to take into account the fluxes of evaporated particles. Therefore, we write (62) in the form In (93), J tM is written by analogy with (92): where for metal M: γ M is the potential ion-induced secondary electron emission yield; Pumping out the vacuum chamber (65).
The system of algebraic equations describing the model of reactive sputtering of a hot target includes Equations (63)-(65), (91) and (93). If target heating initiates insignificant thermal electron emission and evaporation, then the result of solving system (63)-(65), (91) and (93) will not differ from that predicted by the reactive sputtering model of a cold target (61)-(65) taking into account (59) and (60).
As an example, the system (63)-(65), (91) and (93) was used in modeling reactive sputtering of a hot titanium target in an Ar + N 2 mixture [165]. The solution is fulfilled for the function p = f (j, Q 0 ). The chemical reaction in this problem is described by Equation (63). In the calculations, the parameters indicated in Table 3 were used. The measurement of the target temperature is of great importance. Most often, noncontact infrared pyrometric methods are used in this case [189][190][191]. In this work, the dependence of the target temperature on the discharge current density was determined by a technique based on measuring the optical spectra of the discharge and separating from them a component associated with the thermal radiation of a heated target [192,193]. The results of this computational procedure are shown with dots in Figure 21.
Next, we perform an analysis of the features of the hot target sputtering process, which were revealed by the proposed model. For comparison, a similar problem was solved for a cold target. results of this computational procedure are shown with dots in Figure 21.
Next, we perform an analysis of the features of the hot target sputtering process, which were revealed by the proposed model. For comparison, a similar problem was solved for a cold target. An analysis of the dependences p = f (j, Q0) obtained using the proposed model (given in [165]) showed that for a hot target: • Hysteresis persists over a wide temperature range. In accordance with expression (95), current densities of 50 and 500 A/m 2 corresponded to temperatures of 720 and 1560 K; • Changes in the operating mode of the target occur at lower values of Q1 and Q2 compared to a cold target ( Figure 22); • The width of the hysteresis loop also has smaller values ( Figure 23).  An analysis of the dependences p = f (j, Q 0 ) obtained using the proposed model (given in [165]) showed that for a hot target:

•
Hysteresis persists over a wide temperature range. In accordance with expression (95) Next, we perform an analysis of the features of the hot target sputtering process, which were revealed by the proposed model. For comparison, a similar problem was solved for a cold target. An analysis of the dependences p = f (j, Q0) obtained using the proposed model (given in [165]) showed that for a hot target: • Hysteresis persists over a wide temperature range. In accordance with expression (95), current densities of 50 and 500 A/m 2 corresponded to temperatures of 720 and 1560 K; • Changes in the operating mode of the target occur at lower values of Q1 and Q2 compared to a cold target ( Figure 22); • The width of the hysteresis loop also has smaller values ( Figure 23).  The results of a more detailed analysis are shown in Figure 23, which presents the dependences ΔQ = Q1 − Q2 = f (j) for cold (CT curve) and hot (HT curve) targets. In addition, Figure 23 shows a similar dependence (dashed curve) obtained by solving Equations (63)- (65), (91), and (93), from which the terms taking into account evaporation were removed. Figure 23 demonstrates that at j < 400 A/m 2 (corresponding to a target temperature of 1450 K), evaporation does not influence the sputtering process.
The presented results confirm that the proposed model can predict the reactive sputtering process of a single hot target.
In [163], we published the results of modeling the reactive sputtering process of a single hot target in a gas mixture of Ar +O2 + N2. The equations describing this model are The results of a more detailed analysis are shown in Figure 23, which presents the dependences ∆Q = Q 1 − Q 2 = f (j) for cold (CT curve) and hot (HT curve) targets. In addition, Figure 23 shows a similar dependence (dashed curve) obtained by solving Equations (63)-(65), (91) and (93), from which the terms taking into account evaporation were removed. Figure 23 demonstrates that at j < 400 A/m 2 (corresponding to a target temperature of 1450 K), evaporation does not influence the sputtering process.
The presented results confirm that the proposed model can predict the reactive sputtering process of a single hot target.
In [163], we published the results of modeling the reactive sputtering process of a single hot target in a gas mixture of Ar + O 2 + N 2 . The equations describing this model are similar to (73)- (80), taking into account the features of sputtering of a hot target expressed by formulas (81)-(86). We will not delve into the details of this model. Let us present only some results of calculations for the model of hot titanium target sputtering.
The system of equations describing the model was solved for incoming oxygen flows Q 01 from 0 to 10 sccm, constant nitrogen flows Q 02 0.1, 0.2, and 0.3 sccm, and current densities from 200 to 400 A/m 2 . All dependencies, p = f (Q 01 ), were S-shaped, characteristic of the hysteresis effect. For comparison, solid lines in Figure 24 show the dependencies p = f (Q 01 ) Q02 = const , calculated from the cold target sputtering equations (see expressions (73)-(80)). It can be seen from Figure 24 that when a hot target is sputtered, the addition of nitrogen shifts to the left the points of change in the operating modes of the target, at which the derivative dp/dQ 01 undergoes discontinuities. In addition, an increase in nitrogen incoming flow leads to a decrease in the width of the hysteresis region.
dependences ΔQ = Q1 − Q2 = f (j) for cold (CT curve) and hot (HT curve) targets. In addition, Figure 23 shows a similar dependence (dashed curve) obtained by solving Equations (63)-(65), (91), and (93), from which the terms taking into account evaporation were removed. Figure 23 demonstrates that at j < 400 A/m 2 (corresponding to a target temperature of 1450 K), evaporation does not influence the sputtering process.
The presented results confirm that the proposed model can predict the reactive sputtering process of a single hot target.
In [163], we published the results of modeling the reactive sputtering process of a single hot target in a gas mixture of Ar +O2 + N2. The equations describing this model are similar to (73)- (80), taking into account the features of sputtering of a hot target expressed by formulas (81)-(86). We will not delve into the details of this model. Let us present only some results of calculations for the model of hot titanium target sputtering.
The system of equations describing the model was solved for incoming oxygen flows Q01 from 0 to 10 sccm, constant nitrogen flows Q02 0.1, 0.2, and 0.3 sccm, and current densities from 200 to 400 A/m 2 . All dependencies, p = f(Q01), were S-shaped, characteristic of the hysteresis effect. For comparison, solid lines in Figure 24 show the dependencies p = f (Q01)Q 02 = const, calculated from the cold target sputtering equations (see expressions (73)- (80)). It can be seen from Figure 24 that when a hot target is sputtered, the addition of nitrogen shifts to the left the points of change in the operating modes of the target, at which the derivative dp/dQ01 undergoes discontinuities. In addition, an increase in nitrogen incoming flow leads to a decrease in the width of the hysteresis region.

Sandwich Target in Ar + X 2
The sputtering unit of the magnetron containing several metal plates with cut-outs on the same axis is called the "sandwich target". Let us construct a physical model for reactive sputtering of the simplest sandwich target with two plates using several assumptions [194].

1.
The vacuum chamber includes four elements: the internal and external plates located in the sputtering unit on the same axis [194], the substrate, and the wall. The areas of these elements are denoted by A i , i = t 1 , t 2 , s, w. Cut-outs with a total area of A co are made in the external plate. The values A t 1 and A t 2 specify the areas of only the sputtered fractions of the plates. On the external plate, this fraction is ring-shaped with geometric dimensions determined by the magnetic field. If the area of the ring is A, then the actual sputtered area of the external plate is A t 2 = A − A co or in relative units A t 2 /A = (A − A co )/A = 1 − A co /A = 1 − δ co . Furthermore, the area of the sputtered region of the internal plate A t 1 in this design is determined by the equality A t 1 = A co , and its relative value is A t 1 /A = A co /A = δ co .

2.
The temperatures of the surfaces T i , i = t 1 , t 2 , s, w are different. The temperature of the gas environment T 0 is equal to the wall temperature T w . The internal plate works in cold mode, and the external one in hot mode.

3.
Two reactions can take place on each surface: where m j and n j (j = 1, 2) are the stoichiometric coefficients; k j (T i ) (j = 1, 2) are the Arrhenius equation for the reaction rate constants, which have the dimension of the flux density by analogy with (56); 4.
The discharge current density where j t 1 is the component related to the internal plate, and it is determined by the fraction of the cut-out area δ co = A co /A in the sputtered region of the external plate: Current share j t 2 in (98) is equal to

5.
Two processes occur on the surface of the external plate, which should be taken into account in the model: • Thermal electron emission (the Richardson-Dushman Equation (86)); • Evaporation (the Hertz-Knudsen equation by analogy with (90)).

6.
At any moment in time: • Fraction θ t 1 of the sputtered region of the internal plate is covered by the M 1m 1 X n 1 film. The rest (1 − θ t 1 ) is pure metal M 1 ; • In the sputtered region of the external plate, the similar regions for M 2m 2 X n 2 and M 2 will be denoted by θ t 2 and (1 − θ t 2 ), respectively; • On the i-th surface (i = s, w) due to reactions (96) and (97) and fluxes generated by the plates, a solid solution M 1 + M 2 + M 1m 1 X n 1 + M 2m 2 X n 2 is formed. On each surface, the solid solution is represented in the form of four regions with relative areas θ i 1 , θ i 1 , θ in 1 , and θ in 2 containing the corresponding components of the solution. In this case, it is obvious that θ i 1 + θ i 2 +θ in 1 + θ in 2 = 1.

7.
Each i-th surface (i = t 1 , t 2 , s, w) consumes the reactive gas to support reactions (96) and (97). Let Q ij denote the fluxes incident on the i-th surface (i = t 1 , t 2 , s, w) and participating in the formation of the j-th compound (j = 1, 2).
The independent variables in this problem are the flow of reactive gas introduced into the chamber Q 0 , the discharge current density j, and the relative area of cut-outs in the external plate δ co . The main dependent variable, as in all previous models, is the partial pressure p of the gas X 2 . Next, we present an analytical description of the physical model as a system of equations.
The system of equations includes two groups: • The first of them contains equations describing the kinetics of processes occurring on all surfaces specified in the model; • The second group consists of equations describing reactive gas flows to all surfaces and equations for pumping and balance of gas flows.
On each surface, the formation of M jm j X n j (j = 1, 2) and their removal proceed according to different schemes. Removal of the film from the internal plate occurs only by sputtering. On the external plate, the film removal process enhances evaporation. The steady state of the plates is described by equations like (91), in which θ t 1 and θ t 2 are unknown. These equations are independent and can be easily solved analytically: θ t 2 = f t 2 [k 2 (T t 2 ), θ 0t 2 , J tC 2 ].
Similar values of θ in 1 and θ in 2 are used in the steady state equations of the i-th (i = s, w) surface of type (93). These equations are easily solved for unknowns: θ i 2 = f in 2 [k 2 (T i ), θ 0i 2 , θ i 2 , A, A co , A s , A w , J tM 1 , J tM 2 , J tC 1 , J tC 2 , θ t 1 , θ t 2 ].
In (101), (103), (104), the fluxes of the compound M 1m 1 X n 1 and metal M 1 with densities J tC 1 and J tM 1 , respectively, are described by Equations (92) and (94) for the case with a single cold target (without the second term). In (102)-(104), the fluxes of the compound M 2m 2 X n 2 and metal M 2 with densities J tC 2 and J tM 2 , respectively, are described by Equations (92) and (94) in full form.
As follows from (103) and (104), the values of θ i 1 and θ i 2 are related to (101) and (102) for the values of θ t 1 and θ t 2 . Equations (101)-(104) allow you to find all six fluxes, which are calculated by analogy with (78). In addition, the pumped-out Q p and the total flow Q 0 should be specified in the form (77) and (80), respectively.
Thus, the physicochemical model of reactive sputtering of a sandwich target is described by a system of 14 algebraic equations. Six of them set the steady state of the surfaces of two plates, the substrate, and the wall. The remaining eight set the reactive gas flows in the vacuum chamber. More details about this system can be found in [194].
The described model was used to study the process of film deposition of a binary solid solution of TiO 2 and Ta 2 O 5 oxides. Its chemical formula was written as Ti x Ta 1-x O. These films are of interest in many areas of technology [195][196][197][198][199]. In this case, the sputtering unit consists of an internal titanium plate and an external tantalum plate.
In this problem, reactions (96) and (97) When performing calculations on the optical spectra of the discharge, the dependence of the effective temperature T t 2 (j + t 2 ) on the ion current density of the external plate j +t 2 was determined [200]: T t 2 (j +t 2 ) = 293 + 1520 1 − e −0.00439j +t 2 .
The results of calculations that were performed with the values of the parameters from Tables 1, 2 and 4 are shown in Figure 6 of [194]. Table 4. Additional parameters of the model for reactive magnetron sputtering of a hot tantalum target in Ar + O 2 .
Parameter γ 2 γ n 2 ϕ 2 , eV ϕ n 2 , eV T t 2 Value 0.13 [186] 0.057 [186] 4.5 [201] 5.15 [202] f (j + t 2 ) Parameter A 2 B 2 , K A n 2 B n 2 , K Value 9.2 [187] 17260 [187] 15.3 [203] 20,000 [203] Materials 2023, 16, 3258 42 of 50 The dependences p = f (Q 0 ) in [194] are given on a linear and logarithmic scale, which made it possible to reveal the influence of the internal titanium plate on the process. This effect is observed in the sections of the curves in the region of low gas incoming flow Q 0 = 10 −3 -10 −1 sccm, where the internal plate makes the transition from the metallic to the oxide mode. The change in the operating mode of both plates at different values of Q 0 is associated with a different rate of removal of oxide films from their surfaces. Tables 1 and 2 show that S Ta 2 O 5 > S TiO 2 . Second, it can be shown that, at the discharge current density of more than 700 A/m 2 , the sputtering of the Ta 2 O 5 film is supplemented by its evaporation. Therefore, the change in the steady state of the surface of the internal plate occurs at a lower rate of the chemical reaction proceeding on it and, hence, at a lower incoming oxygen flow.
The proposed model made it possible to calculate the dependence of the stoichiometric coefficient x = f (Q 0 , δ co )| j=100A/m 2 in the formula Ti x Ta 1-x O (see Figure 12 of [194]).
In [194], it is shown that there is an area Q 0 with a negative derivative and a complete transition of the process to the steady mode with a constant value of x.
Concluding the review of works devoted to the modeling of reactive sputtering based on chemical reactions occurring under nonisothermal conditions, we note the main thing in this class of models.
Let us start with the disadvantage. The mathematical description of the model is complicated for application by the presence of three unknown parameters in the description of the kinetics of the chemical reaction. Moreover, if the value of the adsorption coefficient can be taken equal to unity, then the parameters of the chemical reaction constant can be determined only from experimental results by solving an optimization problem.
An advantage of the nonisothermal physicochemical model is a more correct description of the process of reactive sputtering and, as a result, the adequacy of the model, which does not require the involvement of additional mechanisms occurring on the target surface.

Summary
The study of articles published over the past 50 years on reactive sputtering and modeling of this process allowed us to identify the main stages and directions of the research development.
Reactive sputtering is a well-known technological method, the most common in various forms of magnetron sputtering. It is used for the deposition of films of simple metal compounds (nitrides, oxides, carbides, etc.) and their diverse double, triple, etc. solid solutions.
In early experimental works, a significant influence of the reactive gas concentration in the gas mixture or its partial pressure on the film growth rate, discharge voltage, and film composition was found. In addition, nonlinear effects were discovered, and it was found that with an increase in the partial pressure of the reactive gas, there is a critical value. At this value, during the sputtering process, an abrupt change in the film deposition rate was observed, in some cases, by an order of magnitude.
In later works, it was found that the partial pressure of the reactive gas was not an independent variable. It reflected only the state of the process with given values of other parameters that could be changed independently. The main ones were the incoming reactive gas flow and the discharge current (or the power released on the target). In a large number of experiments, when changing the flow or discharge current (power), the hysteresis was discovered.
For a more detailed study of reactive sputtering, many experts worked on the development of their models. All the variety of known models for reactive sputtering was actually based on two assumptions: there are two processes competing on the excited target surface (the formation of a thin layer of metal compound with reactive gas and its removal by accelerated argon ions); on the substrate and walls of the vacuum chamber, the sputtered target material is deposited, and the reactive gas molecules are chemisorbed.
At the initial stage, when constructing models, which are called "specific" in this work, the partial pressure of the reactive gas was taken as an independent variable, and only the processes occurring on the target were considered without taking into account its temperature. The difference between the models of different authors was not fundamental. The main unifying element of all cited publications was the kinetic equation for the target surface.
Further, specific models of reactive sputtering were developed. They included the surfaces of the vacuum chamber wall and the substrate. Moreover, this entailed a complication of the analytical description, in which equations for gas flows to each surface and an equation for the balance of gas flows appeared.
Specific models based on chemisorption received the most consistent development in the work of Professor S. Berg and co-workers. These models eventually became known as the Berg model, and in this work, they are called "general". In the initial Berg model, he described the simplest process of reactive sputtering of a single metal target in a mixture containing one reactive gas. Subsequently, models of more complex processes were proposed. They included, for example, sputtering of a single target in a mixture containing two reactive gases, two targets in a mixture containing one reactive gas, etc. In the latter version, the authors included implantation and the knock-on effect in the model, excluding sputtering of the compound from the target surface in molecular form, which changed the kinetic equations for all surfaces.
Ease of study and use led to great interest in the Berg isothermal chemisorption model over the last twenty years. Specialists used it to describe practical problems. Some of them offered options for its development. In the series of publications on modeling of reactive sputtering, there were works describing other models based, for example, on the laws of thermodynamics, statistics, etc.
The second group of models was proposed by the group of Professor D. Depla. This team has consistently studied reactive sputtering and its modeling over the past twenty years.
A new reactive sputtering model called the RSD model, has been introduced stepwise. In the first step, the authors took into account only the process of implantation of reactive gas ions, which is accompanied by a bulk chemical reaction and sputtering. In this problem, the main one was the kinetic equation describing the change in the state of the target surface. In the next step, the RSD model was fully formed by including chemisorption and the knock-on effect. Subsequently, the authors aimed to improve the computational procedures for solving problems of reactive sputtering modeling. Thus, the RSD2009 and RSD2013 models, etc., appeared.
In fact, the modified Berg model has become an analog of the RSD model in its full form, including implantation and knock-on effect. The exception was the assumption about the features of compound sputtering from the target surface. In the RSD model, it took the molecular form. This greatly distinguished it from the Berg extended model.
Finally, the models developed by our scientific group formed the third approach to the problem. These models differ from all previous ones as they replace chemisorption with surface chemical reactions and remove restrictions on surface temperatures inside the vacuum chamber. Our model was called "nonisothermal physicochemical" or the Barybin model. The initial model described sputtering of a single metal target in a mixture with one reactive gas. The chemical reaction was represented using the Langmuir monomolecular adsorption isotherm and the law of mass action. In subsequent works, the Barybin model was extended to describe more complex sputtering processes in an environment containing two reactive gases, sputtering of a single hot target, and a sandwich target. Here it is necessary to pay attention to the fact that chemisorption models cannot be applied to the last two versions of the magnetron.
Concluding the article, I would like to express the hope that the modeling of reactive sputtering is far from complete. There are complex processes involving more than two magnetrons. The sandwich target sputtering unit may contain more than two plates with variations in their chemical composition. The hot target of the magnetron can operate in a liquid state. There are many practical problems of this kind. To solve them with the help of chemisorption models, first of all, it is necessary to abandon isothermality. In addition, none of the currently known models has studied the effect of diffusion of reactive gas atoms into the target. This process can become important, especially in hot targets. The possibility of including polymolecular adsorption in the model is also not ruled out.
Funding: This research received no external funding. Data Availability Statement: Data is contained within the article.