Gradient chromatography has been widely applied in analytical and preparative chromatography since it provides better peak shapes and the ability to elute analytes in a shorter time frame. Apart from changes in the composition of a mobile phase also alteration of process temperature can be applied during sep-aration procedures to improve efficiency. However, proper mathematical modeling of the gradient chro-matography and further correct prediction of solutes' retention behavior have become a serious challenge as it involves the need to develop computational procedures that accurately account for the time and spa-tial gradients of crucial parameters. In this work, a computational procedure including the equilibrium- dispersive two-dimensional mass transfer model, the two-dimensional (2D) heat transfer model together with Darcy's law and the continuity equation have been proposed. Additionally, the calculation procedure was simplified by replacing the 2D model with the one-dimensional (1D) mass transfer model in order to speed up the computations. Both proposed solutions were validated employing external experimental data of temperature gradient HPLC [1] as well as with predictions based on the linear elution strength (LES) model available therein. The proposed procedures made it possible to efficiently predict the concen-tration profiles with average relative errors of calculated retention times not exceeding 3.22%. Moreover, the effect of the axial dispersion coefficient determination method on the obtained peak shapes was ex-amined involving the Gunn, the Wen-Fan, and the Chung-Wen correlations, indicating that the latter pro-duces gives the most accurate results. Finally, the proposed mathematical procedures were tested under UHPLC conditions, and due to a significant difference in retention times found the 2D model is strongly advised. (c) 2021 Elsevier B.V. All rights reserved.