The predictions from two formulations, the finite element method (FEM) and the fast Fourier transform-based (FFT) method, were compared in terms of the elastic field distribution in crystalline solids under thermal loading. Firstly, the thermoelastic responses of heterogeneous materials system from both methods are verified against the theoretical solution of Eshelby and Thin film case. Then, the spatial distributions of elastic fields of three-dimensional polycrystalline material under thermal loading were simulated. Both similarity and discrepancy of the responses from two methods were examined. In particular, the distributions of the extreme values in hot spots of stress, strain and elastic energy density (EED) were examined.