Thermochemical conversion of larger biomass particles (thermally thick regime) toward high-end products still suffers from an unrevealed quantitative relationship between process and product parameters. The main issue relates to the influence of heating rate within the particle, critical conversion-wise but difficult to assess experimentally. Computational fluid dynamics (CFD) modelling may help, but first the model must prove its reliability to prevent error transfer to the results. This study aimed to provide an unbiased, state-of-the-art model constructed in a stepwise mode to investigate the heating rate's distribution. Several datasets with broadly varying parameters from the literature were used for the development and validation since the reproduction of datasets would not bring novelty to solving the problem. Instead of the model's calibration to fit to the data, the parameters for each step-model were meticulously selected to match the experimental conditions. The stepwise development showed the best accuracy when the anisotropy and the heat sink drying sub-model were implemented. Moreover, using the Ranzi-Anca-Couce (RAC) scheme led to more accurate results than the Ranzi scheme. The comprehensive model was positively validated against a broad range of production parameters (pyrolysis temperature: 500 °C - 840 °C, diameter of particles: 10 mm - 20 mm, shapes: cylinders and spheres). Investigation showed a pattern in volatiles release profiles and homogeneous heating rate distribution when particle size is below 4 mm. Despite basing the models on the literature's data, the study includes novel and valuable insights for biomass conversion and constitutes a solid foundation for future development.