Elucidating the Molecular Mechanism of CO2 Capture by Amino Acid Ionic Liquids

Amino acid ionic liquids have received increasing attention as ideal candidates for the CO2 chemisorption process. However, the underlying molecular mechanisms, especially those involving proton transfer, remain unclear. In this work, we elucidate the atomistic-level reaction mechanisms responsible for carbamate formation during CO2 capture by amino acid ionic liquids through explicit ab initio molecular dynamics augmented by well-tempered metadynamics. The resulting ab initio free-energy sampling reveals a two-step reaction pathway in which a zwitterion, initially formed from the reaction between the anion of serine and CO2, undergoes a kinetically facile intermolecular proton transfer to the O atom of the COO– moiety in the nearby serine. Further analysis reveals that the significantly reduced free-energy barriers are attributed to enhanced intermolecular interaction between the zwitterion and serine, thus facilitating the kinetic favorability of the proton transfer, which governs the overall CO2 capture mechanism. This work provides valuable insight into the important mechanistic and kinetic features of these reactions from explicit condensed phase ab initio MD free-energy sampling of the CO2 capture process.


Detailed Description of Computational Methods and Simulation Setup
Ab initio molecular dynamics (AIMD) simulations based on density functional theory were performed with the CP2K software. 1 To prepare the AIMD simulations, all-atom classical molecular dynamics simulations were first computed using the LAMMPS MD software. 2 The simulation box containing 15 amino acid ionic liquids with cation and anion pairs of cholinium and serine with five CO2 molecules was prepared using PACKMOL. 3 The cubic simulation box was run for 10 ns in the constant NPT ensemble, from which the equilibrated density is obtained. Thereafter, the simulation box was equilibrated under the NVT ensemble for 5 ns, followed by production runs for 15 ns. The final configuration from the production runs is used as the initial structure for AIMD runs. All AIMD simulations are computed at 300 K. The generalized gradient approximation (GGA) functional of revised Perdew, Berke, and Ernzerhof (revPBE) is used. 4,5 The conclusions drawn from this work are unaffected by choice of functionals (see below). Normconserving Goedecker-Teter-Hutter (GTH) pseudopotentials are employed to describe the interactions between ionic cores and valence electrons. 6 A hybrid Gaussian and plane-waves (GPW) method 7 is used in the QUICKSTEP molecule 8 , in which atom-centered Gaussian-type orbitals are used to describe the wave functions and an auxiliary plane wave basis set for expansion of the electron density. A triple-zeta Gaussian basis set with two polarization function sets (TZV2P) is employed; a plane-wave kinetic energy cutoff of 450 Ry is used. The Brillouin zone was sampled using only the gamma point. The semi-empirical dispersion corrections developed by Grimme (DFT-D3) with Becke-Johnson (BJ) damping is used to treat the long-range van der Waals interaction. 9 A timestep of 0.5 fs is used to integrate the equations of motion. AIMD simulations were run for 15 ps under the constant NVT ensemble, followed by production runs for 90 ps. Thereafter, well-tempered metadynamics simulations 10,11 in conjunction with PLUMED plug-in 12 were carried out for the free-energy sampling of reaction pathways. More details on metadynamics simulation setup, convergence tests, and transition state verifications can be found below.

Sensitivity of Functional Choices (Levels)
All simulations are prepared with ab initio molecular dynamics (AIMD) augmented with welltempered metadynamics. To check the sensitivity of the density functional theory (DFT) exchangecorrelation functionals on our simulations and the corresponding results, we also performed additional free-energy sampling calculations for the two-step reaction pathways discussed in the main text. Here, the generalized gradient approximation (GGA) level of theory of the revised Perdew, Burke, Ernzerhof (revPBE) functional with is compared with that of Becke, Lee-Yang-Parr (BLYP) functional 13,14 , metaGGA level of theory of revised Tao, Perdew, Staroverov, Scuseria (revTPSS) 15,16 and Strongly Constrained and Appropriately Normed (SCAN) 17 semilocal functionals. The corresponding Helmholtz free-energy barrier (under NVT ensemble) between the initial and the transition states for the two-step reaction pathways during CO2 capture process are reported in Table S1. As illustrated in Table S1, the relative free-energy barrier (FEB) is reported for the reactions of case (a), (b), and (c), as described in the main text. The conclusions made from our simulation results on the FEB are unaffected by the higher level of DFT functionals (metaGGA), and revPBE (GGA) functional used in our simulations is suitable. Table S1. Free energy barrier (FEB) for the reactions of case (a), (b), and (c) (i.e., FEBa, FEBb, FEBc) based on various functionals of BLYP, revPBE, revTPSS, and SCAN with the level of theory of GGA and metaGGA from AIMD simulations for the two-step reaction pathways for CO2 capture by amino acid ionic liquids (AAIL) studied in this work.

