This paper is concerned with the energy consistent simulation of motions of a viscoelastic continuum body, under inclusion of the coupling of thermal and mechanical fields. The corresponding algorithm is based on a four-field formulation in the Lagrangian description, in which the deformation mapping, the velocity field, the temperature field and a strain-like viscous internal variable field are independent unknowns. Hence, the equations of motion are formulated in first-order form. The Lagrangian temperature field is determined by the first-order entropy evolution equation, associated with Fourier's law of heat conduction. The first-order viscous evolution equation is derived from an internal dissipation being quadratic in a nonlinear viscous strain-rate tensor. This coupled system of nonlinear differential equations is discretised by a new space-time finite element method, consisting of continuous as well as discontinuous finite element approximations in time. Owing to particular time approximations in the constitutive laws, beside the total linear momentum as well as the total angular momentum balance, a nonlinear stability estimate with respect to a relative energy function is exactly fulfilled in the fully discrete case as well.
Hence, the resulting time integration algorithm is long-time nonlinear stable also when changing the time step size. The obtained coupled system of nonlinear algebraic equations is solved by a monolithic solution strategy. The corresponding Newton-Raphson methods on the global and the element level are based on a consistent linearisation. The new convergence criteria used for these iterative solution procedures take the energy consistency into account, and is free of the scaling in the independent variables. Representative numerical simulations with various boundary conditions show the higher-order accuracy and the superior stability of the new time integration algorithm.