Simulation of Ionic Resistivity of Ultra-Thin Passivating Films on Metals During their Galvanostatic Transient Growth

Para o caso de crescimentos transientes de filmes de passivação em condições galvanostáticas foram formuladas e resolvidas, por método numérico, equações diferenciais descrevendo a variação das concentrações de defeitos em filmes ultra-finos. Nas equações levou-se em conta a injeção de defeitos, a recombinação de alguns destes e o aumento na espessura do filme. Nos resultados das simulações as concentrações passam por um máximo durante o crescimento, com este máximo aumentando com a corrente e tornando-se independente da mesma para altos valores desta. Para que se obtenham valores razoáveis da concentração de defeitos nos filmes, o trabalho mostra que nem toda a corrente deve ser utilizada para a geração de defeitos. Nessas condições, os cálculos mostram que a resistividade iônica específica nos filmes passa por um mínimo com o aumento de sua espessura e o sobrepotencial neles varia quase linearmente com a mesma. A análise dos transientes leva à idéia de que a constante de velocidade de recombinação de defeitos poderia depender da densidade de corrente utilizada. Todos estes dados são analisados e discutidos do ponto de vista físico.


Introduction
There are practically no works published in the literature with respect to simulation of differential equations related to the growth of passivating films under transient conditions, before coming to stationary or quasistationary situations. 1,2 One of the problems is that, often, the differential equations are non-linear and, as a consequence, have no analytical solutions. The other problem is that under transient conditions there are sound reasons to think that the film parameters are not constant. Perhaps these are the two problems wich explain why no differential equations have been proposed for this phenomenon, except those of D'Alkaine and co-workers. [3][4][5] The present work is the first of a series in the theoretical area making simulations through numerical calculations. It will be dedicated to show the evolution of the defect concentrations, ionic specific resistivity and overpotential inside a passivating film growing under transient galvanostatic conditions. The fact that some defects can recombine will be considered.
For transient conditions, which have not been discussed theoretically in the literature, two cases 2 are possible: a) the ultra-thin film model, in which it is not conceivable that there are concentration gradients inside the film, due to its very thin thickness and high field; b) the thin film model, in which the thickness of the film makes it necessary to consider that the injection of defects in both interfaces (the metal/film and the film/solution interfaces) generates concentration gradients in opposite directions inside the film. The present paper is related to the first model.

The Equations
Considering a metal/initial film/solution interface through which an anodic galvanostatic current density i g will pass, with the initial film homogeneously covering the metal with an initial constant thickness l f 0 and having an initial defect concentration c i 0 , the i g will introduce into the film during growth, for example, interstitial metal ions (M 2+ i,(ox) ) at the metal/film interface and metal ion vacancies ) at the film/solution interface, this latter due to the introduction of O 2-(ox) in the anionic lattice of the oxide. The related reactions at the film interfaces can be described as: where the second one can be considered at equilibrium, if the superficial pH of the solution in front of the growing oxide is constant, 3 M 2+ (m) is a metal ion inside the metal, V (m) is a metal vacancy inside the metal and H 2 O and H + are considered at the solution side in front of the oxide surface.
Under these circumstances the variation of any defect concentration (c i ) inside the film, considering the case of ultra-thin films, will be done by: where V f is a constant giving the volume per unit charge of the film, q f is the charge density of the film, t is the elapsed time from the beginning of the experiment and n i is the number of moles of the defect i inside the film per unit area (c i = n i /V f q f ). The first term on the right hand side of equation 3 is related to the increase in defect concentration due, in first approximation, to the passage of i g . The second term is related to the decrease of i concentration with the increase in film thickness. Under galvanostatic conditions, q f is given by: where q 0 corresponds to l f 0 which is equal to V f q 0 . Inside the film, the recombination reaction can be considered, which can be given by: where k r and k d are the recombination and dissociation rate constants, respectively. k d corresponds to the formation of a Frenkel pair. The reactions in both directions will contribute to the variation of c i .
As a consequence, (dn i /dt) will be given by the injection of defects, the formation of Frenkel pairs and the recombination of them. Therefore, this term must be expressed by: where z is the charge of the i defect, n 0 is the moles of nodes of cations or anions in the lattice per unit film area and n i and n v are the concentrations of interstitial and vacancy defects, also per unit film area. From equations 3, 4 and 6 it is possible to deduce the final differential equation for the variation of the defect concentrations c i with q f : where c 0 is the lattice concentration of the cation or anion nodes of the film (equal to d/M, with d the density of the film and M its molar mass). In equation 7, c i has been made equal to c v and is represented as c i . Equation 7 is nonlinear and has no analytical solution, but its resolution would give the evolution of c i with film growth, expressed as q f . The integration will be done numerically.

