Modelling diffusion controlled electro-deoxidation of metal oxide to metal in molten salt

Diffusion is a fundamental irreversible process intervening in the evolution of many out-of-equilibrium systems and is successfully described by Fick’s law obtained from non-equilibrium thermodynamics. Despite this, numerical simulations of solid state electro-deoxidation in the diffusion-controlled regime in molten salt remain elusive. Here, a new model for diffusion controlled three-phase interline (3PI) penetration in a porous cathode during electro-deoxidation is validated against experimental observation. This penetrating 3PI model is applied at high overpotential and benchmarked using the oxygen ionisation TiO2(s) + 4e− → Ti(s) + 2O2− at the 3PI. The model couples slow diffusive transport and fast oxygen ionisation while assuming a negligible ohmic potential drop in bulk molten CaCl2 electrolyte. The 14 nm s−1 penetration rate of the 3PI and the order of magnitude of 3PI currents (derived from an exchange current density and cathodic transfer coefficient of 0.32 A cm−2 and 0.01, respectively) in the chronoamperometric data for porous cathodes are in good agreement with experimental observation.


Introduction
Evolution of the FFC Cambridge process is necessary since kinetic barriers persist in solid-state electro-deoxidation, despite an acceptable level of energy consumption, including the low-current efficiency for Ti extraction relative to the extractions for Cr and Zr [1].Electro-deoxidation in molten salt can reduce a variety of metal compounds, particularly oxides, to their respective metals, alloys, and intermetallic compounds [2][3][4][5].Although thermodynamics is the best predictor of reaction and phase equilibria, kinetic effects diminish at high molten-salt temperatures [6,7] while the innate complexity of electro-deoxidation leads to factors that complicate thermodynamic calculations.Moreover, threephase interline (3PI) models derived at low overpotential [8,9] are unfavourable for fast oxygen ionisation rates [3] and assume that charge transfer kinetics at the cathode is dominant.Thus, a Tafel equation, which is a limiting case of the Butler-Volmer equation at large overpotentials, is needed to model diffusion-limited current in electro-deoxidation.Also, existing models of electro-deoxidation, based on a phase field model [10,11], strong solution theory [12], and the Nernst-Planck equation [13,14], have speciously modelled electroactive minority oxygen ion transport by migration instead of diffusion while assuming a significant ohmic potential drop in bulk molten salt electrolyte [15].Critically, although relevant to process scale-up, very little systematic exploration is available on enhancing the performance of electro-deoxidation cells based on rate-limiting oxygen ion transport inside a porous cathode [15].
In this work, a new model for diffusion-controlled 3PI penetration in a porous cathode during electro-deoxidation is validated against experimental observation.Unlike other cell designs, our penetrating 3PI model in a 2D axisymmetric electro-deoxidation cell permits the inclusion of alternative anode designs with low current density [16] while allowing unprecedented observation of 3PI movement.During electro-deoxidation, the model assumes that the potential drop in the bulk molten CaCl 2 electrolyte is negligible, since non-electroactive chloride ions provide a high conductivity.The unique formalism of combining non-equilibrium thermodynamics [17,18] with phenomenological 3PI theory [2,19], alongside fitting parameters (i.e.overlap function and exchange current density), empirical parameters and model assumptions, is calibrated against experimental data.The model simulates chronoamperometric data by boundary integration of the 3PI current.To this end, fitting using the minimisation of least squares reveals good agreement between numerical and experimental chronoamperometric data for porous cathodes with inferential support from a statistical chi-squared goodness-of-fit test.The lack of agreement for a dense cathode is explained by the fundamental difference between the electro-deoxidation in a porous cathode and a dense cathode, including the different mechanism for metal nucleation with diffusion-controlled growth.This work paves the way to developing experimentally validated models of electro-deoxidation in metal powder production.

