Calculation of the space charge distribution in poled soda-lime glass

An analytical model of electric field assisted diffusion of ions into a multi-component glass is extended to calculate the space charge that forms between the poled layer and the potassium peak in a poled soda-lime glass. The model is compared with numerical solutions to the drift-diffusion equations and Poisson’s equation and shows good agreement. Some recent experimental results in corona poled soda-lime glass are also discussed using this model.


Introduction
A recent paper presented an analytical model of field assisted diffusion in a multicomponent glass containing mobile Na + and K + ions [1]. The model is applicable to the analysis of poling experiments when a non-blocking anode is used which results in H + /H 3 O + ion injection into the glass and predicts the pile-up of K + ion below the poled glass layer. This model has recently been compared in detail to concentration profiles obtained by energy dispersive x-ray analysis measurements of corona poled, soda-lime glass slides [2] and with computer simulations based on drift-diffusion models. The authors observed a giant pile-up of K + ions below the poled layer, which initially had a very low K + concentration in the original glass, and noted that this was in good correspondence to the analytical model [1,2]. They also visualized a large space charge region between the H + poled layer and the K + pile-up layer using scanning electron microscope (SEM) measurements [2]. This was also observed in their simulations as was a much smaller space charge at the leading * Author to whom any correspondence should be addressed.
Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. edge of the K + pile-up layer, next to the un-poled sodium rich substrate [2].
Although concentration profiles are well explained, one criticism of the analytical model was that since it adopted a space charge neutral approximation, it was unable to model the space charge that exists between the poled glass and the trailing edge of the K + ion pile-up layer [2]. As pointed out space charge effects are important in explaining certain phenomena including second order non-linearity effects, Pockel effects and nano-printing applications that have been observed in various poled glasses [3][4][5].
The purpose of this paper is to show that using the analytical results of [1] coupled with previous analytical results [6], it is possible to include space charge effects accurately in an analytical model. This has utility in that it reveals the parameters that the space charge distribution is theoretically dependent upon in a typical poling experiment, which up to now has been lacking. The improved analytical model is shown to be in good agreement with space charge profiles obtained from computer simulations based on drift-diffusion equations and Poisson's equation with parameters relevant to poled soda-lime glass, and with simulated space charge profiles reported in [2] at both the trailing and leading edges of the K + ion pile-up region. This addresses the criticism of the analytical model [2].
The authors in [2] also observed that the diffusional smearing at the H + -K + interface obtained experimental profile was larger than expected and offered a suggestion why this might be so. We also discuss this aspect of their results and compare previous optical reflectivity data on poled samples in the same glass in the light of the space charge model presented here.

