Joint probability density function (PDF) calculations are reported of the bluff-body stabilized flames (HM1, HM2, and HM3) and the results are compared with the available experimental data. The calculations are based on the modeled transport equation for the joint PDF of velocity, turbulence frequency, and composition (species mass fractions and enthalpy) using the interaction by exchange with the mean and Euclidean minimum spanning tree mixing models. The methane chemistry is described by a 19-species augmented reduced mechanism, and is implemented using in situ adaptive tabulation. The numerical accuracy of the calculations is carefully studied, and the associated errors are quantified. For flame HM1 (which has the least local extinction), there is generally good agreement between calculations and measurements, although (for all flames) the quality of the agreement deteriorates at downstream locations. The calculations correctly show essentially inert mixing in the shear layer between the recirculation zone and the coflow in flame HM1, but not in flames HM2 and HM3. In general, the calculations of flames HM2 and HM3 are not in good agreement with the experimental data and do not exhibit the observed local extinction. This deficiency is attributed to the inaccurate calculations of the mean mixture fraction in the recirculation zone (for flames HM2 and HM3). The sensitivity of the calculation to the mixing model constant is investigated, and the mean scalar dissipation is reported.