Field distribution in a 3PI model
In the penetrating 3PI model, solid boundaries (except for the anode) and the symmetry axis in contact with bulk electrolyte as well as the electrolyte surface are assumed insulating to current flow and ion diffusion in a 2D axisymmetric electro-deoxidation cell domain [13,14].The domain comprises active (i.e.porous metal, solid metal, porous metal oxide and bulk electrolyte) and inactive (i.e.solid anode, solid current collector and solid metal oxide) subdomains in which a boundary is assumed planar since the cell domain size is much larger than the Debye length [18].Assuming a thin, porous metal layer at time t = 0 (Fig. 1a), the porosity of porous metal oxide ε o and the cathode thickness h th of the porous cathode is ε o ≥ 0.32 and h th = 2.6 mm whereas for the dense cathode, it is ε o = 0.22 and h th ≥ 3.3 mm.The Maxwell equation curl(E) = 0, where E = − ∇ φ is the elec- tric field, imposed on a boundary in the cathode implies that the electric potential φ at the boundary is spatially constant and independent of adjacent subdomains (i.e.φ| i − = φ| i + in Fig. 1a).The potential at the current collector/metal boundary and at the electrolyte/anode boundary is − 3.2 V and zero, respectively.
The ionic flux vector in the cathode is expressed by the Nernst-Planck equation for the superposition of diffusion, migration and convection: where c is the ion concentration, z is the charge number, F is the Faraday constant, u is the electrical mobility, D is the diffusion coefficient, u is the velocity vector and κ is the electrical conductivity.The diffusion-migration equation assumes steady-state conditions and a stagnant electrolyte in the pores [15]: where κ o and κ m is respectively the effective electrical conductivity of porous metal oxide and porous metal, and D o and D m is respectively the effective diffusion coefficient of oxygen ions in porous metal oxide and porous metal.The direct conversion of porous metal oxide to porous metal at t > 0 leads to a significant impact of the porosity of reductiongenerated metal ε m = 1 − [V r (1 − ε o )/(1 − S)] on the tortuosity in porous metal (i.e.Bruggeman exponent > 1.5) [20].Thus, Eqs. 2 and 3 are given by where V r = V m / V o where V m and V o is respectively the molar volume of metal and metal oxide, D e and κ e is respectively the diffusion coefficient of oxygen ions in and electrical conductivity of bulk molten salt electrolyte, and S1) where m o , d o , h o and ρ o are respectively the mass, diameter, thickness and absolute density of a porous metal oxide precursor.The temperature dependence of κ e can be described by an Arrhenius equation [21].Setting the highest (equilibrium) oxygen ion concentration c| i at the 3PI (Table 1) and the lowest oxygen ion concentration (c| an = 0) at the anode ensures electro-neutrality in the bulk electrolyte.Non-electroactive majority charge carriers eliminate the electric field in the bulk electrolyte so that electroactive minority oxygen ion charge carriers undergo Fickian diffusion: It is important to note that electron transfer across the bulk electrolyte is of a diffusional nature and does not provide an ohmic potential drop [15].

