Here, we address the turbulent dynamics of the gravitational instability in accretion disks, retaining both radiative cooling and irradiation. Due to radiative cooling, the disk is unstable for all values of the Toomre parameter, and an accurate estimate of the maximum growth rate is derived analytically. A detailed study of the turbulent spectra shows a rapid decay with an azimuthal wave number stronger than ky−3, whereas the spectrum is more broad in the radial direction and shows a scaling in the range kx−3 to kx−2. The radial component of the radial velocity profile consists of a superposition of shocks of different heights, and is similar to that found in Burgers’ turbulence. Assuming saturation occurs through nonlinear wave steepening leading to shock formation, we developed a mixing-length model in which the typical length scale is related to the average radial distance between shocks. Furthermore, since the numerical simulations show that linear drive is necessary in order to sustain turbulence, we used the growth rate of the most unstable mode to estimate the typical timescale. The mixing-length model that was obtained agrees well with numerical simulations. The model gives an analytic expression for the turbulent viscosity as a function of the Toomre parameter and cooling time. It predicts that relevant values of α = 10−3 can be obtained in disks that have a Toomre parameter as high as Q ≈ 10.