Metadynamics Simulations Details -I) Setup
Chemisorption of CO2 by AAILs can be initiated from nucleophilic attack by a basic N atom in serine (anion) at an electrophilic C atom in CO2, forming zwitterion. Thereafter, the zwitterion may undergo three possible routes involving proton transfer for its conversion into carbamate products. The overall reaction pathways are described in Figure 1 Figure 1 of the main text). To investigate the thermodynamic and kinetic favorability of the reaction paths described herein, we first computed ab initio metadynamics. A representative simulation setup with 15 AAIL pairs of cholinium (cation) and serine (anion) with five CO2 molecules in a cubic periodic box is illustrated in Figure 2 of the main text.

Metadynamics Simulations Details -II) Convergence
The free-energy barriers for the reaction of CO2 by AAILs are calculated from AIMD with welltempered metadynamics simulations. The three reaction cases of (a), (b), and (c) are the followings: (a) intermolecular proton transfer to the O atom in serine, (b) intermolecular proton transfer to the N atom in serine, and (c) intramolecular proton transfer to the O atom within the zwitterion (as shown in Figure 1 of the main text). As shown in Figure S1, the convergence of free-energy barriers is obtained after 20 (x 10 5 ) number of timestep for AIMD-metadynamics for all cases of (a), (b), and (c). Figure S1. Free-energy barrier distributions during AIMD-metadynamics simulations for the twostep reaction pathways for CO2 capture by amino acid ionic liquids (AAILs) studied in this work. The three reaction pathways involving proton transfer (a) intermolecularly to a nitrogen atom in serine (Ninter, in blue star), (b) intermolecularly to an oxygen atom in serine (Ointer, in orange diamond), and (c) intramolecularly to an oxygen atom within zwitterion (Ointra, in green square) are indicated.
Thereafter, another convergence is confirmed for the AIMD well-tempered metadynamics by investigating the deposited Gaussian hill heights. For the well-tempered algorithm, the Gaussian hill height should decay to zero once fully converged, which is also shown and confirmed with Figure S2 around 20 (x 10 5 ) number of timestep. Figure S2. Gaussian hill height deposited distributions during AIMD-metadynamics with welltempered algorithim for the two-step reaction pathways for CO2 capture by amino acid ionic liquids (AAILs) studied in this work. The simulations are fully converged once Gaussian hill decays to zero.

Metadynamics Simulations Details -III) Transition state verifications
A large set of unbiased AIMD simulations were computed, in which at least 17 different molecular configurations near the identified transition state region during the CO2 reactions with AAILs were simulated with random Boltzmann distributed velocities. It was subsequently determined, through these simulation runs, that the unbiased MD simulations fell into reactant or product states as represented along the chosen CV space, after a true reactive trajectory had been identified.

Metadynamics Simulations Details -IV) CV specifications
A free-energy surface (FES) was obtained from AIMD coupled with well-tempered metadynamics for the carbamate formation during CO2 chemisorption process by AAILs. Here, the collective variables (CVs) chosen employ two collective variables. The first CV (x3) is a bond distance between the N atom in serine molecule and the C atom in CO2 molecule. A relatively larger value of x3 indicates that the serine and CO2 molecules are separated with no significant interaction between the two, while a comparatively smaller value of x3 indicates that the zwitterion is formed. The second CV is composed of two coordination numbers (CNs): (i) a bond distance-dependent coordination number between a proton and the base N atom in zwitterion (x1), and (ii) a bond distance-dependent coordination number between the proton and the base N atom in nearby serine molecule (x2). Here, a linear combination of the two CVs (x1-x2) can conveniently describe S6 protonation states of both the zwitterion and nearby serine molecule. A positive value (specifically a value of unity, 1, due to x1-x2 = 1-0 = 0) indicates that the proton is fully bound to the zwitterion, while there is no significant interaction between the proton and serine molecule; a negative value (specifically a value of -1, due to x1-x2 = 1-0 = 0) denotes that the proton is successfully transferred from the zwitterion to the nearby serine. The CNs were defined using a rational switching function for the purpose of ensuring differentiability and the values of the CVs were set to unity or negative unity. The switching function description of the second CV employing two CNs are the followings.  Figure S3, also confirming the convergence of AIMD well-tempered metadynamics.
S7 Figure S3. Collective variables (CVs) distributions during the AIMD-metadynamics with welltempered algorithm. CV1 is x3, mimicking the bond distance between N and C atoms, while CV2 is x1-x2, denoting a difference between two coordination numbers for protonation states.