  On the example of the Poynting-Thomson-Zener rheological model for solids,
which exhibits both dissipation and wave propagation - with nonlinear
dispersion relation -, we introduce and investigate a finite difference
numerical scheme. Our goal is to demonstrate its properties and to ease the
computations in later applications for continuum thermodynamical problems. The
key element is the positioning of the discretized quantities with shifts by
half space and time steps with respect to each other. The arrangement is chosen
according to the spacetime properties of the quantities and of the equations
governing them. Numerical stability, dissipative error and dispersive error are
analysed in detail. With the best settings found, the scheme is capable of
making precise and fast predictions. Finally, the proposed scheme is compared
to a commercial finite element software, COMSOL, which demonstrates essential
differences even on the simplest - elastic - level of modelling.
