Abstract

The standard Poisson-Boltzmann (PB) model for molecular electrostatics assumes a sharp variation of the permittivity and salt concentration along the solute-solvent interface. The discontinuous field parameters are not only difficult numerically, but also are not a realistic physical picture, as it forces the dielectric constant and ionic strength of bulk in the near-solute region. An alternative to alleviate some of these issues is to represent the molecular surface as a diffuse interface, however, this also presents challenges. In this work we analyzed the impact of the shape of the interfacial variation of the field parameters in solvation and binding energy. However we used a hyperbolic tangent function tanhkpx$$ \left(\tanh \left({k}_px\right)\right) $$ to couple the internal and external regions, our analysis is valid for other definitions. Our methodology, restricted to the linear PB, was based on a coupled finite element (FEM) and boundary element (BEM) scheme that allowed us to have a special treatment of the permittivity and ionic strength in a bounded FEM region near the interface, while maintaining BEM elsewhere. Our results suggest that the shape of the function (represented by kp$$ {k}_p $$) has a large impact on solvation and binding energy. We saw that high values of kp$$ {k}_p $$ induce a high gradient on the interface, to the limit of recovering the sharp jump when kp ->infinity$$ {k}_p\to \infty $$, presenting a numerical challenge where careful meshing is key. Using the FreeSolv database to compare with molecular dynamics, our calculations indicate that an optimal value of kp$$ {k}_p $$ for solvation energies was around 3. However, more challenging binding free energy tests make this conclusion more difficult, as binding showed to be very sensitive to small variations of kp$$ {k}_p $$. In that case, optimal values of kp$$ {k}_p $$ ranged from 2 to 20.