A new SPN theory valid for piecewise homogeneous regions is rigorously derived without mathematical approximations. The resulting SPN equations are the same as the conventional ones, but both the interface and the boundary conditions are different from the ad hoc conditions conventionally assumed. The new theory does not need the assumption of locally 1D behavior near a surface.Gelbard pointed out in 1960 that SPN is rigorously valid for an infinite domain 3D problem with constant total cross-section and arbitrarily distributed isotropic sources, and speculated that SPN could be a good approximation for practical applications. The Gelbard’s example does not involve any interface or boundary condition. For a generic problem with piecewise homogeneous regions, even though the SPN equations can still be rigorously derived for each homogeneous region, the interface or boundary condition poses the problem for generalizing the rigorous SPN theory. In this paper we point out that in each of the homogeneous regions, the solution to the SPN equations can be used to rigorously reconstruct an angular flux distribution which is a particular solution to the PN equations. The complete PN angular flux solution is the sum of this particular SPN angular flux solution and a set of additional functions that are solution to the homogeneous part of the PN equations. This particular SPN angular flux solution alone does not have enough independent functions to fulfill all the PN interface or boundary conditions. In this case, one can only satisfy a subset of the conditions. Because this subset of conditions cannot be decoupled from the full set of conditions, the whole core SPN solution so obtained cannot provide exactly the same scalar flux solution as that given by the PN equations. Nevertheless if this particular SPN angular flux solution is to be used alone, its corresponding interface and boundary conditions can be explicitly derived without introducing mathematical approximations. The resulting conditions are different from the conventionally used ad hoc conditions for SPN. The new interface and boundary conditions involve additional terms of higher order derivatives of the SPN solution and are consistent with what Selengut found for SP3 in 1970. The arguments and techniques needed to establish the above theory can in fact be all deduced from the detailed PN theory formulation for piecewise homogeneous regions that was worked out by Davison in 1957 and 1958.This new SPN theory is valid for the multi-group case as well and can be generalized to the case of anisotropic scattering. It may be a good candidate for practical applications to 3D pin-cell meshed whole core calculations, if the pin-cell meshes can be properly homogenized. The new interface and boundary conditions can be implemented in conventional SPN codes but require more work to revise the coupling equations, resulting in new SPN nodal equations. The discontinuity factors needed for using homogenized meshes can be generated without ambiguity by using the particular SPN angular flux solution as well. A physical interpretation for the SPN theory is suggested. It is believed that the missing steps are now filled to have a sound and complete SPN theory.