Concentration and space charge profile model
We first present the analytical results from the space charge neutral model since these will used in the analysis. On the assumption of a voltage applied to the sample during poling at an elevated temperature, the ion flux j(t) through the sample decreases with time due to the increasing glass resistance associated mainly with the hydrogen ion drifting and diffusion into the glass at the anode. The hydrogen ion concentration profile is given approximately by [1] (1) once a so called quasi-stationary state has been established. This assumes an unlimited hydrogen ion supply at the surface so C H (0, t) = C 0 . The corresponding Na + ion concentration profile is given by In these equations j(t) is the ion flux(=J(t)/e where J(t) is the current density), C NaO and C KO are the uniform Na + and K + ion concentrations in the original glass, M K is the K + to Na + mobility ratio, M H is the H + to Na + mobility ratio and C 0 = C NaO + C KO is the total mobile ion concentration in the glass. An analytical expression for j(t) is also given in [1] but is not reproduced here. The results were obtained by extending the so called quasi-stationary approach, which was first applied to ion exchange into a glass with a single ion species [7], to a multi-component glass.
The K + ion profile C K (x, t) in the model was calculated on the assumption of space charge neutrality hence Previously, the author has also analysed the drift and diffusion of a single ion species drift-diffusing into a glass containing just a single mobile ion species [6]. This analysis however, was based on the solution of the drift-diffusion equation on the assumption of a constant ion injection rate, j, rather than a defined voltage boundary condition but did include space charge effects [6]. The relevant results of that analysis was that the normalised invasive ion concentration profile C * = (C/C 0 ) is given by and In these equations C 0 is the concentration of ions in the glass, D A and μ A are the diffusion coefficient and mobility of the invasive ions and M is the invasive ion to sodium ion mobility ratio.
Unfortunately, an exact closed form solution of equation (4) could not be found. However, the parameter B is normally 1. For example, it was shown that B ∼ 10 −8 for field assisted Ag + -Na + ion exchange and B ∼ 10 −4 for K + -Na + ion exchange in soda-lime glass [6]. So for these cases an approximation was made by putting C * = g in the right-hand side of equation (4) to give the first order approximation and the authors concluded that space charge effects were insignificant for these applications and C * = g, which is the conventional space charge neutral solution for a single ion glass with constant current boundary condition [8]. Now, the analytical model in [1] shows that once K + ion pile-up occurs, the H + ions follow the quasi-stationary concentration profile equation (1) and so effectively are being driven into a glass region where all sodium ions are replaced by K + ions. Further from equation (1), the slope of the leading edge of the H + ion profile is determined by the instantaneous ion flux j(t) at the end of the poling process and the time integral of the ion flux simply determines the depth of the poled glass region. This implies that equations (4)-(6) can also be used to determine the H + concentration profile with space charge effects included at the H + -K + interface provided M in equation (6) is replaced by M H /M K (=D H /D K ) and In order to determine the space charge profile, we consider the expression for the electric field, E, from [6] where we change notation from μ A to μ H for this particular application Substituting the spatial derivative dE/dx from equation (10) into Poisson's equation and using dη/dx = j(t)/C O D H , from equation (9) gives the space charge ρ which is the last term on the right-hand side of equation (4). So for example, if we sub (8) into (11) it gives a first order formula for the negative space charge distribution This may also be deduced from (4) directly: since C * = g is the electrical neutral solution then any deviation from that implies the development of a space charge.
Using the experimental data from [2] for corona poling, it can be estimated that the average poling current density is J ∼ 0.8-0.9 A m −2 for their H + ion injection and J ∼ 0.5 A m −2 at the end of the poling process, so using μ Na = 1.2 × 10 −15 m 2 V −1 s −1 , then for 300 • C, B ∼ 0.12. So for poling experiments the parameter B is much larger than that for conventional ion exchanges considered in [6] and it is no longer clear that a first order approximation, equation (12) is valid. We have thus extended the approximate solution of equation (4) to include higher order terms to extend its range of validity. The analytical details are presented in the appendix and we present the final result here for the hydrogen concentration profile and the normalised space charge and so In figure 1 we plot this space charge function for a number of values for B < 1. We compare the accuracy of this function by including only the linear term in B, including the B 2 term and including the B 3 term. It can be seen that including only the linear term in B (=equation (12)), is adequate for B < 0.1 whilst second or third order terms in B are needed for B ∼ 0.2. Higher order terms, an alternative approximate solution or a numerical integration of equation (4) would be necessary for larger B. A feature of equation (14) is that the space charge is a function of g(1 − g), which is an even function of η. From the first order approximation, equation (12), it can be seen that the peak magnitude of the space charge distribution is given when g = 0.5 and has a value-B(1 − M) 2 /4. For B = 0.12 and M = M H /M K = 1 × 10 −3 /0.0125 = 0.08, then the space charge at its peak magnitude is ρ/(eC 0 ) = −0.024 (based on including third order terms in equation (14)) is in excellent agreement with the computer results from figure 2 of [2]. Figure 2 shows the two concentration profiles C H * , C K * using the space charge neutral approximation and those including space charge with terms including B 3 in equation (13) for B = 0.2. It can be seen that for the concentration profile there is minimal difference as observed if space charge effects are included. This conclusion was reached previously by numerical computation [1,2]. When plotted as a function of the normalised variable η it can be seen that the space charge width and concentration distribution transitions appear wide. However, using the parameters in [2] and equation (9) a change of η of 10 corresponds to approximately 1.5 nm.

