Testing a heterogeneous polarizable continuum model against exact Poisson boundary conditions

The polarizable continuum model (PCM) is a computationally efficient way to incorporate dielectric boundary conditions into electronic structure calculations, via a boundary-element reformulation of Poisson's equation. This transformation is only rigorously valid for an isotropic dielectric medium. To simulate anisotropic solvation, as encountered at an interface or when parts of a system are solvent-exposed while other parts are assumed to be in a nonpolar environment, ad hoc modifications to the PCM formalism have been suggested, in which a dielectric constant is assigned separately to each atomic sphere that contributes to the solute cavity. The accuracy of this "heterogeneous" PCM (HetPCM) method is tested here for the first time, by comparison to results from a generalized Poisson equation solver. The latter is a more expensive and cumbersome approach to incorporate arbitrary dielectric boundary conditions, but one that corresponds to a well-defined scalar permittivity function, epsilon(r). We examine simple model systems for which a model function epsilon(r) can be constructed in a manner that maps reasonably well onto a dielectric constant for each atomic sphere, using a solvent-exposed dielectric constant epsilon_solv = 78 and a variety of smaller values to represent a hydrophobic environment. For nonpolar dielectric constants epsilon_nonp <= 2, differences between the HetPCM and Poisson solvation energies are large compared to the effect of anisotropy on the solvation energy. For epsilon_nonp = 4 and epsilon_nonp = 10, however, HetPCM and anisotropic Poisson solvation energies agree to within 2 kcal/mol in most cases. As a realistic use case, we apply the HetPCM method to predict solvation energies and pKa values for a blue copper protein. The HetPCM method affords pKa values that are more in line with experimental results as compared to either gas-phase calculations or homogeneous (isotropic) PCM results.