Quenching is one of the most important heat-treatment processes for improving the mechanical properties of steel components in manufacture industry. The quenching stress is a source of cracking, which is frequently detrimental to steel properties. Therefore, the investigation of quenching stress is very important for the control of distortion, cracking and residual stress distributions of components. In the study of quenching stress, the measurement of stress distribution is necessary to the stress analysis and design of quenching process. However, in most cases, the cracking of a quenched component is caused by transient stress during quenching, while experiment can only measures the final internal stress (residual stress), rather than transient stress. As a result, the measurement of residual stress associated with finite element simulation (FES) has been a mainstream direction in the investigation of quenching stress. In this work, a full through-hardened 40CrNiMo cylinder with 60 mm diameter was water-quenched, and cooling curves at three positions along the radius of cylinder were measured. Then, an optimized heat transfer coefficient as a function of surface temperature was obtained by fitting with the measured cooling curves using the trial and error method. Based on an exponent-modified (Ex-Modified) normalized function describing transformation plasticity kinetics proposed, the thermo-elasto-plastic constitutive equations were deduced. The commercial finite element software, Abaqus/Standard, was used to solve the coupled temperature field, microstructure field and stress (strain) field. The results indicate that the quenching stress and its distribution predicted by FES is well consistent with those measured by XRD, which verified that the models employed in coupling of thermal field, phase transformation field and stress field including transformation plasticity function proposed are correct. Meanwhile, the features of residual stress distribution were revealed that compressive stress exists in the core and surface of cylinder and the maximum tensile stress exists at subsurface. The separated calculation of thermal stress and phase transformation stress by FES reveals the origin of residual stress distribution feature in quenched cylinders, that is, the relative higher phase transformation compressive stress and lower tensile thermal stress at the core of cylinder make the residual stress to be compressive, while at the surface of cylinder the compressive stress is predominantly from thermal stress, because it is much larger than the tensile stress caused by phase transformation stress. The tangential residual stress distributions in cylinders with several diameters from 3 mm to 100 mm were predicted by FES, and the results indicate that when diameter is less than 5 mm, the tensile stress at the surface increases with increasing diameter until to 5 mm, then decreases with increasing diameter to 20 mm, finally the tensile stress becomes compressive stress. Besides, with the increase of diameter, maximum tensile stress shifts from the surface to the location of 0.6 radius. The effects of different quenching media on quenching stress were also investigated by FES. The results demonstrated that although there is compressive stress at surface of cylinder quenched in water or salt solution, the maximum stress locates the subsurface, meaning that cracking easily occurs at the subsurface, which is consistent with cracking in practical components. This work is helpful for the analysis of cracking from quenching stress in components with different sizes and under different quenching media.