Biomolecules, such as DNA, RNA, and proteins, are the building blocks of living systems. Each of these molecules is an assembly of atoms, with many carrying electric charges. Under normal circumstances, biomolecules are found in a solvent such as salted water. Highly charged parts of the biomolecules will polarize the solvent to produce ions that are mobile charges. The electrostatic interactions, together with covalent bonding and short- ranged van der Waals repulsion, give rise to the dominant forces that determine the dynamics and function of underlying biological systems. Understanding and simulating electrostatic interactions in biomolecular systems is the main goal of this dissertation. Both molecular dynamics and implicit-solvent descriptions have been the major modeling concepts for biomolecular interactions. Recent years have seen the development of a new class of efficient and accurate approaches, called variational implicit solvent modeling. This dissertation begins with this new framework and develops continuum models, mathematical analysis, and computational methods for simulating biomolecular interactions at different scales. First, the interfacial motion of dielectric boundaries is characterized and studied in the implicit- solvent description of biomolecules. Such motion is driven by the competition of mean curvature flow and the electrostatic free energy governed by the Poisson or Poisson-Boltzmann equation. Second, a new dynamic implicit solvent model is developed that introduces solvent fluid motion with fluctuations via the Landau-Lifshitz Navier- Stokes equations. This model significantly improves the existing variational implicit solvent model. And third, several multi-scale schemes combining a Brownian dynamics simulation of charged molecules with a continuum-based finite element description of coarse-grained concentrations are examined. Robust numerical algorithms are developed and extensive numerical computations are performed. They demonstrate the initial success of the theory. New directions in the mathematical modeling and computation of biomolecular systems are discussed