Supporting Information Kinetic Analysis on the Role of Bicarbonate in Carbon Dioxide Electroreduction at Immobilized Cobalt Phthalocyanine

The mechanism for carbon dioxide reduction (CO2RR) to carbon monoxide (CO) at immobilized cobalt phthalocyanine (CoPc) in aqueous electrolytes has been widely debated. In this work, we investigated the mechanism of CO2RR to CO on CoPc via experimental reaction kinetics coupled with model fitting. Unexpectedly, reactant order dependences and Tafel slopes deviate from commonly expected values and change depending on the testing conditions. For example, (1) the effect of bicarbonate deviates from power law kinetics and transitions from inhibitory to promotional with increasingly reductive potential, and (2) the CO2 order dependence deviates from unity at more-reductive potentials. We propose a kinetic model, chosen from more than 15 candidate models, that is able to quantitatively fit all of the experimental data. The model invokes (1) catalyst poisoning via bicarbonate electrosorption, (2) mixed control between concerted proton− electron transfer (CPET) and sequential electron transfer-proton transfer (ET-PT), and (3) both water and bicarbonate as kinetically relevant proton donors. The proposed model also predicts that the relative importance of the above factors changes depending on the reaction conditions, highlighting the potential downfalls of broadly applying reaction mechanisms that were inferred from kinetic data collected in a narrow range of testing conditions. This study emphasizes the importance of cohesively using kinetic data collected over a wide range of operating conditions to test and formulate kinetic models of electrocatalytic reactions.

Disc shaped electrodes with a diameter of 12 millimeters were punched from carbon paper (Toray, TGP-H-060, Fuel Cell Earth LLC) and inserted into a preheated tube furnace at 800 °C in static air for 10 min to generate oxygen-functionalized carbon paper (OxCP).
Catalyst ink was prepared by dispersing 10 mg of Cobalt(II) Phthalocyanine (CoPc, Strem Chemicals, Lot 24745700) in a mixture of 25 μL of 20 wt.% Nafion solution (D2021 Nafion dispersion, Fuel Cell Store) and 9974.5 μL N,N-dimethylformamide (DMF, Sigma Aldrich, 99.8%) with 15 minutes of sonication to yield a final solution containing a CoPc concentration of 1.75 x 10 -3 mol/L. This solution was then serially diluted in DMF to obtain catalyst inks containing CoPc concentrations ranging from 2 x 10 -7 to 4.38 x 10 -6 mol/L CoPc, with the corresponding Nafion content ranging from 6.25 to 12.5 × 10 -4 mg/ml (for the 4.38 x 10 -6 mol/L CoPc ink, the color of the solution was just barely blue). The specific loading of the catalyst, corresponding to the ink with CoPc concentration of 4.38 x 10 -6 mol/L, did not affect the measurement of turnover frequencies ( Figure S8), though the CoPc loading used for every measurement is tabulated in the supplemental excel document (SupplementalKineticData.xlsx).
The CoPc/OxCP working electrode was then prepared by dropcasting 15 μL of freshly sonicated ink onto the center of OxCP. The droplet instantaneously spread throughout the entire carbon paper without penetrating onto the underlying aluminum foil. The electrode was then oven-dried in air at 80 °C for 10 min. Fresh ink was prepared every time a new batch of electrodes was made, and electrodes were used within 24 hours of preparation. This is because there was a noticeable decrease in activity when electrodes were prepared from ink that was more than several days old or using electrodes that had been prepared even just one or two days prior. Catalyst loading was kept below 5.8 x 10 -11 mol/cm 2 .