Comparison with numerical computations
We compare the analytical expressions derived in the previous section with numerical solutions of the drift-diffusion equations and Poisson's equations in order to confirm accuracy of the space charge function. Following Petrov et al [9] we numerically solve We solve these numerically by using the MATLAB PDE toolbox with model parameters C O = 7.27 × 10 27 m −3 and C NaO /C O = 0.9477 and C KO /C O = 0.0523, which are appropriate for Menzel soda-lime glass [2]. We also assume a processing temperature of 300 • C, μ Na = 1.3 × 10 −15 m 2 V −1 s −1 , M K = 0.0125, M H = 0.001 and a Havens ratio of 1. Figure 3 shows concentration profiles calculated from simulations corresponding to a poling time of 0.25 h. Excellent agreement between simulations and the analytical model is demonstrated [1,2]. For this poling time the pile-up region of K + ions and hence the quasi-stationary profiles for Na + and H + ions are well established.
The space charge between the H + and K + rich regions is shown in detail in figure 4 for the same simulation parameters. It can be seen that the analytical model of the space charge region is in good agreement with that obtained from the computer simulations both in width and depth. The value of B for this example is 0.084. Although the space charge width is small, the total charge per unit area is −0.013 Cm −2 for this example due to the large value of C 0 . This value was obtained by integration of ρ but can also be determined by applying Gauss's law and using the difference in electric fields from either side of the space charge in the H + and K + rich regions where the concentrations are constants. The large value of C 0 in soda-lime and borosilicate glasses explains why concentration profiles from the quasi neutral approximation are  in good agreement with calculations (analytic and simulations) that include space charge effects. In such a case a sizable charge can develop for only a small deviation in the concentration profile.
The model can also be used to predict the space charge profile at the leading edge of the K + pile-up region. Here K + ions are being driven into the sodium rich glass bulk so the relevant value for B, which we call B 2 , will be much smaller since μ A and D A in equation (5)    order approximation in equation (14) is more than adequate in calculating the magnitude of the space charge peak and this gives ρ 2 /(eC 0 ) = −B 2 (1 − M K ) 2 /4 = −1.87 × 10 −4 which is close to the value observed in previous computer simulations [2].
The analytical approximations, equations (13) and (14), for the H + /K + interface presented here would not be usable for glasses with very low mobile ion concentrations, for example for poling silica glass where C 0 ∼ 10 22 m −3 , since for these B > 1. In glasses with mobile ion concentrations of this magnitude a complete depletion region is formed with ρ = −eC 0 [9,10].
Finally we note that SEM measurements have shown a transition between the poled layer and the K + peak that was wider (∼100 nm) than the analytical model or computer simulations predicted [2]. Although the potential for convolving effects due to beam size and electron scatter in SEM measurements were acknowledged, it was also suggested that since the space charge represents a concentrated region of unoccupied non-bridging oxygen ion sites these could result in enhanced inter-diffusion at the H + -K + interface due to a lower activation energy. Hence a wider transition in the concentration profiles and a wider space charge region could result due to this effect. However, recent optical reflectivity measurements of leaky modes supported by poled layers in the same glass indicated the need for a more rapid change than 100 nm between the low refractive index of the poled glass and the higher refractive index of the glass in order to obtain the observed sharp reflectivity peaks [11]. Since the model in [1] and the profiles in [2] demonstrate a significant pile-up of K + ions (close to a value C 0 ), a re-interpretation of the reflectivity measurements suggests that the rapid transition in refractive index responsible for the reflectivity spectra corresponds to the H + -K + transition in the concentration profile. However, the values for current density were at least a factor of 2 smaller and B from equation (5) is about a factor of 4 smaller in these experiments compared to those in [2] and thus figure 1 suggests a less significant space charge in these samples. So the optical reflectivity data is not in contradiction with the enhanced diffusion coefficient suggestion [2].

Conclusion
The analytical model of poling has been extended to calculate the space charge between the poled layer and the trailing edge of the K + pile-up region in a soda-lime glass. The analytical model presented shows good agreement when compared to those obtained from computer simulations based on the full drift-diffusion equation coupled with Poisson's equation. A symmetrical space charge region between the H + and K + regions has been predicted. Parameters that the space charge is dependent upon have been deduced.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.