Lorena A. Barba group

Accepted: Biomolecular electrostatics solver using Python and GPUs

Surface charge plot on a lysozyme molecule, obtained using the PyGBe code.

Surface charge plot on a lysozyme molecule, obtained using the PyGBe code.

Submitted: 20 September, 2013. Accepted: 24 October, 2013.

This paper presents a study of the effect of solvent-filled cavities and Stern layers in a biomolecular electrostatics solver based on a boundary integral formulation. The tool for this study was the PyGBe code: a solver for biomolecular electrostatics using Python, GPUs and boundary elements. To determine the impact on accuracy of including or not features like cavities and Stern layers, we compared with a community code based on finite-difference solution of the Poisson-Boltzmann equation. Rigorous comparisons like the one we offer are scarce and we bring  to bear in this work our experience with verification and convergence analysis of numerical methods.


The continuum theory applied to bimolecular electrostatics leads to an implicit-solvent model governed by the Poisson-Boltzmann equation. Solvers relying on a boundary integral representation typically do not consider features like solvent-filled cavities or ion-exclusion (Stern) layers, due to the added difficulty of treating multiple boundary surfaces. This has hindered meaningful comparisons with volume-based methods, and the effects on accuracy of including these features has remained unknown. This work presents a solver called PyGBe that uses a boundary-element formulation and can handle multiple interacting surfaces. It was used to study the effects of solvent-filled cavities and Stern layers on the accuracy of calculating solvation energy and binding energy of proteins, using the well-known APBS finite-difference code for comparison. The results suggest that if required accuracy for an application allows errors larger than about 2%, then the simpler, single-surface model can be used. When calculating binding energies, the need for a multi-surface model is problem-dependent, becoming more critical when ligand and receptor are of comparable size. Comparing with the APBS solver, the boundary-element solver is faster when the accuracy requirements are higher. The cross-over point for the PyGBe code is in the order of 1-2% error, when running on one GPU card (NVIDIA Tesla C2075), compared with APBS running on six Intel Xeon CPU cores. PyGBe achieves algorithmic acceleration of the boundary element method using a treecode, and hardware acceleration using GPUs via PyCuda from a user-visible code that is all Python. The code is open-source under MIT license.

Open data



This research is made possible by support from the Office of Naval Research, Applied Computational Analysis Program, N00014-11-1-0356. LAB also acknowledges support from NSF CAREER award OCI-1149784.