Models and Their Limitations in the Voltammmetric Parameterization of the Six‐Electron Surface‐Confined Reduction of [PMo
 12
 O
 40
 ]
 3−
 at Glassy Carbon and Boron‐Doped Diamond Electrodes

In acidic media, three well-defined surface confined overall two electron transfer reactions derived from closely spaced oneelectron steps are observed in the DC and AC voltammetry of 1⁄2PMo12O40� 3 . In this study, parameterization of these processes by sinusoidal AC voltammetry at glassy carbon (GC) and boron doped diamond (BDD) electrodes has been undertaken in 1 M HClO4. The experiment-simulation comparison of even a simplified reaction scheme required estimation of 21 parameters, which is an exceptionally large number in voltammetry. The six reversible potentials were robustly recovered from data obtained at low (9.02 Hz) and higher (60.05 Hz) frequencies at both GC and BDD electrodes. Examination of the harmonic content by the Fourier transformed version of AC voltammetry enabled imperfections in the model associated with faradaic and non-faradaic aspects of the voltammetry to be readily detected. After consideration of these imperfections, all processes were considered to be reversible within experimental uncertainty. Subtle differences in reversible potentials observed at both electrodes account for electrode dependence in shapes of the overall two electron steps. Fits of experimental data were superior using BDD at both frequencies and for the higher frequency GC case. Non-compliance with an assumed potential independent model of the double layer capacitance was considered to be a significant contributor to discrepancies found with the GC electrode.