1.2) Electrochemical testing
Electrochemical measurements were conducted in a previously reported, 1 customized threecompartment cell ( Figure S1) fabricated from polycarbonate, containing a counter electrode compartment, working electrode compartment, and gas compartment. The working electrode compartment was separated from the counter electrode compartment by an activated 2 Nafion membrane (Nafion117, Fuel Cell Store). The counter and reference electrodes were platinum foil (99.9% metals basis, Alfa Aesar) and Ag/AgCl leak-free reference (LF-2, Innovative Instrument Inc.), respectively. The leak-free reference was calibrated every day before use against a saturated calomel electrode, (CHI150, CH Instruments) in a S3 saturated potassium chloride solution (The potential at the SCE was assumed to be, according the manufacturer specification, +0.241V vs SHE at 25 o C). CoPc/OxCP electrode was placed between the working electrode compartment and gas compartment.
Prior to experiments, 1.75 mL of CO 2 -saturated electrolyte of specified composition was added into the working electrode compartment and counter electrode compartment, respectively. High-purity CO 2 gas (Airgas, 99.999%, 10 sccm) was controlled by an Alicat mass flow controller and introduced into the cell at atmospheric pressure; CO 2 gas entered the cell through the gas compartment, traversed the working electrode, and exited through the working electrode compartment and flowed directly to an online gas chromatograph (SRI Instruments, Inc., MG #5, Model 8610C). The anode compartment was not purged and was left open to atmosphere. The working compartment was purged with CO 2 for 5 min before electrochemical polarization. For CO 2 order dependence tests, the working compartment was purged with CO 2 for 30 min before polarization. Cells were disassembled and rinsed with MilliQ water, and a fresh electrode was used for each test.
The electrochemical measurements were controlled with a VMP3 Multi-channel potentiostat. Resistance between the reference and working electrodes was measured with Potential Electrochemical Impedance Spectroscopy (PEIS) and automatically compensated by 85%. The remaining 15% was not manually compensated, but because the cell resistance was always around 10 ohms, and operating currents were usually less than 1 mA (and were at most 4 mA), 15% of the IR drop should have been at most 1 mV, if not significantly less. For all experiments, the reaction was run for a total of 20 minutes, with gas products analyzed by an on-line gas chromatograph (SRI Instruments) every 5 minutes. The rates reported are calculated from average of the products detected at the 10, 15, and 20 minute marks. Data collection began after 10 minutes to ensure that the sampled gas headspace composition reflected steady state conditions. The gas chromatograph was equipped with a thermal conductivity detector (TCD) and flame ionization detector (FID), a methanizer, and a Hayesep-D column. Every 5 minutes, sample flowing through the sample loop was flowed to the Hayesep-D column, where it flowed to the FID and TCD. CO 2 eluted after CO, and flow to the FID was diverted during the times at which CO 2 was expected to elute. Column temperature was maintained at 100 o C. Each run was 3.5 minutes long. Sample chromatograms from the FID and TCD are provided in Figure S2.
For the experiments where CO 2 partial pressure was varied, the total flow rate and pressure were kept constant by adding in a diluent stream of nitrogen. 1M sodium bicarbonate electrolytes were prepared by bubbling CO 2 (Airgas, 99.999% research grade) through 0.5M sodium carbonate (Sigma Aldrich, 99.999% trace metal basis) solution (prepared with MilliQ ultrapure water) overnight. To test the order dependence on bicarbonate at constant ionic strength, sodium perchlorate (Sigma Aldrich, ≥ 98.0%) was used as supporting electrolyte. For order dependence tests, sodium bicarbonate concentration was varied by mixing the 1M sodium bicarbonate solution in an appropriate ratio with the 1M sodium perchlorate solution. For all kinetic tests, constant voltage vs the Ag/AgCl leak free reference (i.e., a pHindependent reference) was used.
Kinetic isotope effect experiments were performed by preparing electrolytes of Na 2 CO 3 and NaClO 4 identically as mentioned above, but using D 2 O (Oakwood Chemical, Lot No. 044626K20K) instead of H 2 O as the solvent.
Both the sodium bicarbonate and sodium perchlorate solutions were treated with Chelex resin (Chelex 100 sodium form, Sigma Aldrich, prepared according to previously reported procedures 3 ) before use as a precaution against electrolyte impurities, which are known to cause issues for bulk metal foil S4 catalysts. Notably, this treatment appeared to have very little effect, probably due to the mostly carbonaceous nature of the electrode used in this work ( Figure S16). Turnover frequency, defined as the reaction rate at each site, was calculated by assuming all deposited catalyst molecules were active. Thus, the formula for turnover frequency was: is the concentration of CO detected in the product stream in mol/cm 3 , is the flow rate of [ ] gas exiting the cell in cm 3 /s (notably, the amount of H 2 evolution is small and any water vapor from the electrolyte was also assumed to be small, so this is assumed to be the input flow rate, 10 sccm), is the loading of CoPc deposited on the electrode in mol/cm 2 and is the active (exposed) area of the electrode, assumed to reflect the size of the electrode hole in the working compartment (specifically, 1 cm 2 , though notably, 1 cm 2 is slightly smaller than the area of the entire physical electrode).

1) Details of equilibrium calculations
Modeling of equilibrium and transport in CO 2 RR systems was based on a previous report in the literature. 4 Equilibrium species concentrations were numerically calculated in Matlab. The following relevant constants were used:

Physical context
Value Eq constant for CO 2 (aq) + H 2 O --> H + + HCO 3 -10 -6.37 *1000 mol/m 3 Eq constant for HCO 3 ---> H + + CO 3 Where are the Sechenov coefficients for each ionic species (listed in Table S1 above), and is the , s i h G h specific gas constant of CO 2 , (in m 3 /mol) given by: 5 7 1.72 10 3.38 10 *( 298.15) Concentrations of all species were solved for by enforcing that ratios of activities equaled equilibrium constants for the reactions listed in Table S1

2.2) Notes on pH changes within order dependence tests
In general, order dependence tests can only be rigorously interpreted when all conditions except the desired species concentration are held constant. It is worth noting that in both the bicarbonate order dependence test at constant P CO2 (main text Figure 1a) and CO 2 order dependence at constant bicarbonate concentration (main text Figure 1c), the pH technically also changes when the concentration of the independent variable is being changed. This is because the equilibrium between CO 2 and carbonic acid that is described above in SI Section 2.1 implies that out of the three variables P CO2 , [H 3 O + ], and [HCO 3 -], it is only possible to independently manipulate 2. Therefore it is not physically possible to vary one of these parameters while holding the other two constant.
In this work, we exclude the kinetic relevance of [H 3 O + ] due to literature precedent that suggests proton concentrations near an electrode during CO 2 RR are generally too low to expect [H 3 O + ] to be a kinetically relevant species (full discussion in main text section 4.1.1). Therefore, in our kinetic modeling, the variation of pH in the bicarbonate and CO 2 order dependencies mentioned would only be significant for proton transfer steps that were assumed to be in thermodynamic equilibrium (where no matter which proton donor one assumes, pH will fall out of the equilibrium expression). However, from the results of the model fitting, the proposed model, M2, only has proton transfer steps that are assumed to be irreversible (i.e., kinetically limited and not equilibrated). Therefore, the variation of electrolyte pH with different [NaHCO 3 ] and P CO2 did not actually have any physical significance in the final proposed mechanism. We acknowledge that we cannot definitively exclude the possibility that H 3 O + is actually kinetically relevant. If it were, in fact, kinetically relevant, the bicarbonate and CO 2 order dependence tests could have different physical interpretations from the ones we have provided.
To illustrate the pH change more quantitatively, in the CO 2 order dependence test at costant [NaHCO 3 ] added , when P CO2 is changed by an order of magnitude, the pH of the bulk electrolyte also changes by about 1 pH unit. The value of [HCO 3 -] calc also changes, but to a much smaller extent. Thus, off unity CO 2 order dependencies could also be explained by some sort of pH dependence. Example calculated values are shown below in Table S3 For the bicarbonate order dependence tests at 1 atm of CO 2 , please see Table S2 to see how the solution pH changes.

S8
3) Mass transport analyses 3.1) Exclusion of mass transport limitations from invariance of TOF to loading Mass transport limitations were not considered (alternatively framed, bulk concentrations for kinetically relevant species were assumed to reflect interfacial species concentrations) mainly due to the experimental observation that turnover frequency (TOF) per molecule of CoPc deposited was, within certain limits of loading, invariant to the loading of CoPc used ( Figure S8). That means that if 5x more catalyst was deposited, then 5x more current was drawn. If mass transport limitations were at play, if 5x more current was drawn, we would expect steeper concentration gradients near the electrode (because the boundary flux of reactants is defined by the reaction rate at the electrode), which should change (most likely decrease) the experimentally measured TOF. Because no difference in rate per active site was observed, even at different overall geometrical current densities, we assumed that interfacial species gradients due to transport limitations were not kinetically relevant, meaning that for kinetically relevant species, the bulk concentrations were likely to be a good reflection of the interfacial concentrations.