Electron transfer at a 3PI
By applying multiple-scale formal asymptotic analysis [22], u = u λ separates domain space into electronically conducting metal M(s), insulating metal oxide MO x (s) and ionicconducting electrolyte phases.This enables the oxygen ionisation MO x (s) + 2xe − → M(s) + xO 2− at a discrete M(s)/ MO x (s)/electrolyte 3PI (Fig. 1a) where electron transfer at the M(s)/MO x (s) interface results in the conversion between two solid phases, accompanied by oxygen ion transfer at the MO x (s)/electrolyte interface to maintain electro-neutrality in the solid phases [19].Substituting a series expansion for u λ far from the discrete 3PI and near  Oxygen ion concentration at the 3PI 4200 mol m −3 [20] Ginzburg-Landau equations gives Eq. 7 for λ ≪ 1 and Eq. 8 for the asymptotic limit λ → 0, respectively [22]: where W(u 0 ) is a non-negative potential metric with minima at the M(s), MO x (s) or electrolyte phase and 1a) where λ represents the interface thickness between phases, r represents a position in space, and τ = λ2 t where t is time.Interface generation occurs as the diffusion term (Eq.7) becomes large indicating that u 0 ** (r, ∞ ) tends to M(s), MO x (s) or electrolyte phase.The inte- ger number m i of discrete 3PIs within a porous cathode is equivalent to a 3PI having a continuous reaction zone of length (1/m i )h th (Fig. 1a) [19] where h th is the thickness of the porous cathode.
The 3PI is a thermodynamic boundary in local equilibrium in which the Nernst equation for the oxygen ionisation reaction determines its thermodynamic feasibility.Thus, the minimum potential required to ionise oxygen is derived from the standard Gibbs free energy change [2].The entropy production [17] in the 3PI is positive along a reaction γ-coordinate according to the second law of thermodynamics: where μ(γ) = RTIn[f(γ)Γ(γ)] is the effective chemical potential where Γ(γ) is the excess oxygen ion adsorption at the 3PI, R is the universal gas constant, and f(γ) = exp{[C(γ) + Fφ(γ)]/RT} is the activity coefficient where C(γ) is an intrinsic energy barrier and T is the temperature of the electrolyte, and T i = T is an isothermal approximation where T i is the temperature of the 3PI.The normal current density for electron transfer j e n| i across the 3PI is from the porous M(s) subdomain Ω − to the porous MO x (s) subdomain Ω + (i.e.quasi-1D approximation n = ± 1) (Fig. 1a).Based on Eq. 9, the thermodynamic force [17] conjugate to the flux j e n| i is given by The linear flux-force relation with an assumed constant j e n| i along the γ-coordinate includes an integral of inverse mobility weighted with a Boltzmann factor: where L(γ) is a phenomenological coefficient, μ eq is the effective chemical potential at equilibrium, b = L(γ)/FΓ(γ) is the transport coefficient for electric mobility, and μ(1) and μ(0) is respectively the effective chemical potential for the state of the products (M(s) and xO 2− ) and reactants (MO x (s) and 2xe − ).Since a large overpotential η i at the 3PI is favourable for fast oxygen ionisation rates, introducing μ(0) − μ eq → 0 and μ(1) − μ eq = − αFφ| i + (assuming φ| i + ≡ η i ) into Eq.12 gives a Tafel equation: where α is the cathodic transfer coefficient for the electric potential that alters the activation energy barrier.By applying Marcus theory, charge transfer alters nuclear configurations from M-O bonds in MO x (s) to M-M bonds in M(s) such that the reorganization energy necessary for this transformation is comparable to polarization (i.e.α < 0.5) [19].The exchange area current density through the porous cathode cross section is where j 0i is the exchange line current density with reference to the 3PI [19].The j 0 is an intrinsic property of cathode phases that encompasses structure and composition effects at the 3PI.

Diffusion-controlled penetration of a 3PI
During 3PI penetration into the cathode interior, removal of oxygen ions from the metal oxide surface occurs while metal adatoms tend to form metal nuclei rather than a continuous layer thereby forming physical gaps between metal oxide and metal.Newly formed metal adatoms diffuse along solid surfaces (interfaces) to reach preformed M(s) nuclei [23].The continuous formation of these nuclei after the application of an electric potential [14,24] is given by where J m and t m is respectively the peak current and time taken to reach peak current, and J is the current.The growth of M(s) nuclei and the interference of their diffusion zones [14] gradually results in the development of a planar ( 12) diffusion current at the 3PI.Current density for multiple nucleation and diffusion-controlled growth [25] at the 3PI is given by where J(t, u) is the current to a single metal nucleus born at time u, and Θ(t) = θ(t)/θ * (t) is the overlap function describing the effect of the interference of diffusion zones where θ * (t) is the non-overlapped fraction of the 3PI area covered by planar diffusion zones associated with M(s) nuclei (Fig. S2) and θ(t) is the overlapped fraction of the 3PI area covered by planar diffusion zones [14].Preferential M(s) growth sites close to a growing M(s) nucleus can be interrupted by oxygen as it prefers to adsorb onto these sites [23].Thus, the appearance rate of M(s) nuclei is where N 0 is the growth site density, N m is the nucleus number density, N ad is the number density of adsorbed sites and A is the nucleation rate constant.Charge transfer across the 3PI is proportional to the amount of porous metal formed (or porous metal oxide consumed) according to Faraday's law: where M o is molar mass of MO x (s) and n is the number of moles of electrons transferred.The circumferential inward-3PI penetration nξ(t)| i π from the current collector is assumed to apply over the entire 3PI.Differentiating both sides of Eq. 18 with respect to t, applying chain rule differentiation and rearranging give Equating electron transfer flux j * (t)/nF to Fickian oxygen ion diffusive flux − ε m 2 D e (∂c/∂l) in Eq.21 gives the electrodeoxidation equation for 3PI penetration: where l d is a linear distance parameter between the current collector and anode boundary.Since the oxygen atom is usually larger than the metal atom, oxygen ion diffusion out of the porous cathode leads to a solid volume decrease of the cathode [26].However, assuming that h th and the cathode diameter are constant during electro-deoxidation, the volume shrinkage of the cathode is assumed negligible (i.e. S = 0).Hence, Eq. 22 simplifies to Although oxygen diffusion occurs in solid cathodes [27], performance enhancement of electro-deoxidation cells is based on rate-limiting oxygen ion diffusion through the pores of cathodes [15].

