The numerical modeling enables us to reduce the risk related to the selection of best localization of wells. Moreover, at the stage of production, modeling is a suitable tool for optimization of well operational parameters, which guarantees the long life of doublets. The thorough selection of software together with relevant methodology applied to generation of numerical models significantly improve the quality of obtained results. In the following paper, we discuss the impact of density of calculation grid on the results of geothermal doublet simulation with the TOUGH2 code, which applies the finite-difference method. The study area is located between the Szczecin Trough and the Fore-sudetic Monocline, where the Choszczno IG-1 well has been completed. Our research was divided into the two stages. At the first stage, we examined the changes of density of polygon calculation grids used in computations of operational parameters of geothermal doublets. At the second stage, we analyzed the influence of distance between the production and the injection wells on variability in time of operational parameters. The results demonstrated that in both studied cases, the largest differences occurred in pressures measured in production and injection wells whereas the differences in temperatures were less pronounced.