We formulate a notion of very weak solution for the Poisson-Nernst-Planck system. The stationary system possesses a local monotonicity formula. Iterative application of the formula reveals improvement in estimates for ion density and potential, and leads to a local boundedness result. Local boundedness extends to steady-state systems for multiple ions and variable coefficients. The formulation applies to the related Keller-Segel system where stationary very weak solutions in two dimensions are regular. Examples illustrate how structure influences this regularity in higher dimensions. (C) 2014 Elsevier Ltd. All rights reserved.