Experimental validation of a 3PI model
The preparation of reduced cathode cross sections and in situ chronoamperometric measurement (see Osarinmwian et al. [13] for details) for the experimental validation of the penetrating 3PI model is briefly described.Cylindrical porous TiO 2 cathode precursors were prepared and then characterised using Helium gas pycnometry.The cathode system was assembled by threading a precursor onto a stainless steel current collector rod through a central 5-mm hole while the anode was assembled by attaching a cylindrical graphite annulus to a specially designed stainless steel current collector.Thus, the 2D asymmetry in the electro-deoxidation cell was formed by enclosing the cathode system within the anode at a 4-cm immersion depth in molten CaCl 2 electrolyte.A controlled potential ramp up to − 3.2 V was applied to initiate electro-deoxidation while eliminating excessive J m during in situ chronoamperometric measurement.Retrieval of the reduced cathode from the cell after electro-deoxidation was followed by sectioning of the cathode into a segment, immersing it in epoxy resin and then polishing it before ex situ cross-sectional analysis using optical microscopy (Fig. S3).
A previously reported electron-transfer-controlled 3PI model failed to show complete electro-deoxidation due to severe numerical instability during 3PI penetration (Fig. S4).Also, existing simulations of boundary-integrated 3PI current in a penetrating 3PI model were an order of magnitude higher than J m while those in a thin layer 3PI model showed good agreement [13].However, the latter was based on erroneous model assumptions and an inapt application of kinetic parameters (i.e.j 0 = 289 A m −2 and α = 1.99) [28] to the oxygen ionisation TiO 2 (s) + 4e − → Ti(s) + 2O 2− .In this work, the analysis of experimental data with respect to the penetrating 3PI model requires fitting Eqs.23 and 13 to experimental data to obtain the parameters Θ(t) and j 0 .To obtain a best fit using the minimization of least squares [29], the model function ŷ(x; p) of the variable x = Θ(t) or j 0 and the vector of the number of model parameters p to a set of m data points (x z , y z ) is fitted to minimize the sum of the weighted squares of the errors between the experimental data y z and the model function: where σ yz is the measurement error for y(x z ).Since ŷ(Θ(t); p) is non-linear, minimization of χ 2 (p) with respect to p is performed iteratively by changing Θ(t).An initial guess of Θ(t) = 2 × 10 −3 in Eq. 23 results in dξ(t)| i /dt = 4 nm s −1 and visual observation of partial electro-deoxidation at t = 24 h (Fig. S5).By carefully increasing the value of Θ(t), the guess Θ(t) = 7.163 × 10 −3 results in dξ(t)| i /dt = 14 nm s −1 and visual observation of complete electro-deoxidation at t = 24 h.The fitting of ŷ(j 0 ; p) to y(j 0 ) is achieved by adjusting j 0 such that the boundary-integrated 3PI current at t = 2 h equates to the arithmetic mean in experimental current at t = 2 h.Numerical simulation at t ≥ 2 h is assessed inferentially against y(j 0 ) using the chi-squared χ 2 goodness-of-fit test.Given that most statistical tests require normality assumptions, goodness-of-fit tests for normality are useful procedures in determining whether a random sample of size m originates from a normal (Gaussian) population [30].The χ 2 test is closely related to least square fits with the exception that the hypothesis is fixed.In this work, the χ 2 test is designed to find whether the observed set O z of mean experimental current is in agreement with the expected distribution E z of boundary-integrated 3PI current.This is achieved using the test statistic: The test rejects the null hypothesis (i.e.no significant difference between numerical and experimental chronoamperometric data) at an α significance level if χ 2 is larger than the (1 − α) 100% quantile of the χ 2 distribution.In other words, χ 2 > χ 2 (m − 1) where χ 2 (m − 1) is the critical chi-squared test statistic and m − 1 is the degrees of freedom.

