This paper presents a three-dimensional mathematical model for steamflooding in heavy oil saturated naturally fractured carbonate reservoirs where the governing equations are written for the matrix and fracture media, separately. Uncertainties associated with the shape factor for the communication between the matrix and fracture was eliminated through setting a continuity boundary condition at the interface. Using this boundary condition, the solution method employed differs from the most of the modeling simulations reported in the literature. The Newton-Raphson approach was also used for solving mass and energy balance equations. The recovery factor and cumulative steam to oil ratio were obtained as a function of steam injection rate, temperature, and also characteristics of the fractured media such as matrix size and permeability. The numerical solution clearly shows that fractures play an important role in improving the conduction of heat into the matrix part It is also concluded that the matrix size and total permeability are the most important parameters affecting the dependent variables involved in steamflooding in this type of reservoir. The results obtained from the mathematical modeling were compared with the real data, exhibiting a reasonable agreement.