Introduction
Voltammetric theory for electrode processes comprising a series of heterogeneous electron transfer steps and coupled chemical reactions is now very well established (see [1][2][3][4][5] for example). Generation of voltammetric data derived from a model selected to mimic the experimental response is known as the forward problem. However, for a complex electrochemical reaction, determining the large number of parameters that have to be deduced by comparison of experimental and theoretical data, in what is termed the inverse problem, often remains unmanageable with respect to achieving a complete and unique solution and in deciding when acceptable agreement between experimental and simulated data has been achieved. In essence, solving the inverse problem requires capturing substantial amounts of very high quality experimental voltammetric data and repetitively comparing with theoretical data deduced from a simulated model until satisfactory agreement is achieved.
In recent work, the Oxford Computing and Monash University Electrochemistry Groups have been developing protocols that reliably and efficiently parameterize complex electrochemical processes (see [6][7][8][9] for example). In summary, very large data sets containing current, potential, and time domain information are collected [6][7][8][9] at high resolution using instrumentation having 18 bit DAC and ADC converters. [10] Experiment-theory comparisons are then undertaken at levels ranging from fully heuristic to highly automated data optimization multi-parameter fitting exercises in attempts to define uniquely the thermodynamic, kinetic, capacitance and resistance and other parameters that are included in the model. In one of the more complex reaction schemes we have quantified, [9] values of thermodynamic, electrode kinetic and other variables present in a model were elucidated for a series of proton-coupled two-electron transfer processes associated with the reduction of the polyoxometalate ½PMo 12 O 40 � 3À (structure included in Figure 1a) adsorbed onto a glassy carbon (GC) electrode in 1 M H 2 SO 4 . In this particular example, integration of experimenter-controlled heuristic and automated data optimization was employed to estimate up to 17 of the parameters employed in the model. ½PMo 12 O 40 � 3À spontaneously adsorbs onto many electrode surfaces and has been used widely in this state in many practical applications (see for example [9,[11][12][13][14][15][16][17] ). In 1 M HClO 4 , the electrolyte of interest in this study, three well-defined surface confined reduction processes are found under conditions of DC cyclic voltammetry at glassy carbon (GC) and boron doped diamond (BDD) electrodes (Figures 1b and 1c respectively) within the potential range employed in these experiments. Under AC conditions (Figures 1d and 1e) three well defined responses are also detected at both types of carbon electrode. However, at the GC electrode (Figure 1d), processes I and II are sharper and have larger peak currents than process III, whereas at BDD (Figure 1e), the middle process II is sharper than either process I or III and has the largest peak current. Each process, despite commonly having a significantly different shape, corresponds to the overall two electron-two proton processes represented by Eq. 1-3. [9,11,18] Clearly, the background non-faradaic current also contributes more extensively to the total current (faradaic plus non-faradaic) in the DC and AC voltammetry derived from the GC electrode than for the BDD one (compare Figures 1b and 1d and 1c and 1e), and this also is expected to have an impact in the degree of difficulty and indeed level of success encountered in the the modelling.  Scanning to more negative potentials than used in Figure 1 reveals a fourth process IV which leads to rapid dissolution of solid from the electrode surface.
Each process summarized in Eq. 1-3, consists of two unresolved one electron transfer steps as in Eq. 4-5, where for convenience charges on species are neglected and redox levels are written in their non-protonated forms.  Figure 1. In terms of generating a model to describe this voltammetry, each electron transfer process can be assigned its own reversible potential (E rev ) which is a thermodynamic parameter derived by combining acid-base or ion-pairing equilibrium constants with the formal reversible potential (E 0 ). If the model assumes all six electron transfer processes are quasi-reversible, as is fundamentally the case, then modelling the electrode kinetics requires the inclusion of six heterogeneous electron transfer rate constants (k 0 values at E rev ) along with six charge transfer coefficients (α values) if the Butler-Volmer relationship is used. Additionally, the double layer capacitance (C dl ), uncompensated resistance (R u ), and the surface coverage (Γ) need to be modelled assuming they are not independently known. Under these circumstances it therefore follows that the simulation of the voltammetry shown in Figure 1 requires the inclusion of at least 21 parameters.
While simulation of the forward multi-step problem relevant to this study is inherently complex, it can be readily achieved nowadays [5] with the aid of commercially available DigiElch or KISSA or web downloadable software packages such as MECSim. Alternatively, in-house coding of the algorithms needed to generate the model can be undertaken as in our recent work on the reduction of ½PMo 12 O 40 � 3À at a GC electrode in 1 M H 2 SO 4 . [9] The major challenge now is to solve the computationally intense inverse problem in a manner that provides a reliable and physically meaningful understanding of the values and significance of a large number of parameters. In reality, none of the models in this or any other quantitative electrochemical study can be expected to mimic exactly the experimental data because they incorporate inexact theory or are devoid of one or more unknown reaction steps that are not accommodated in the simulation model. That is, there is an uncertainty associated with the underlying model itself as well as in parameter values estimated by data optimization. As a consequence, in a complex parameterization exercise it is wise to see if any systematic departure of experimental data from that predicted by the model are present, understand the origin and if possible revise the model. Analyses of AC voltammetry in its total current and Fourier transformed versions (FTACV) is ideal for both multi-parameter estimation and for detection of systematic deviation from the model used, as will be demonstrated in this study.
In this extension of studies at a GC electrode in 1 M H 2 SO 4 , we have compared simulated and experimental AC data for the reduction of ½PMo 12 O 40 � 3À at both GC and BDD electrodes in 1 M HClO 4 . Differences in density of states (DOS values), [19] adsorption capacity and other characteristics of the two carbon based surfaces [20,21] could impact on the values of parameters such as E 0 , k 0 , R u , C dl and Γ. Data analysis was initially conducted within the framework of knowledge gained and reported [9] when applying the mechanistically sensitive AC voltammetric technique to the analysis of the reduction of ½PMo 12 O 40 � 3À at a GC electrode with 1 M H 2 SO 4 as the electrolyte, but now with simultaneous parameterization of all 21 unknown variables from the total AC current data obtained at 9.02 and 60.05 Hz. In the 1 M HClO 4 environment, non-conformance of experimental data obtained at both GC and BDD electrodes to that from simulated models was searched for by examining the individual aperiodic DC and AC harmonics as resolved by Fourier transformation methods in what is termed FTACV. After assessing the impact of non-ideality, it was confirmed that all k 0 values at BDD remain extremely fast as is the case at GC. Interestingly, in studies where both the oxidized and reduced species undergoing electron transfer are soluble in the electrolyte medium used, k 0 values at BDD electrodes have been reported to be orders of magnitude slower than at GC. [19,22] However, E rev values, when deconvoluted from the overall two-electron processes, differ subtly and differences also were found in parameters such as R u C dl , and Γ which describe the physical chemistry of the GC and BDD electrodes. Voltammetric behavior compared to that reported [9] at a GC electrode in 1 M H 2 SO 4 is also considered. Data are also compared with those deduced from a quantitative square wave study at BDD [11,18] and GC [18] electrodes in 1 M HClO 4 , where the modelling was significantly simplified by including a number of variables parameterized in our study as knowns rather than unknowns. In reference, [11] it was assumed that all processes are reversible, whereas in reference, [18] quasi-reversibility was introduced into the modelling.

Chemicals and Reagents
Perchloric acid (HClO 4 ) and phosphomolybdic acid (H 3 PMo 12 O 40 , 99.99 %) were used as received from Sigma-Aldrich. Water employed for preparation of all aqueous solutions was obtained from a Milli-Q water purification system.

Instrumentation and Procedures
AC voltammetric measurements were undertaken using in-house instrumentation [10] with an applied sine wave having an amplitude (ΔE) of 20 mV and frequency (f) of 9.02 or 60.05 Hz being superimposed onto the DC cyclic ramp which had a scan rate (v). DC cyclic voltammetric data was also obtained with this instrument. A conventional three-electrode potentiostated system was used to collect experimental data at 25°C. GC (diameter = 3.0 mm, CH Instruments) or BDD (diameter = 1.0 mm) were used as the working electrode, Ag/AgCl (3 M NaCl) as the reference electrode and platinum wire as the auxiliary electrode with 1.0 M HClO 4 being the supporting electrolyte. The BDD disc electrodes were kindly supplied by Professor Julie Macpherson and prepared by the method described in reference [23] using materials provided by Element 6. In the preparation of the POM modified electrodes, the GC or BDD surface was polished with a 0.3 μm alumina aqueous slurry on a polishing cloth, rinsed with water, sonicated and rinsed again with water. The cleaned GC and BDD electrodes were then dried under a stream of nitrogen gas and placed in a 1 mM phosphomolybdic acid (H 3 PMo 12 O 40 ) aqueous solution for 1 min. The resultant POM-modified electrodes were rinsed several times with water to remove unbound solid before being used in electrochemical experiments.

Simulations and Data Analysis
In the automated parameterization of the 21 variables present in the model, possible values available for k 0 , α R u , C dl , and Γ were constrained to confine the extent of parameter space search and hence computational time needed to complete each data optimization exercise. Sensible range restrictions for these parameters were easily identified based on knowledge gained in the study in 1 M H 2 SO 4 at a GC electrode. [9] Significantly, this earlier study revealed that much tighter and carefully considered restrictions on E rev ranges are needed to avoid physically unrealistic values been obtained via data optimization. To achieve appropriate ranges, simulations of FTACV voltammograms assuming all processes were reversible were initially undertaken with MECSim software and heuristic comparison with experimental data used to estimate approximate E rev values for each of the six electron transfer reactions using methodology described in. [9] Parameter estimation with each data set was undertaken five times with a randomly selected starting point in the parameter space selected for each variable. Values for the 21 parameters reported for each experiment are those obtained from the data optimization having the minimum objective function. Full details of the parameterization modus operandi employed with the automated date optimization methods are available elsewhere. [9] In FTACV, the experimental and simulated total current time domain data were converted to the frequency domain to generate the power spectrum. [10,24] Frequencies corresponding to the AC harmonics and the aperiodic DC component were then selected from the power spectrum. Band filtering followed by inverse Fourier transformation provided access to resolved aperiodic DC and first, second, and third harmonic AC components which could then each be analyzed independently.

The Model
In order to generate the simulated data to mimic the AC voltammetry of the 3 surface confined processes summarized in Eq. 1-3, a combination of relationships accommodating six quasi-reversible electron transfer steps, an adsorption isotherm, background current and uncompensated resistance containing a total of 21 unknown parameter values are required. However, there are alternatives that could be used in almost every theoretical relationship employed in the model. In this sense, the model itself represents an additional or twenty second parameter containing uncertainty. The relationships employed in this study along with plausible alternatives that might rationalize systematic differences identified in experimental and simulated data are as follows: * Electron transfer model -the Butler-Volmer relationship was used for each of the six electron transfer steps which requires the introduction of six sets of E rev , k 0 and α parameters. Alternately, the mathematically more sophisticated Marcus-Hush relationship could have been used. [3] Neither thermodynamic nor kinetic dispersion that can occur with surface confined reactions [25,26] were included in the electron transfer model.
* Acid-base and ion-pairing chemistry -all proton transfer and ion-pairing reactions coupled to electron transfer were assumed to be diffusion controlled and hence reversible and subject to a solely thermodynamic description as is E 0 . On this basis, the unknown acid equilibrium constants (pKa) and ion-pairing constants have been combined with the E 0 values. Effectively this is a simplification which means only six E rev thermodynamically relevant parameters with their six k 0 and α values need to be parameterized. This implies that the individual E 0 , k 0 , α, pKa and ion pairing constant values remain unknown and that reported E 0 , k 0 , α, actually only represent apparent or E 0 app , k 0 app and α app values. No double layer correction was applied to the electrode kinetics.
* Adsorption isotherm -the Langmuir model was used, which means that a surface coverage parameter Γ is needed. This implies that all surface confined polyoxometalate moieties present in monolayer or sub-monolayer coverage are electrochemically identical in their behavior and act independently of each other. Other isotherms are available such as Frumkin or Temkin. [3] The simulation also assues that no dissolution of surface bound POM occurs during the course of the voltammetric experiment.
* Uncompensated resistance -Ohms law used (uncompensated resistance R u parameter introduced).
* Background current was modelled assuming it arises solely from double layer capacitance. Modelling undertaken assumes a simple R u C dl time constant applies at each potential and that C dl is independent of potential. In fact C dl depends on potential. Furthermore, other imperfections in the use of the simple R u C dl time constant model have been pointed out theoretically by Myland and Oldham. [27] Additionally, background current at a polyoxometalate modified electrode could exhibit far from simple double layer behavior for other reasons. If the polyoxometalate is at sub-monolayer levels as will emerge to be the case, uncovered heterogeneous GC or BDD surface as well as the modified part of the electrode surface which exists in six redox states will contribute to the background current to provide a potentially complex scenario relative to the modelled one. Finally, there may be another special feature of this particular polyoxometalate. When in contact with carbon materials, H 3 PMo 12 O 40 is said to form an electrochemical supercapacitor, [28] which in principle implies the presence of pseudo-capacitance via integration of double  1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57 layer and redox activity. However, caution should be applied in interpreting the significance of this phenomenon as ambiguities in the terminology and the implications of pseudocapacitance are widespread in the literature. [29] * Other parameters present in the model (see theory below) such as AC amplitude, DC scan rate, electrode area, temperature etc., are assumed to be accurately known. Parameters such as AC frequency, could have been included in the simulations and parameterized as part of the data optimization exercise to confirm their fidelity, but were not.
In summary, the model chosen for simulation of AC voltammograms requires a total of 21 parameters to be estimated by automated data optimization under circumstances where there are possible imperfections and uncertainties in almost every aspect of the model. In the extensive field of electrochemistry of surface confined enzymes, coupled electron and proton transfer reactions with similar characteristics to the POM reactions considered in this paper occur. Not surprisingly, imperfections in models used to describe the electrochemistry of surface confined enzymes have also been shown to be important as summarized in a recent review. [26] The mathematics underpinning the theoretical relationships employed in the model used to simulate the AC voltammetry in order to mimic the overall six electron reduction of surface confined ½PMo 12 O 40 � 3À at GC and BDD electrodes along with features of the data optimization protocols used are described in detail in. [9]

Preliminary Estimate of Reversible Potential and Other Parameter Values via Heuristic Data Analysis
In order to obtain an approximate estimate of the six E rev values and limit the range of values permitted in the fully automated data optimization procedure (see above and [9]), a heuristic form of data analysis was employed initially in which the experimentalist varied the unknown parameter combinations in simulations and decided which set of E rev values gave satisfactory agreement between experimental and simulated data obtained with frequencies of 9.02 and 60.05 Hz. In brief, to make this very tedious heuristic method tractable and only require modelling of five parameters rather than 21, all six electron transfer processes in MECSim generated simulated AC voltammograms were assumed to be fully reversible (all k 0 values used were very large and set at 1 × 10 4 s À 1 and all α values set at 0.500). Additionally, for further simplicity, the two electron transfer processes making up the significantly sharper middle process II were assumed to have identical E rev values for reasons given in reference [9]. To estimate C dl it was noted that in the fundamental AC harmonic derived from FT analysis that a potential (time) region exists prior to onset of faradaic current for process I that was assumed to be purely capacitive. By matching experiment and theory for this faradaic current free potential (time) region in the fundamental AC harmonic, C dl was estimated as a constant that was assumed to be valid at all potentials, although in reality it is strongly potential dependent. Analysis of the R u C dl time constant in this faradaic free potential region showed that the uncompensated resistance was sufficiently small to have negligible impact on the voltammetry. Finally, the surface coverage parameter Γ was estimated initially by integration of DC current time-data to give the charge associated with reduction of the POM and use of Faraday's Law [3] for a series of chemically modified electrodes prepared as described above. The charge consumed in these estimates of Γ was calculated from the sum of that derived from processes I and II (n = 4) or process III (n = 2) after background (capacitance current) correction of a DC linear sweep voltammogram. The value of Γ calculated in this way varied with each preparation of POM modified electrodes but was always in the sub monolayer range of 1.5 to 3.5 × 10 À 11 mole cm À 2 . However, since all processes are assumed to be reversible in this initial exercise and the impact from iR u drop assumed to be minimal, Γ was subsequently treated as a scaling factor which could be varied until simulated faradaic current magnitudes matched closely to the experimental ones in this preliminary heuristic form of data analysis. The 21 parameter model is now reduced in the preliminary heuristic data analysis exercise to a manageable one with only five E rev parameters remaining that require estimation. Approximate E rev values in mV versus Ag/AgCl obtained via this preliminary data analysis exercise at a BDD electrode in 1 M HClO 4 with a frequency of 60.05 Hz were E 0 1 ¼ 364 and E 0 2 ¼ 331 mV (for process I), E 0 3 ¼ E 0 4 ¼ 230 mV (for process II) and E 0 5 ¼ À 11 and E 0 6 ¼ À 29 mV (for process III). The corresponding values at the GC electrode were E 0 It is important to ensure that parameters finally reported by fully automated computational means are physically reasonable. Parameter values used in a heuristic form of data analysis by an experienced experimentalist should inevitably be chemically rational. A value of R u ¼ 50 Ohm for both GC and BDD modified electrodes was used arbitrarily in simulations selected heuristically although the voltammetry is not very sensitive to the value of this parameter provided it is less than about 500 Ohm, as confirmed by subsequent automated data optimization results presented below. The relevant preliminary BDD simulations employed C dl ¼ 6 � 10 À 6 Farad cm À 2 along with G ¼ 2:9 � 10 À 11 mole cm À 2 (9.02 Hz) and 2:3 � 10 À 11 mole cm À 2 (60.05 Hz) case with those at GC being C dl ¼ 8 � 10 À 6 (9.02 Hz) or 6:5 � 10 À 6 F cm À 2 (60.05 Hz) with G ¼ 3:4 � 10 À 11 (9.02 Hz) or 1:5 � 10 À 11 mole cm À 2 (60.05 Hz) based on the geometric area. These physically sensible parameter estimates employed as knowns in the preliminary analysis by the heuristic method provide a valuable reference point for verification of the fidelity of the far more sophisticated data optimization method, where the existence of local minima can play havoc with the reported outcome. While the automated data optimization method used to estimate all 21 parameters simultaneously will be shown below to generate a superior set of E rev values and other parameters than those estimated heuristically, the total AC current and fundamental, second and third AC harmonics still match very well with experimental data. This outcome highlights the problem of establishing a unique solution to an electrochemical mechanism as complex as the one being quantified in this study.

Automated Date Optimization using Low (9.02 Hz) Frequency Data and the Quasi-Reversible Electron Transfer Model
Fully automated data optimization was initially applied at a low frequency of 9.02 Hz to the total alternating current data obtained for the overall six electron reduction of surface confined ½PMo 12 O 40 � 3À at both GC and BDD electrodes in 1 M HClO 4 assuming each electron transfer process is quasireversible. Parameters deduced from this exercise for the parameter set providing the minimum objective function are summarized in Table 1. Figures 3a and 3b for GC and BDD electrodes respectively, provide a comparison of experimental total current AC voltammograms and those simulated with the estimated parameters provided in Table 1. Clearly, the fit of experimental data is substantially superior in the case of the BDD electrode. In particular, the GC background current clearly is not independent of potential and differs in magnitude at the start and end of the experiment, which is not predicted by the model. However, at both GC and BDD electrodes, all k 0 values are greater than about 5000 s À 1 and α values were unreliable and often were placed at the bounds that were given to the optimisation algorithm. This implies that all processes at both electrode surfaces lie very close to or are at the reversible limit, with the simulations being totally insensitive to α and k 0 values at best reflecting lower limit rather than an absolute value.
As expected, the value of R u estimated for the semiconducting BDD electrode obtained by data optimization is larger than at the more highly conducting GC. However, examination of all five data optimization parameter sets as well as the tabulated data one for the set with the minimum objective function include examples where R u is close to zero. The factor limiting the impact of R u is not its magnitude, but that under conditions prevailing in 1 M HClO4 the iR u (Ohmic) drop in this highly conducting electrolyte medium only has a minor impact on the voltammetric data. R u represents contributions for both the electrolyte and electrode resistances, with the latter being larger for BDD as correctly reflected in values for this parameter in Table 1. Additionally, errors in modeling background current noted above, particularly with the GC electrode will have an impact on value of R u reported via incorrect modelling of the R u C dl time constant. The primary conclusion reached for the R u value is that it is less than about 600 Ohm at BDD and less than 300 Ohm at GC with the small contribution from iR u readily allowing estimates of all six E rev to be highly meaningful.
In summary, parameterization of total current AC voltammetry at low frequency by automatic data optimization using a model having 21 unknown variables provides E rev values at a high level of confidence. However, only lower limits are available for k 0 , upper limits for R u , and uncertainties exist in C dl and R u due to inadequacies in the model used to describe the background current. The parameter α cannot be determined as the theory is totally insensitive to its magnitude.
The similar coverage at GC versus BDD is consistent with DC voltammetric data (Figure 1) where the area of the faradaic responses (charge per unit area) should accurately reflect the value of this parameter. The higher value of C dl at GC versus BDD also is reflected in the background current present in the DC voltammetry displayed in Figures 1a and 1b after taking into account differences in electrode areas.
Interesting differences in E rev and in pairs of potentials (E 0 1 and E 0 2 or DE 1 , E 0 3 and E 0 4 or DE 2 , E 0 5 and E 0 6 or DE 3 ) are detected and account for the shape variations in AC voltammograms obtained at GC and BDD electrodes. At the GC electrode differences in E 0 1 and E 0 2 are small (10 to 15 mV) for process I, potential inversion occurs for process II (E 0 3 less positive than E 0 4 by about 25 mV) and largest (about 35 mV) for process III. These results are consistent with processes II being sharper with a larger peak current than processes I and III at GC. At BDD, the shapes of processes I and III are similar to each other, with process II being sharper and having a larger peak current. At this electrode, the potential inversion deduced for process II (negative DE 2 value of about 20 mV) is consistent with its distinctly different shape and current magnitude relative to processes I and III where differences are about positive 30 mV. That is, as expected, the similar shapes of both processes I and III translate to similar DE 1 and DE 3 values. The shape variation is almost independent of frequency as are estimated values of E rev at each electrode surface as shown below.

Automated Data Optimization using Low Frequency (9.02 Hz) Data and the Reversible Electron Transfer Model
The above considerations suggest that low frequency data at both GC and BDD electrode could be modeled, plausibly equally well, assuming all electron transfer processes are reversible. To achieve this outcome all k 0 and α values were set at 1:0 � 10 4 s À 1 , and 0.500, respectively. Automated data optimization of the total AC current (f ¼ 9:02 Hz) voltammogram was then applied with only 9 variables (6 E 0 values plus R u , C dl and Γ) treated as unknowns. The values estimated in this manner at both electrodes are also included in Table 1. Theory and experimental data are also included in Figure 2 for the parameter sets obtained with the minimum objective functions. Again, the agreement between theory and experiment is significantly superior for the BDD data with E rev , C dl and Γ values being minimally dependent on whether the electron transfer process is modelled as reversible or quasi-reversible electron transfer. Some variation in R u is detected but as noted above when using the quasireversible model, the data are only weakly sensitive to this parameter. The substantial non-ideality of the background current at the GC electrode now translates into apparently model dependent R u , C dl and Γ values, but robustness in consistency of model independent E rev values is achieved.

Automated Data Optimization using Higher Frequency (60.05 Hz) Data
In principle, use of a much higher frequency of 60.05 Hz AC experiments should increase the prospects for quantifying the electrode kinetics parameters for fast electron transfer processes because departures from the reversible limit should be enhanced and hence the AC voltammetry is significantly more sensitive to variation in k 0 and α. The larger currents obtained at higher frequency can mean that iR u drop now becomes considerably more significant. The parameter values derived from models based on reversible and quasi-reversible electron transfer models at higher frequency are again included in Table 1 and simulated and experimental comparisons are displayed in Figure 3. As was the case with the lower frequency 9.02 Hz data, robustness in E rev parameter estimation is retained at 60.05 Hz with values differing slightly at GC and BDD electrodes but at both surfaces being essentially independent of frequency or as to whether electron transfer is assumed to be reversible or quasi-reversible. This outcome further establishes the fidelity of the recovery of all six E rev values by fully automated data optimization. At the higher frequency, success in matching the background current on the basis that this arises solely from that predicted by a potential independent double layer capacitance model is improved but is still imperfect and contributes to apparently frequency dependent scatter in C dl and R u parameter estimations. Again α remains non-determinable and k 0 remains too fast to quantify at 60.05 Hz except  Table 1.
that it now apparently lies in excess of about 1,000 s À 1 for all processes. On the basis of this higher frequency AC study, R u still cannot be reliably estimated but can only be confirmed to be sufficiently small so as not to induce significant iR u drop and hence R u remains of only low sensitivity with respect to impacting on the shape of the AC voltammograms.

Detection of Imperfections in the Modeling by FTACV
The FTACV method resolves the AC response into the aperiodic DC and fundamental and higher order harmonics, with the magnitude of the second and higher order AC harmonics being substantially enhanced with use of larger sine wave amplitudes. Studies from our laboratories commonly have used an amplitude of around 100 mV which for fast processes typically provides access at least 8 harmonics with good signal to noise ratios. In the present study, the sinewave magnitude was restricted to 20 mV to minimize overlap of processes I and II (large amplitudes broaden response) which limits data with acceptable signal to noise to just the first four harmonics. Analysis using the total current, which equates to the sum of the aperiodic DC and all AC harmonic components weights the outcome to the fundamental harmonic content which provides by far the largest contribution to the total signal. Thus, examination of outcomes from data optimization undertaken on the individual AC harmonic components or data subsets can allow systematic harmonic dependent deviations in experiment and theory to be detected which in turn can facilitate detection and understanding of implications of model inadequacies.
The aperiodic DC and individual AC harmonic content in envelope form are displayed in Figures 4 (GC) and 5 (BDD) at low 9.02 Hz frequency and in Figures 6 (GC) and 7 (BDD) at a higher frequency of 60.05 Hz along with theoretical predictions. This alternative format of data presentation makes visualization of differences between theory and experiment easier to detect. It is now apparent that simulated data predict a fully symmetrical faradaic response in each AC harmonic. Thus, faradaic peak heights observed from the negative potential sweep of potential (reduction) and potential (oxidation) sweep directions are predicted to be identical as indeed are the peak potentials. In contrast, experimental data for each harmonic at both GC and BDD electrodes show a systematic decrease in peak heights of faradaic current derived from oxidation versus those from reduction. Oxidation and reduction peak potentials also differ slightly, rather than being identical as predicted theoretically for a reversible process. Both GC and BDD electrodes are heterogeneous allowing for the presence of thermodynamic and kinetic dispersion which is not included in the model, although of course kinetic dispersion is not relevant for reversible electron transfer, as considered to be the case in this study. The Langmuir model for adsorption also may be imperfect with, for example, some interaction of adsorbed polyoxometalate anions occurring.  Table 1. A very obvious departure in theory versus experiment also is evident now with respect to the background current potential dependence in the aperiodic DC component. The theory used to model the total background current assumes an ideal capacitor which leads to both DC and AC responses (see Eq. 25 in [9] and related equations). The deviation is significantly larger in data obtained at GC electrodes than at BDD ones (Figures 4a  and 5a). In comparison, the fundamental harmonic background experimental current data are much closer in agreement with a model based on assuming the background current is fully capacitive in nature and independent of potential (Figures 4b  and 5b). Furthermore, background current data agreement at 60.05 Hz are superior to 9.02 Hz in the fundamental harmonic AC data (compare Figures 6a and b with 4a and b). The second and third harmonic AC background components at both kinds of carbon electrode, almost fully comply with the model in the sense that these higher order background currents are essentially zero as predicted from use of a model based on potential independent capacitance. [10] The DC scan rate used to collect data at both frequencies is similar so the predicted DC background current contribution is similar. However the magnitude of the background AC component in the fundamental harmonic, theoretically predicted to depend on the square route of frequency, is much larger at 60.05 Hz than at 9.02 Hz. Effectively, examination of the resolved DC and AC content in the FTACV display of data allows us to hypothesize that departures in the modeling and experiment background contributions are predominantly associated with poor compliance of the double layer model which is substantial on long DC time scales, but less so at shorter time domain AC timescales. Plausibly, interaction of quinones or other carbon surface functional groups present particularly on GC surfaces [30][31][32] with   Table 1. reduced forms of the polyoxometalate may generate pseudocapacitance that is more significant at the long DC time scale or at low AC frequencies. Many polyoxometalates, [28,[33][34][35][36] including H 3 PMo 12 O 40 [28] when mixed with carbon materials are also said to form supercapacitors. This feature may have implications in the present study, although caution is needed in the interpretation of the so called pseudo-capacitance.
We also note that a small amount of dissolution of POM, which increases with the extent of reduction may provide a contribution to mass transport by diffusion and contribute to some mismatch with our diffusionless model implied by use of the Langmuir adsorption isotherm. As noted in the Introduction, dissolution accompanying process IV was so rapid that this process was excluded from the parameterization exercise. Apparently, the extent of dissolution increases with an increase in the level of reduction of the polyoxometalate. On this basis it was considered likely that parameterization from just the reduction or negative potential sweep part of the data set could be more coincident with the model and indeed this is the case (data not shown). Nevertheless, reversible potential E rev data derived from use of this kind of restricted data set remain almost invariant with respect to the data set used.
The above observations imply that even superior agreement could be achieved between theory and experiment by removal of the DC and fundamental AC harmonic components that contain a large component of the charging current removed, particularly for the 9.02 Hz data set case. That is, if data evaluation is confined to just the second and higher order AC harmonics where no background current is present. This alternative use of restricted data sets for parameterization has been investigated using the reversible model for electron transfer. Significantly, all forms of data analysis show that the   Table 1. six estimated E rev values remain robustly independent of the data set used for parameterization. That is, yet again, confidence in the fidelity of estimates of reversible potentials is retained using FTACV. In summary, estimated E rev parameter values are almost independent of whether reversible or quasireversible electron transfer is used in modelling applied to either the full data set or a selected part of the data set in the data optimization exercise at either GC or BDD electrodes.
It is of interest to compare our results with those of Molina et al. [11,18] who used a square wave technique at this electrode to analyze the same three processes described in this paper in 1 M HClO 4 . The square wave data they provide in 1 M HClO 4 for BDD [11,18] have very similar characteristics to our ac voltammograms. In all reports, process II at BDD is sharper and has a larger peak response than processes I and III. However, instead of using simulation-experiment comparisons with data optimization for parameter estimation, as in our study, Molina et al in their initial study [11] used an analytical solution to the reversible theory to generate working curves that accommodated variation of the peak parameters (peak potentials, half-peak widths and peak heights) as the means for estimation of the six E rev values. Different reference electrodes [SCE (saturated calomel electrode) versus Ag/AgCl (3 M NaCl) with complex K + ClO 4À junction in SCE case] were used so that absolute values of E rev are difficult to compare. However differences in pairs of potentials (E 0 1 and E 0 2 or ΔE 1 , E 0 3 and E 0 4 or DE 2 , E 0 5 and E 0 6 or DE 3 ) should be directly comparable. Differences reported are À 17, À 85 and 10 mV for processes I, II and III respectively in the work of Molina et al whereas we estimate differences of about 30, À 30 and 30 mV respectively. We find shapes of processes I and III to be very similar with process II being sharper and having a larger peak current. In our case, the potential inversion deduced for process II (negative DE 2 value) is consistent with its distinctly different shape and current magnitude relative to processes I and III. Furthermore, and as expected, the similar shapes of both processes I and III translate to similar DE 1 and DE 3 values. The analytical data analysis method is rather complex with respect to accommodating overlapping processes and background current (iR u drop neglected) and it is difficult to know how the non-ideality we detected versus the reversible model will play out in this scenario. The aspect of agreement in both studies at BDD in 1 M HClO 4 is that differences in ΔE1, DE 2 and DE 3 are small. Furthermore, we note that simulated data when potentials of E 0 1 and E 0 2 are similar [9,[37][38][39] are not highly sensitive to ΔE variation which leads to uncertainties for experimentalists in heuristic forms of data analysis, as noted in this study. It is also probable that greater uncertainties are present in the method used by Molina et al. [11] compared to our data optimization approach. Clearly, the method of data analysis can play a role in parameter values reported for complex problems of the kind addressed in this study. Variation in the nature of the BDD in the work described by Molina et al who sourced their electrode from Windsor Scientific and our BDD electrodes prepared as in reference [23] using materials supplied by Element 6 also may contribute to differences. Boron concentration, homogeneity of boron, presence of graphitic sp 2 carbon, grain structure and whether BDD is oxygen or hydrogen terminated are all well known to play a critical role in electrode kinetic aspects of fully solution phase electrochemistry and in the background current. [23,40,41] An analogous situation applies to the GC electrode. In the case of surface confined material, the actual nature of the BDD, or GC surface, can impact on the reversible potentials for the ½PM 12 O 40 � 3À reduction processes, although this is not fundamentally possible for purely solution based reactions. In the second paper by Molina et al., [18] all processes are assumed to be quasi-reversible with k 0 1 ¼ k 0 2 , k 0 3 ¼ k 0 4 and k 0 5 ¼ k 0 6 , with all k 0 values being about 10 s À 1 at both electrodes. However, they described larger difference in reversible potentials and other features at BDD and GC than found in this work.
Finally, it is noted that data reported at GC in 1 M H 2 SO 4 we published previously [11] using the same data optimization method provided data consistent with reversible processes with E rev values of E 0 1 ¼ 0:362, E 0 2 ¼ 0:358, E 0 3 ¼ 0:215, E 0 4 ¼ 0; 240, E 0 5 ¼ 0:000 and E 0 6 ¼ À 0:016 V vs Ag/AgCl. Thus only small variations in reversible potentials are encountered in the slightly more acidic conditions, but less non-ideality in the modelling attributed to dissolution of reduced ½PMo 12 O 40 � 3À was evident in 1 M H 2 SO 4 than in 1 M HClO 4 .

Conclusions
Three well-defined surface confined overall two electron transfer reactions derived from closely spaced one-electron steps are observed in the DC and AC voltammetry of ½PMo 12 O 40 � 3À in 1 M HClO 4 . Parameterization of these processes has been achived by sinusoidal AC voltammetry at glassy carbon and boron doped diamond electrodes. The experiment-simulation comparison by a simplified reaction scheme required estimation of 21 parameters, which is an exceptionally large number in voltammetry. A heuristic form of data evaluation initially was undertaken to ensure physically realistic parameters were reported from computationally supported experiment-simulation comparisons. Automated data optimization approaches applied using the total current form of AC voltammetry allowed the six reversible potentials to be robustly recovered from data obtained at low (9.02 Hz) and higher (60.05 Hz) frequencies at both GC and BDD electrodes, irrespective of whether the electron transfer reactions were modelled as reversible or quasi-reversible. Examination of the harmonic content by the Fourier transformed version of AC voltammetry enabled imperfections in the model associated with faradaic and non-faradaic aspects of the voltammetry to be readily detected. After consideration of the implications of imperfections in the models, all processes were considered to be reversible within experimental uncertainty. Subtle differences in reversible potentials observed at GC and BDD electrodes account for electrode dependence in shapes of the overall two electron steps. Fits of experimental data were superior using BDD at both frequencies and for the higher frequency GC case than for low frequency data obtained with the GC electrode. Non-compliance with an assumed potential independent model of the double layer capacitance was considered to be a significant contributor to discrepancies 1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57 found with the GC electrode. It is concluded that the method of data analysis can play a role in parameter values reported for complex processes of the kind addressed in this study when impact of model non-ideality is parameterisation method dependent.

Data and Software Availability
The experimental data and software used for the automated data optimisation is available at https://github.com/martinjrobins/pom_project.