Methods
For the numerical approach the Maple applicative (Maple V, Release 4) and the Runge-Kutta method was used. From this applicative the possibility to represent direction fields and vector fields of the differential equation to detect qualitatively the behavior of the solutions was also considered.
For the calculations the corresponding values of the constant parameters were: q 0 of 0.44 mC cm -2 ; c 0 of 0.069 mol cm -3 ; c i 0 of 6 x 10 -3 mol cm -3 (about 10% of c 0 ); z equal 537 Simulation of Ionic Resistivity of Ultra-Thin Passivating Films Vol. 13, No. 4, 2002 to 2; V f of 7.52 x 10 -5 cm 3 C -1 . The theoretical value of ZnO and, consequently, d and M were selected 6 as 5.61 g cm -3 and 81.38 g mol -1 , respectively. In this way, the initial film thickness was 0.33 nm and the simulation experiments were run up to a maximum thickness of about 3.0 nm (q f of 4 mC cm -2 ) so that the calculations correspond to ultrathin film situations.
The values of i g , k r and k d will be varied to study their influence, as will be pointed out in each simulation experiment.

Results and Discussion
When equation 7 is integrated, the c i values always pass through a peak due to the fact of injection and recombination of defects, together with the increase of the film thickness 7 (see Figure 1). The problem is that for normally acceptable values of the parameters the results give a very high c i peak, practically of the order of c 0 . This is shown in Figure 1 curve (a) for k d 0.3 s -1 , k r 10 4 cm 2 mol -1 s -1 and i g 0.1 mA cm -2 . Reasonable variations of i g or of c i 0 do not permit to resolve the appearance of these high values of c i at the peak. Mathematical analysis shows that the high values come from the high value of c 0 which influences the first term of equation 7, and which can not be changed to resolve the problem without loosing the physical meaning of the description (the true value of c 0 must be of the order of d/M for the bulk values). Even if in the first term of equation 7 k d is reduced to the point at which it has no more effect ( Figure  1 curve (b), k d reduced to 3 x 10 -3 s -1 ), the problem can not be eliminated. Another way to resolve the problem is to increase the recombination rate, which will increase the third term of equation 7, but this case quickly arrives at a situation in which the maximum happens at very low thicknesses or times (see Figure 1 curve (c), for k d equal to 3 x 10 -3 s -1 and k r equal now to 10 6 cm 2 mol -1 s -1 ), being contrary to general normal experimental results 7 for different systems. This means that the only possible way to come out from these contradictions is to consider that not all the current must inject defects inside the film. Some of the current produces point film defects and the rest directly gives the normal constituents of the film. This means introducing a factor K in equation 7, corresponding to the fraction of i g which is used to produce the point defects. This constant will multiply the term on 1/q f in the first term on the right of equation 7. In the following a hypothetical value of K of 0.1 will be used to explore this idea.
In Figure 2 the evolutions of the defect concentrations given by equation 7 are plotted with the introduction of a K factor equal to 0.1, for different i g . There are several interesting theoretical points which can be indicated. The first one is that now the defect concentrations are acceptable values even at the maximum. The second one is related to the fact that the increase of i g has, as a consequence, an increase of the c i values as must be qualitatively expected. The third one is the fact that higher the i g , the higher the thickness for the recombination to have an effect producing the maximum. This is an unexpected result. In general, experimental galvanostatic results 7 have shown that the increase of i g always has as a consequence a decrease of the thickness for which this maximum or its experimental equivalent is obtained. The equivalent experimental fact 7 is the inflection point in the E/t or E/q f transient.  This result shows the importance to go from qualitative explanations to simulations to see if the "qualitative" thinking is followed by the mathematical formulations. In relation to the topic under analysis, the problem is to make a conclusion between evolution of defect concentrations and evolution of potential vs. charge density, or ionic specific resistivity vs. charge density curves, as will be discussed further. The fourth interesting point is that there is always a region of high current densities for which the c i results become independent of i g . For a long time attention 3,5 has been called to the existence of at least two kinds of stationary or quasi-stationary states in relation to passivating film growths. One would exist for long time experiments (the stationary situation related in most of published papers). The other would exist for high rates of growth, where aging has no time to happen and, consequently, the film is always practically the same. The latter appears in the simulations at very high i g , as practically an independent curve for c i vs. q f (Figure 2 curves (c), (d) and (e)).
From the calculations of c i given above it is possible now to obtain the evolution of the specific ionic resistivity r f inside the film for different q f . r f will be given by where m j corresponds to the ionic mobility of the defect j. For the present calculations all c j and all m j will be taken the same for all j's.
In Figure 3 are plotted the results of r f versus q f for the different transients of Figure 2 using a m i value 2 of 8 x 10 -12 cm 2 V -1 s -1 . This Figure can be compared with experimental results and even the parameter values can be selected to make that the theoretical and experimental plots become compatible, but this will be the objective of a forthcoming work. 8 Here it is important to point out that the r f values pass through a minimum and that this minimum decreases with the increase of i g . Both results are in agreement with what is found experimentally. 3,7,8 Nevertheless, in Figure  3, as for the case of the c i vs. q f plots, the increase of i g shifts the minimum to higher q f , against the experimental results. 7,8 From the point of view of the present description this is a clear indication that the variation of i g in equation 7 must be followed by a variation in k r , the only mathematical way in the equations to make the experimental facts found possible (the minimum must shift to lower q f values for higher i g ). Variations of m i or k d can not lead to compatibility of the theoretical and experimental results. This will be shown in detail in a forthcoming work. 8 The variation of k r with i g is possibly related with a variation of the film characteristics with i g , for example, the hydration level. This new result shows for the second time (the first was the introduction of K) the importance of the simulation technique for the development of theoretical analysis and not only for the interpretation of experimental results. A final fact which is possible to be theoretically obtained is the overpotential through the film (h f ). This dependent variable corresponds to the potential drop inside the film and is a component of the total potential (E) since this is a sum 9 of the overpotential at the metal/film interface (h m/f ) plus h f and E F (the Flade potential) for cases in which the pH of the solution is constant and the system is in the passivity region. 5 The h f is given in an Ohmic model by 3,4 h The results calculated from equation 10 for the variation of h f with q f are shown in Figure 4 for the same conditions as for Figure 2. The first point is that all the curves extrapolate to zero, as must be expected physically and from equation 10. The second point is that with the increase of i g , for the same q f (the same film thickness) there are higher fields (higher h f ). The last point is that for low q f there is some curvature which decreases with the increase of i g . This is related with the weighting of variation of r f , as shown in Figure 3. These facts are in agreement with a general understanding of the Ohmic model, 5 considering that the model defining parameters (in this case only r f ), at low i g , will have time to present aging, before coming to the quasistationary situation. At high i g there is no time for aging. The aging here must be understood fundamentally in terms of the recombination process.

