A three-dimensional reactive multi-phase particle-in-cell (MP-PIC) model is employed to investigate biomass combustion and gasification in fluidized bed furnaces. The MP-PIC model considered here is based on a coarse grain method (CGM) which clusters fuel and sand particles into parcels. CGM is computationally efficient, however, it can cause numerical instability if the clustered parcels are passing through small computational cells, resulting in over-loading of solid particles in the cells. To overcome this problem, in this study, a distribution kernel method (DKM) is proposed and implemented in an open-source CFD code, OpenFOAM. In DKM, a redistribution procedure is employed to spread the solid volume and source terms of the particles in the parcel to the domain in which the particles are clustered. The numerical stiffness problem caused by the CGM clustering can be remedied by this method. Validation of the model was performed using data from different lab-scale reactors. The model was shown to be able to capture the transient heat transfer process in a lab-scale bubbling fluidized bed reactor under varying fluidization velocities and loads of sand. Then, the model was used to study the combustion/gasification process in a bubbling fluidized bed reactor under varying ambient temperatures, equivalent air ratios, and steam-to-biomass ratios. The performance of DKM was shown to improve the accuracy and the robustness of the model. © 2021 The Author(s)