  We present and analyze a discontinuous Galerkin method for the numerical
modelling of the non-linear fully-coupled thermo-poroelastic problem. For the
spatial discretization, we design a high-order discontinuous Galerkin method on
polygonal and polyhedral grids based on a novel four-field formulation of the
problem. To handle the non-linear convective transport term in the energy
conservation equation we adopt a fixed-point linearization strategy. We perform
a robust stability analysis for the linearized semi-discrete problem under mild
requirements on the problem data. A priori hp-version error estimates in
suitable energy norms are also derived. A complete set of numerical simulations
is presented in order to validate the theoretical analysis, to inspect
numerically the robustness properties, and to test the capability of the
proposed method in a practical scenario inspired by a geothermal problem.
