A mathematical model for simulation of refinery furnaces is proposed. It consists of two different submodels, one for the process side and another for the flue gas side. The process side is appropriately modeled as a plug flow due to the high velocity of the fluid inside the tubes. The flue gas side is composed by a radiative chamber and a convective section both connected by a shield tube zone. Both models are connected by the tube surface temperature. As the flue gas side model uses this temperature as input data, the process side model recalculates this temperature. The procedure is executed until certain tolerance is achieved. This mathematical model has proved to be a useful tool for furnace analysis and simulation.