We present the implementation of two grand canonical ensemble approaches in the open-source computational chemistry software CP2K that go beyond the existing canonical ensemble paradigm (Chai, Z.; Luber, S. Under Review). The first approach is based on implicit solvent models and explicit atomistic solute (electrode with/without adsorbed species) models, and includes two recent developments: Grand canonical self-consistent field (GC-SCF) method (J. Chem. Phys. 146, 114104 (2017)) allowing the electron number of the system to fluctuate naturally and accordingly with the experimental electrode potential. Planar counter charge (PCC) (J. Chem. Phys. 150, 041722 (2019), Phy. Rev. B 68, 245416 (2003)) salt model completely screening the net charge of the electrode model. In contrast with previous studies, in our implementation, the work function (absolute electrode potential if the potential drop at the electrolyte-vacuum interface is omitted) is the constrained quantity during a self-consistent field (SCF) optimization instead of the Fermi energy. The second approach (referred to as the two-surface method and the numerical titration method) (Phys. Rev. Lett. 88, 213002 (2002), J. Chem. Phys. 122, 234505 (2005), J. Am. Chem. Soc. 126(12), 3928–3938 (2004)) is based on fully explicit modeling of solvent molecules and ions. It is used to calculate the electron chemical potential corresponding to an equilibrium electrochemical half-reaction (M^(n+m)+ + n e^(-) ⇌ M^m+) which involves DFT molecular dynamics. This opens the way for forefront electrochemical calculations in CP2K for a broad range of systems.

In the self-consistent continuum solvation (SCCS) approach (J. Chem. Phys. 136, 064102 (2012)), the analytical expressions of the local solute-solvent interface functions determine the interface function and dielectric function values at a given real space position based solely on the electron density at that position, completely disregarding the surrounding electron density distribution. Therefore, the low electron density areas inside the solute will be identified by the algorithm as regions where implicit solvent exists, thereby resulting in the emergence of non-physical implicit solvent regions within the solute and even potentially leading to the divergence catastrophe of Kohn-Sham SCF calculations. We present a new and efficient SCCS implementation based on the solvent-aware interface (J. Chem. Theory Comput. 15, 3, 1996–2009 (2019)) which addresses this issue by utilizing a solute-solvent interface function based on convolution of electron density in the CP2K software package, which is based on the mixed Gaussian and plane waves (GPW) approach (Chai, Z.; Luber, S. To Be Submitted). Our implementation has been tested to successfully eliminate non-physical implicit solvent regions within the solute and achieve good SCF convergence, as demonstrated by test results for both bulk and surface models of liquid H₂O, titanium dioxide, and platinum.