3.2) 1D Cartesian transport model
To more rigorously test our assumption that bulk species concentrations were a good reflection of interfacial species concentrations, we implemented a 1D reaction-diffusion model. We assumed a flat, planar electrode with a static boundary layer of δ=100 microns. This is an approximation, because the actual system in this work is a porous carbon paper electrode with gas bubbling through it.
Within the static boundary layer, the following governing equation applies: Where the first term corresponds to diffusion, the second to migration, and the last to homogeneous reaction (in this case, buffering reactions and autoionization of water). When simplifying to one dimension, the above is simplified to: The above equation is additionally paired with an electroneutrality requirement: ) ( CO f j These equations were discretized using finite differences and numerically solved in Matlab. Diffusion and rate constants were obtained from previous literature. 4 The 1D transport model explicitly simulated the following species: K + , HCO 3 -, CO 3 2-, H + , OH -, CO 2 (aq). We found that the concentrations of CO 3 2-, H + , and OHdeviated the most from bulk concentrations, particularly at higher currents and lower buffer concentrations. This may serve as further justification for not explicitly modeling (CO 3 2-, H + , OH -) in the kinetic model, because these interfacial species concentrations appear to be very current-sensitive ( Figure S3).

3.3) 1D Spherical transport model
An additional explanation for the invariance of TOF to loading could be that the isolated sites are distant enough from each other that the transport problem should be modeled locally as spherical (or hemispherical, since the catalyst is on a solid support) to the catalyst. If local spherical transport to one catalyst does not affect transport to another (i.e., the boundary layers do not coincide), then even if there were transport limitations, one would still expect TOF to be invariant with loading. We considered this possibility and excluded it, the physical rationale being mainly because the radius of the catalytic center in CoPc is very small. When transport to an object is spherically symmetric, (in the absence of homogeneous reactions), an intrinsic mass transport boundary layer will form, and at steady state the concentration profile will take the following form: Where is a constant that can be resolved from some boundary condition at the catalyst surface 1 A and is the spherical distance coordinate. r We also more rigorously modeled 1D spherical diffusion-reaction problem by using the same procedure as SI Section 3.1 to numerically discretize and solve the following spherical conservation equation: Using the boundary conditions: On CoPc, we considered the radius of the active site ( ) to be 1 nm, which is a bit of an cat R overestimation, but at much smaller length scales, the continuum limit for these bulk transport equations breaks down. We found that at TOF values of 200 s -1 , we did not get any significant depletion of reactants at the interface. S11 4) Considered mechanisms 4

4.4) Consideration of RDS steps involving the second electron transfer
We also briefly considered the possibility of a second electron transfer being rate determining. As such, we define R6 and R7 according to Figure S4. We show the model fitting results for some mechanisms that considered these RDS possibilities in Table  S7.

4.5) Consideration and exclusion of alternative anionic poisoning species
We also considered the possibility that another anion other than bicarbonate could be poisoning the catalyst, and excluded this possibility due to the following analysis: the only other anions present in the electrolyte would be perchlorate (ClO 4 -), hydroxide (OH -), and carbonate (CO 3 2-). We were able to exclude perchlorate as a poisoning anion from some control perchlorate order dependence tests ( Figure  S14). We find it unlikely that hydroxide is poisoning the catalyst, because some studies have shown that at lower external buffer concentrations, the interfacial pH during CO2RR is more basic, 5 which means that a hydroxide poison should have a trend opposite to that of a bicarbonate poison.
It is possible that carbonate is a poisoning agent. However, it is difficult for us to experimentally discern this possibility from bicarbonate because (1) This effect could also be convoluted with a pH effect, because the only way to vary carbonate while keeping bicarbonate constant is to change the pH S14 and (2) Carbonate is a minority species and also carries a -2 charge, so from both concentration and electrostatic repulsion arguments, we think bicarbonate is a more likely poison (main text, Section 4.1.1).

