A thermomechanical analysis of elasto-plastic bodies is a necessary step toward the understanding of tribological behaviors of machine components subjected to both mechanical loading and frictional heating. A three-dimensional thermoelastoplastic contact model for counterformal bodies has been developed, which takes into account steady state heat flux, temperature-dependent strain hardening behavior, and interaction of mechanical and thermal loads. The fast Fourier transform and conjugate gradient method are the underlying numerical algorithms used in this model. Sliding of a half-space over a stationary sphere is simulated with this model. The friction-induced heat is partitioned into two bodies based on surface temperature distributions. In the simulation, the sphere is considered to be fully thermoelastoplastic, while the half-space is treated to be thermoelastic. Simulation results include surface pressure, temperature rise, and subsurface stress and plastic strain fields. The paper also studies the influences of sliding speed and thermal softening on contact behaviors for sliding speed ranging three orders of magnitude.