Continuum electrostatics methods have become increasingly popular due to their ability to provide approximate descriptions of solvation energies and forces without the expensive sampling required by all-atom solvent models. In particular, the Poisson–Boltzmann equation (PBE) provides electrostatic potentials, solvation energies, and forces by modeling the solvent as a featureless dielectric material and the mobile ions as a continuous distribution of charge. Polar solvation forces and energies obtained from the PBE are often supplemented with simple solvent-accessible surface area (SASA) models of nonpolar solvation. However, while polar and nonpolar continuum models have been assessed on their ability to reproduce global properties, such as solvation free energies, their ability to provide accurate representations of local solvation properties such as forces has not previously been adequately studied. We have developed efficient software for describing polar biomolecular solvation by solving the Poisson-Boltzmann equation using multigrid and adaptive finite element methods. Additionally, we have implemented new models to describe nonpolar solvation phenomena. These models have been used to study solvation forces for protein, RNA, and alkane systems. In particular, we have performed comparisons of continuum and all-atom representations of solvation forces for these very different molecular systems in order to assess the performance of continuum models in the presence of widely varying charge densities. The results of these comparisons show that current implementations of the PBE are capable of generating polar solvation forces that correlate well with explicit solvent forces for protein systems but provide significantly less accurate representations of polar solvation forces for RNA systems. Conversely, SASA-based nonpolar forces are found to have no significant correlation with nonpolar explicit solvent forces for either protein or RNA molecules. Good correlation between explicit and continuum nonpolar forces is only obtained when area, volume, and attractive dispersion forces are included in the continuum model. We discuss the implications of these studies in the context of molecular simulation as well as the impact of this work on basic models for understanding experimental observations of biomolecular binding and folding.