Abstract: In order to improve the status that the spherical harmonics method in the ANSYS-Fluent software is limited to its lowest order (P1) form for radiation heat transfer calculations, based on the basic principles of solving the radiative transfer equation using high-order spherical harmonics, by coupling the numerical solution of the governing equations and boundary conditions using the user-defined function (UDF) and user-defined scalar (UDS) interfaces provided in the ANSYS-Fluent software, a high-order spherical harmonics radiation transfer equation solution model suitable for the software was developed. The model accuracy was verified by calculating the radiative heat transfer in one-dimensional, two-dimensional flat plates, and two-dimensional axisymmetric flame, and comparing with the corresponding analytical or photon Monte Carlo (PMC) solutions. Results show that the P3 model developed in this paper is more accurate than the built-in P1 model of ANSYS-Fluent software in calculating both the angular distribution of radiation intensity and the overall radiative heat source.