Photovoltaic laser beam degradation in lithium niobate (LiNbO3) has been investigated through a finite-differences beam propagation method for nonlinear media. The simulations use a two-center model (Fe2+↔Fe3+, NbLi4+↔NbLi5+) that has been recently proved to be necessary to successfully describe the photorefractive effect in nominally pure LiNbO3. Refractive index profiles and intensity and phase beam profiles have been calculated for a wide intensity range and several material distances. A good agreement is obtained on comparing simulations and experimental data. This includes self-defocusing for moderate intensities and the complex beam profile structures appearing for high intensities, providing a complete description of this nonlinear optical phenomenon.