Continuum modeling of electrostatic interactions based upon the numerical solutions of the Poisson-Boltzmann equation has been widely adopted in biomolecular applications. To extend their applications to molecular dynamics and energy minimization, robust and efficient methodologies to compute solvation forces must be developed. In this study, we have first reviewed the theory for the computation of dielectric boundary forces based on the definition of the Maxwell stress tensor. This is followed by a new formulation of the dielectric boundary force suitable for the finite-difference Poisson-Boltzmann methods. We have validated the new formulation with idealized analytical systems and realistic molecular systems.