In this article, we develop a local discontinuous Galerkin (LDG) discretization of the (non)-isothermal Navier–Stokes–Korteweg (NSK) equations in conservative form. These equations are used to model the dynamics of a compressible fluid exhibiting liquid–vapor phase transitions. The NSK-equations are closed with a Van der Waals equation of state and contain third order nonlinear derivative terms. These contributions frequently cause standard numerical methods to violate the energy dissipation relation and require additional stabilization terms to prevent numerical instabilities. In order to address these problems we first develop an LDG method for the isothermal NSK equations using discontinuous finite element spaces combined with a time-implicit Runge–Kutta integration method. Next, we extend the LDG discretization to the non-isothermal NSK equations. An important feature of the LDG discretizations presented in this article is that they are relatively simple, robust and do not require special regularization terms. Finally, computational experiments are provided to demonstrate the capabilities, accuracy and stability of the LDG discretizations.