We propose the construction of a liquid-vapor relaxation model which is able to capture metastable as well as saturation states of the nonisothermal van der Waals model. The Gibbs formalism is recast as a constrained optimization problem, from which we recover physical laws, such as the Gibbs phase rule. From these results we deduce a dynamical system whose long time equilibria are the physically stable states. This behaviour is illustrated by numerical simulations. A tentative coupling with fluid dynamics is then presented: the dynamical system is coupled to a set of compressible Euler equations supplemented by convection equations on volume, mass and energy fractions of one of the phases.