The vibro-impact structures are introduced into the piezoelectric vibration energy harvesting (VEH) systems, which can broaden the operating bandwidth of the piezoelectric VEH systems and thus improve the harvesting efficiency. The process of vibro-impact not only causes energy loss during the impact process but also causes the local deformation between two impact bodies. Therefore, this manuscript proposes a class of stochastic non-classical inelastic impact VEH models, whose non-classical impact processes are described by Hertz contact forces. Firstly, the relationship between the total energy of the unperturbed system and the potential energy of the system at the barrier position is used to determine whether the impact phenomena have occurred or not, then used to derive the period and frequency of the system. Secondly, using the quasi-conservative stochastic averaging method, the averaged drift term and the averaged diffusion term of the piecewise smooth averaged Ito^ equation are derived. The stationary probability density functions (PDFs), the mean square output (MSO) voltage and the root mean square (RMS) voltage of the system response can be obtained from solving the corresponding Fokker-Planck-Kolmogorov (FPK) equation for averaged Ito^ equation. Finally, an example is provided to verify the effectiveness of the proposed method, and the effects of the model parameters on the MSO voltage and power conversion efficiency (PCE) are analyzed. The relationship between the MSO voltage and the PCE is also analyzed. It is of great significance for the design and optimization of stochastic VEH systems.