First-order system least squares (FOSLS) is a recently developed methodology for solving partial differential equations. Among its advantages are that the finite element spaces are not restricted by the inf-sup condition imposed, for example, on mixed methods and that the least-squares functional itself serves as an appropriate error measure. This paper studies the FOSLS approach for scalar second-order elliptic boundary value problems with discontinuous coefficients, irregular boundaries, and mixed boundary conditions. A least-squares functional is defined and ellipticity is established in a natural norm of an appropriately scaled least-squares bilinear form. For some geometries, this ellipticity is independent of the size of the jumps in the coefficients. The occurrence of singularities at interface corners, cross-points, reentrant corners, and irregular boundary points is discussed, and a basis of singular functions with local support around singular points is established. A companion paper shows that the singular basis functions can be added at little extra cost and leads to optimal performance of standard finite element discretization and multilevel solver techniques, also independent of the size of coefficient jumps for some geometries.