5.1) Model fitting procedures
Parameter estimation was performed using Athena Visual software. Models were written in the mode of "Parameter Estimation"  "General form with user defined explicit models", since the expressions for rate could always be written explicitly in terms of the inputs and parameters. Parameter estimation was performed via Bayesian estimation, with a single response variable (which was ln(rate)). Gradient calculation was via central differences scheme, and lack-of-fit analysis was performed based on replicate experiments. Various tolerances defining convergence criterion were left as default values:  Convergence tolerance: 1 x 10 -1  Threshold pivot tolerance: 1 x 10 -4  Parameter tolerance: 1 x 10 -6  Objective function tolerance: 1 x 10 -6 Initial guesses were varied randomly, and several random initial guesses were tried when fitting each candidate model. Specifically, initial guesses for the and values were a random number between 0 and 1, and for the and rate and equilibrium constants, initial guesses were a random number between -3 and 4. For the best fit model, it generally did not matter what the initial guess was. For some rejected models, the initial guess did matter, and the model fit corresponding to the highest goodness-offit parameter was selected and reported.
Goodness-of-fit statistics were obtained from the Athena Visual-generated model discrimination and criticism report, which reports a variance ratio (which is referred to as the F-ratio in this work) and the sampling probability of greater ratio (which is referred to as alpha in this work). A theoretical discussion is provided in a textbook resource. 6 S15 5.2) Statistical graphs of best-fit model Figure S5. Normal probability plot for optimum parameter fitting of the proposed model. Since the plot appears linear, there is a good indication that errors between model and data are normally distributed Figure S6. Lag plot for the optimum parameter fitting of the proposed model. Since the plot appears random, that is a good indication that there is no correlation or drift within the data S16 6) Supplemental data 6.1) Negligible effect of Nafion in the electrode A small amount of Nafion was added to the electrode to serve as an ion conductive binder. However, Nafion has been implicated to block heterogenized transition metal macrocycle active sites in some electrocatalytic contexts. 7 To address this possibility, we performed control experiments to determine if Nafion had a similar detrimental effect in our system and found a negligible difference ( Figure S7).  Figure S8. Turnover frequency to CO as a function of various catalyst loadings. All points collected at the most reductive potential analyzed in this work, -1.12 V vs SHE, which is the potential that should be most susceptible to transport limitations. Total geometric current density is given by multiplying TOF CO with CoPc loading. Each point represents one experiment, solid lines are linear regressions to help guide the eye. . Rate data for hydrogen evolution reaction under some representative reaction conditions, using carbon paper with and without CoPc catalyst. Error bars represent standard deviation with n ≥ 2. Although the TOF towards HER is seen to increase with more reductive potentials in Error! Reference source not found., because the overall catalyst loading was generally lowered at more reductive potentials, the overall current towards HER was roughly constant with potential, as seen in this figure.  Figure S11. All data used to fit experimental Tafel slopes that were reported in Figure 1 of the main text. Unlike Figure 1 of the main text, all points represent one singlicate experiment; multiple points collected at the same conditions are shown as multiple points on the graph (only points that were not repeated in at least duplicate were shown in Figure 1 of the main text, but all points were considered for the fitting in this figure). For all bicarbonate concentrations, the region less reductive of -1 V vs SHE was considered for the Tafel fit. Lines represent linear regressions of the data. S19 6.5) CO 2 order dependence fitting -  Figure S12. All data used to fit experimental CO 2 order dependences that were reported in Figure 1 of the main text. Unlike Figure 1 of the main text, all points represent one single experiment; multiple points collected at the same conditions are shown as multiple points on the graph. Lines represent linear regressions of the data.  4 ] from 0 to 0.95M (thus, the ionic strength was not kept constant). The below data in Figure S14 is not corrected for "salting out" because the overhead partial pressure of CO 2 was kept at 1 atm, which means that the activity of dissolved CO 2 was constant in the electrolyte solution. We also considered the possibility that changing the electrolyte composition could cause changes in the junction potential between the leak-free reference and the electrolyte. This hypothesis seemed unlikely due to the different observed bicarbonate order dependence trends at less vs more reductive cathode potentials. We nonetheless tried calibrating the reference electrode against the ferri/ferrocyanide redox couple in the different electrolyte compositions. This was done by running cyclic voltammograms in solutions consisting of 5mM K 3 [Fe(CN) 6 ] dissolved in the appropriate electrolyte composition, using calcined carbon paper as the working electrode and Pt as the counter electrode. As seen in Figure S15 and Supplemental Table 8, we did not observe shifts large enough to account for the observed kinetic trends. (A 4mV shift for a reaction that has a Tafel slope of approx. 120 mV/dec represents about 8% error in measurement of TOF.) This data should be interpreted with caution because some white precipitate formed upon dissolution of 5mM K 3 [Fe(CN) 6 ] in the electrolyte solutions. Additionally, the ferri/ferrocyanide redox potential may also change depending on electrolyte composition.  Table 8. Tabulated values of reference electrode shift due to junction potential, determined from cyclic voltammograms S22 Finally, we also considered the possibility that the observed kinetic trends could be due to a difference in the amount/types of impurities present in the NaHCO 3 solution vs the NaClO 4 solution. Upon using higher purity sodium carbonate salts (99.999% rather than >99% Na 2 CO 3 , Sigma Aldrich) and also treating the electrolyte with chelating agent according to previously reported procedures 3 , we did not notice any change in the kinetic data ( Figure S16). We hypothesize that due to the mostly carbonaceous and CO2RR-inert nature of the electrode, trace metal impurities would be less likely to plate onto catalyst site and affect activity.

6.8) Addition of exogenous phosphate buffer
We also performed an experiment in which phosphate buffer (0.5 M, pH 7: 0.289 M Na 2 HPO 4 + 0.211 M NaH 2 PO 4 ) was added to the system while maintaining constant total ionic strength of 1 M ( Figure S17). Because H 2 PO 4has a pKa of 7.2 whereas HCO 3has a pKa of 10.33, we expect that when the two ions are present at similar concentrations, that H 2 PO 4should be the preferential proton donor. Therefore, the attenuated positive bicarbonate order dependence in the presence of phosphate buffer at more reductive potentials adds further evidence to support our mechanistic hypothesis.