Results and discussion
Figure 2a shows a decrease in 3PI length during the penetration of porous metal into porous metal oxide in a porous cathode (ε o = 0.32, h th = 2.6 mm and ρ o = 4274 kg m − 3 ).The order of magnitude of j 0 = 0.32 A cm − 2 from TiO 2 deoxidation is in good agreement with recently reported values of j 0 [9] and is assumed constant during 3PI penetration in 2 varies linearly as t → ∞ due to the increasing effects from concentration and ohmic polarizations during electro-deoxidation (Fig. 2b) [8].For metal oxides with relatively small V r (Table S1), electrodeoxidation forms a highly porous metal that improves oxygen ion diffusion pathways leading to higher dξ(t)| i /dt (Fig. 2b).The highest dξ(t)| i /dt = 26 nm s −1 for TaO 2.5 deoxidation, which possesses the lowest V r and highest M o , is in good agreement with a theoretical model [20] whereas the lowest dξ(t)| i /dt ≤ 14 nm s −1 for TiO 2 and ZrO 2 deoxidation is consistent with the interruption of metal growth by oxygen adsorption affecting metals with a large solubility for oxygen [23].In TiO 2 deoxidation, experimental observation of slow inward 3PI penetration from the current collector (Fig. S3) is consistent with lower φ (Fig. 1b) and higher c (Fig. 1c) near the current collector.
The numerical chronoamperometric plot for TiO 2 deoxidation (Fig. 3a, b) shows that the diffusion current to an area of 3PI corresponding to planar diffusion zones, taking into account their overlap, slowly declines at t ≥ 2 h consistent with the decrease in 3PI length (Fig. S3).Similar to an overlap model [25], the Θ(t) < 1 describes the effect of overlapping diffusion zones on the current to a single metal nucleus at time u.The mechanism of progressive nucleation with diffusioncontrolled growth in a porous cathode is slow (i.e.A → 0 and N 0 → ∞ ) and involves a non-uniform size distribution of dif- fusion zones.Moreover, the good agreement between numerical and experimental chronoamperometric plots for porous cathodes (Fig. 3a) may contribute to addressing the challenge of correlating time, current and cathode composition caused by a slow CaTi 2 O 4 decomposition rate and strong electronic background current [31].However, the statistical inference of a good agreement for a dense cathode (Table S2) is contrary to visual observation (Fig. 3b).The violation of the normality assumption underpinning the χ 2 test may contribute to this invalid inference [30].Furthermore, the fundamental difference between the electro-deoxidation in a porous cathode and a dense cathode explains the significant deviation between the experimental and numerical chronoamperometric plots for the dense cathode.
The porous cathode initially undergoes electrochemical reactions in which there is a reduction of TiO 2 to intermediate Ti oxides of decreasing valency and CaTiO 3 throughout the entire cathode, culminating in Ti 2 O 3 formation with little oxygen removal [31] during the formation of a thin Ti layer at t < 2 h (Fig. S3).As the Ti layer thickens, the sharp decrease in the rate of change in current is caused by complete reconstruction of the Ti 2 O 3 lattice [31].The reaction pathway in a porous cathode is more thermodynamically favourable and kinetically facile than in a dense cathode [32].This is demonstrated by a sharp decline in current beyond J m [13] and the presence of current waves during potential ramping [14].In a dense cathode, the alternative reaction pathway is more direct and characterised by an inward moving 3PI since the blocking of cathode porosity, by initial Ca-incorporation reactions, restricts electrodeoxidation to the cathode surface prior to 3PI penetration into the cathode interior [32].Thus, unlike porous cathodes [14], in a dense cathode, the nucleation mechanism cannot be described as progressive nucleation (Fig. 3c).Also, since cathode geometry is related to diffusion distance and the available surface area in contact with the electrolyte [33], the greater available surface area in the porous cathode leads to higher dξ(t)| i /dt (Fig. S3).
Overcoming the perception that a random packed bed of particulate feedstock cannot be reduced due to insufficient electrolyte flow [34] is aided by applying electro-deoxidation modelling in metal powder production.The ability to control dξ(t)| i /dt by changing h th (Fig. S3) is consistent with optimal electro-deoxidation of a Fig. 2 Diffusion-controlled 3PI penetration.Arbitrary Lagrangian-Eulerian method simulates 3PI penetration using the finite element method thereby predicting the location of the 3PI.Penetration via Eq.23 circumvents inverted finite elements by maintaining positive finite element volume (Fig. S4).Physicochemical properties of the metal oxides are shown in Table S1.a The dξ(t)| i /dt for TiO 2 deoxidation is fitted to experimental observation by adjusting Θ(t), which results in good qualitative agreement with experimental observation in terms of dξ(t)| i /dt (Fig. S3).Thus, the fitted Θ(t) is used to simulate 3PI penetration for all metal oxides.The negligible effect of j 0 on the shape and rate of penetration (Fig. S5) indicates that the penetration is diffusion controlled for all metal oxides.b [ξ(t)| i ] 2 transient and dξ(t)| i /dt for different metal oxides random TiO 2 particle packing at low particle size and packing depth [35].Moreover, increasing the volume fraction of a conductive additive in a composite cathode, comprising metal oxide and conductive additive particles, increases both the electrical conductivity κ x = κ A [(f A − f c )/ (1 − f c )] ω [36], where κ A is the electrical conductivity of the conductive additive, f c is the percolation threshold, and ω is the critical exponent of the composite cathode, and volume-specific 3PI length (Fig. S6) in the cathode.These increases are consistent with the increase in the electro-deoxidation rate of Si-SiO 2 composite cathodes with increasing Si content [37].Investigating the effect of volume-specific 3PI length in composite cathodes on dξ(t)| i /dt is important because the size and location of metal oxide particles in a cathode bed strongly affect the internal current distribution [35].Finally, combining an electro-deoxidation cell design comprising cathode stacks and horizontal anodes [13] with the penetrating 3PI model provides a first step to applying electrodeoxidation modelling to the development of bipolar electro-deoxidation cells.

Conclusions
Non-equilibrium thermodynamics provides a physically consistent and structured theory for guiding model development in electro-deoxidation.The local equilibrium hypothesis provides a theoretical framework for irreversible processes during electro-deoxidation over long time scales while mesoscopic phenomena in the 3PI account for shorter time scales since additional degrees of freedom that have not yet equilibrated exert an influence on 3PI penetration.This systematic framework couples slow diffusive transport and fast electron transfer at the 3PI while assuming a negligible ohmic potential drop in the bulk electrolyte.This work provides an important milestone in the development of experimentally validated models of electro-deoxidation while providing a much needed, but missing to date, design tool for aiding process scale-up and commercialisation.

Table 1
Constants from Eqs. 13 and 23 that are used as input data in the penetrating 3PI model e Diffusion coefficient of oxygen ions in molten CaCl 2 2.5 × 10 −9 m 2 s −1 [20] R Universal gas constant 8.314 J mol −1 K −1 l d Linear distance between current collector and anode 0.0255 m c| i