We are reporting development and initial applications of the light-induced and thermally induced spectral diffusion modeling software, covering nonphotochemical spectral hole burning (NPHB), hole recovery, and single-molecule spectroscopy and involving random generation of the multiwell protein energy landscapes. The model includes tunneling and activated barrier-hopping in both ground and excited states of a protein-chromophore system. Evolution of such a system is predicted by solving the system of rate equations. Using the barrier parameters from the range typical for the energy landscapes of the pigment-protein complexes involved in photosynthesis, we (a) show that realistic cooling of the sample must result in proteins quite far from thermodynamic equilibrium, (b) demonstrate hole evolution in the cases of burning, fixed-temperature recovery and thermocycling that mostly agrees with the experiment and modeling based on the NPHB master equation, and (c) explore the effects of different protein energy landscapes on the antihole shape. Introducing the multiwell energy landscapes and starting the hole burning experiments in realistic nonequilibrium conditions are not sufficient to explain all experimental observations even qualitatively. Therefore, for instance, one is required to invoke the modified NPHB mechanism where a complex interplay of several small conformational changes is poising the energy landscape of the pigment-protein system for downhill tunneling.