2009ApJ...704..715V

Motivated by the recent detection of a large number of embedded young stellar objects (YSOs) with mass accretion rates that are inconsistent with the predictions of the standard model of inside-out collapse, we perform a series on numerical hydrodynamic simulations of the gravitational collapse of molecular cloud cores with various initial masses, rotation rates, and sizes. We focus on the early Class I stage of star formation when circumstellar disks are exposed to high rates of mass deposition from infalling envelopes. Our numerical modeling reproduces the large observed spread in accretion rates inferred for embedded YSOs in Perseus, Serpens, and Ophiuchus star-forming regions by Enoch et al., yielding 37%-75% of objects with "sub-Shu" accretion rates M_{dot} ≲ 10^–6 M_☉/yr and 1%-2% of objects with "super-Shu" accretion rates M_{dot} > 10^–5 M_☉/yr. Mass accretion rates in the Class I stage have a lognormal distribution, with its shape controlled by disk viscosity and disk temperature. The spread in M_{dot} is greater in models with lower viscosity and smaller in models with higher viscosity and higher disk temperature, suggesting that gravitational instability may be a dominant cause of the observed diversity in M_{dot} in embedded YSOs. Our modeling predicts a weak dependence between the mean mass accretion rates and stellar masses in the Class I stage, in sharp contrast to the corresponding steep dependence for evolved T Tauri stars and brown dwarfs.