The huge number of elementary flux modes in genome-scale metabolic networks makes analysis based on elementary flux modes intrinsically difficult. However, it has been shown that the elementary flux modes with optimal yield often contain highly redundant information. The set of optimal-yield elementary flux modes can be compressed using modules. Up to now, this compression was only possible by first enumerating the whole set of all optimal-yield elementary flux modes. We present a direct method for computing modules of the thermodynamically constrained optimal flux space of a metabolic network. This method can be used to decompose the set of optimal-yield elementary flux modes in a modular way and to speed up their computation. In addition, it provides a new form of coupling information that is not obtained by classical flux coupling analysis. We illustrate our approach on a set of model organisms.