Conclusions
Using differential equations to describe the evolution of the point defect concentration inside a passivating film growing galvanostatically under ultra-thin film conditions, it was possible to obtain the general physical description of the process, the evolution of the defect concentrations, the ionic specific resistivity and the overpotential through the film.
The results show that during the transients the point defect concentrations pass through a maximum and that this maximum increases with the current density, as must be expected. Nevertheless, to obtain reasonable values for the point defect concentrations at the maximum it is necessary to consider that not all the current is used to produce point defects.. The major portion of the current, for this situation, is used to directly produce the normal lattice of the film. On the other hand, the defect concentration maximum shifts higher thicknesses when the current density is increased, differing from prior experimental results. 7,8 This is interpreted in terms of the model, as a demonstration that the recombination rate constant of some of the point defects considered needs to change with the current density, possibly due to changes in the structure of the film and a future work will specifically discuss this problem. Finally, the theoretical existence of a quasi-stationary situation for the point defect concentrations at very high rates of growth is shown, as must be expected for the case in which there is no time for normal aging.
In relation to the specific ionic resistivity it is shown that it passes through a minimum and it decreases with the increase of the current density, as must be theoretically expected and was previously found experimentally. Due to the fact that, in the present work, the variation of the recombination rate constant with the current density was not considered, for the case of the specific ionic resistivity, and opposed to prior experimental results, the minima shift to higher thicknesses for higher current densities.
Finally, with the mathematical description of the model it was possible to calculate the overpotentials through the film which qualitatively follow the expected theoretical results of the Ohmic model with no constant parameters (transient conditions).