This paper presents a theoretical modeling study of starburst galaxies using two stellar population synthesis codes, PEGASE v2.0 and STARBURST99, to generate spectral energy distributions (SEDs) of young star clusters. The models also incorporate photoionization and shock physics using the MAPPINGS III code, which includes self-consistent treatment of dust physics and chemical depletion. The study uses optical diagnostic diagrams to assess the hardness of the EUV radiation field in these galaxies. The results show that warm infrared starburst galaxies have a relatively hard EUV field in the 1-4 Rydberg region. The PEGASE ionizing stellar continuum is harder in this range than that of STARBURST99, which is attributed to differences in stellar atmosphere models used for Wolf-Rayet stars. The Schmutz, Leitherer & Gruenwald (1992) atmospheres are considered more applicable to the starburst galaxies in the sample, but they do not produce the hard EUV field required by observations. The inclusion of continuum metal blanketing in the models may be a solution. Supernova remnant (SNR) shock modeling shows that the contribution by mechanical energy from SNRs to the photoionization models is <<20%. The models are used to derive a new theoretical classification scheme for starbursts and AGN galaxies based on the optical diagnostic diagrams. The study also discusses the role of Wolf-Rayet stars in determining the EUV spectrum and the effect of continuum metal opacities on the EUV ionizing continuum. The results suggest that a hard EUV field in the 1-4 Ryd region is necessary to model the starburst galaxies on the optical diagnostic diagrams. The study concludes that the inclusion of continuum metal opacities in the models may be a possible solution to the discrepancy between observations and models using Schmutz extended atmospheres.This paper presents a theoretical modeling study of starburst galaxies using two stellar population synthesis codes, PEGASE v2.0 and STARBURST99, to generate spectral energy distributions (SEDs) of young star clusters. The models also incorporate photoionization and shock physics using the MAPPINGS III code, which includes self-consistent treatment of dust physics and chemical depletion. The study uses optical diagnostic diagrams to assess the hardness of the EUV radiation field in these galaxies. The results show that warm infrared starburst galaxies have a relatively hard EUV field in the 1-4 Rydberg region. The PEGASE ionizing stellar continuum is harder in this range than that of STARBURST99, which is attributed to differences in stellar atmosphere models used for Wolf-Rayet stars. The Schmutz, Leitherer & Gruenwald (1992) atmospheres are considered more applicable to the starburst galaxies in the sample, but they do not produce the hard EUV field required by observations. The inclusion of continuum metal blanketing in the models may be a solution. Supernova remnant (SNR) shock modeling shows that the contribution by mechanical energy from SNRs to the photoionization models is <<20%. The models are used to derive a new theoretical classification scheme for starbursts and AGN galaxies based on the optical diagnostic diagrams. The study also discusses the role of Wolf-Rayet stars in determining the EUV spectrum and the effect of continuum metal opacities on the EUV ionizing continuum. The results suggest that a hard EUV field in the 1-4 Ryd region is necessary to model the starburst galaxies on the optical diagnostic diagrams. The study concludes that the inclusion of continuum metal opacities in the models may be a possible solution to the discrepancy between observations and models using Schmutz extended atmospheres.