A numerical two-phase one-dimensional mathematical model of a single channel catalytic monolith reactor was developed in this study. The model is capable of simulating a final steady state starting from an initial distribution of temperatures and concentrations. An algorithm was written in order to solve the conservation equations. In a catalytic reactor these equations are the mass and energy balances, rate equations and physical property relationships. The new method used to solve the equations was the combination of the Euler's method and central finite difference method. The advantage of this combined method is that the percentage error produced by the program code was negligible. The program code in MATLAB environment was executed to describe the profiles of temperatures and concentrations of each individual species in the gas and solid phase along the reactor. The reforming process is complicated for heavy hydrocarbons and it is not well defined. The model was applied, as far as the authors are aware for the first time, for autothermal reforming process with n-hexadecane (C16H34) feed to depict the performance of the monolith channel for hydrogen production at specified operating conditions. In this model, it was also assumed that the reactions only occur on the surface of catalyst. The temperatures and the concentrations of each species of interests (i.e. H-2, CO2, CO, and C16H34) were obtained in both the gas phase and at the wall. The results obtained using reaction kinetics data from literature were successfully validated against experimental results from a different study with respect to reactor temperatures and concentrations of products H-2, CO2 and CO. The influence of thermal conductivity and effective wall thickness of the solid phase, and mass flow rate of gas on wall temperatures and mole fractions of components were also investigated and the results were in agreement with literature data. Copyright (C) 2013, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.