The Poisson-Boltzmann equation (PBE) gives a mean field description of the electrostatic potential in a system of molecules in ionic solution. It is a commonly accepted and widely used approach to the modelling of the electrostatic fields in and around biological macromolecules such as proteins, RNA or DNA. The PBE is a semilinear elliptic equation with a nonlinearity of exponential type, a measure right hand-side, and jump discontinuities of its coefficients across complex surfaces that represent the molecular structures under study. These features of the PBE pose a number of challenges to its rigorous analysis and numerical solution. This thesis is devoted to the existence and uniqueness analysis of the PBE and the derivation of a posteriori error estimates for the distance between its exact solution and any admissible approximation of it, measured either in global energy norms or in terms of a specific goal quantity represented in terms of a linear functional. These error estimates allow for the construction of adaptive finite element methods for the fully reliable and computationally efficient solution of the PBE in large systems with complicated molecular geometries and distribution of charges.