A new level-set model is proposed for simulating immiscible thermocapillary flows with variable fluid-property ratios at dynamically deformable interfaces. The Navier–Stokes equations coupled with the energy conservation equation are solved by means of a finite-volume/level-set approach, adapted to a multiple marker methodology in order to avoid the numerical coalescence of the fluid particles. The temperature field is coupled to the surface tension through an equation of state. Some numerical examples including thermocapillary driven convection in two superimposed fluid layers, and thermocapillary motion of single and multiple fluid particles are computed using the present method. These results are compared against analytical solutions and numerical results from the literature as validations of the proposed model.