The terrestrial heat flux density, an essential information to evaluate the deep geothermal resource potential, is rarely defined over urban areas where energy needs are important. In an effort to fill this gap, the subsurface thermal conductivity estimated during two thermal response tests was coupled with undisturbed temperature profile measurements conducted in the same boreholes to infer terrestrial heat flow near the surface. The undisturbed temperature profiles were reproduced with an inverse numerical model of conductive heat transfer, where the optimization of the model bottom boundary condition allows determining the near-surface heat flow. The inverse numerical simulation approach was previously validated by optimizing a steady-state and synthetic temperature profile calculated with Fourier's Law. Data from two thermal response tests in ground heat exchangers of one hundred meters depth were analyzed with inverse numerical simulations provided as examples for the town of Québec City, Canada, and Orléans, France. The temperature profiles measured at the sites and corrected according to the paleoclimate effects of the quaternary glaciations were reproduced with the model. The approach presented offers an alternative to assess heat flow in the preliminary exploration of deep geothermal resources of urban areas, where thermal response tests may be common while deep wells are sparsely distributed over the area to assess heat flow.