Methods based on Green's functions have long been a workhorse for homogeneous boundary value problems, as they lead to integral equations posed on the boundary of a region of interest, with additional special advantages accruing for exterior problems. Persistent challenges include the need to evaluate singular integrals with kernels that decay slowly and require special care for efficient computations of long-range interactions. But even more fundamentally, the methods have not seen success for inhomogeneous problems for which methods typically involve one or several volume integral operators (VIOs). We will discuss a class of numerical methods for VIOs that rely on analytic regularization: they use Green's identities to regularize the singular kernels in the volume integrals and lead to high-order accuracy even with the use of standard, singularity-oblivious triangle/tetrahedral quadratures. A complete error analysis, over 2D and 3D meshes possibly containing curved elements, for several key VIOs of progressively increasing singularity strength that arise in applications will be discussed. Numerical examples will be presented in the context of (1) IMEX time-stepping methods for nonlinear PDEs and (2) variable-coefficient scattering problems arising in inhomogeneous media. For the latter, a key challenge involves the incorporation of efficient volume preconditioners into iterative solvers, while retaining high-order accuracy in the presence of material discontinuities.