We use a ghost-cell based high-order immersed boundary method (IBM) to study the thermal interaction between entrained solid spherical particles and a turbulent velocity- and temperature-carrying flow. The sphere diameter (D) is about eight times the Kolmogorov scale (eta). The ambient turbulent field is isotropic, and the Taylor microscale Reynolds number is 50. The inflow turbulent velocity intensity varies from 0.05 to 0.1, and the intensity of temperature fluctuation varies from 0.1 to 0.4. The particle volume fractions are 0.01 and 0.02, and particle-to-fluid density ratios are 1.2, 10.0, and 100.0. It is observed that the particle-to-fluid density ratio affects the Nusselt number the most, followed by the solid volume fraction and turbulence intensity, while the effect of the intensity of temperature fluctuation is relatively small. It is also shown that correlations that have been proven to be valid for a single stationary particle can deviate significantly from the exact value obtained by directly integrating the dimensionless temperature gradient over the surface of the particle. Better estimates can only be gained by also taking the local flow and thermal conditions into consideration, in addition to the particle Reynolds number. (C) 2017 Elsevier Ltd. All rights reserved.