Abstract. In this paper we propose a method for the construction of locally conservative flux fields from Generalized Multiscale Finite Element Method (GMsFEM) pressure solutions. The flux values are obtained from an element-based postprocessing procedure in which an independent set of $4 \times 4$ linear systems need to be solved. To test the performance of the method we consider two heterogeneous permeability coefficients and couple the resulting fluxes to a two-phase flow model. The increase in accuracy associated with the computation of the GMsFEM pressure solutions is inherited by the postprocessed flux fields and saturation solutions, and is closely correlated to the size of the reduced-order systems. In particular, the addition of more basis functions to the enriched coarse space yields solutions that more accurately capture the behavior of the fine scale model. A number of numerical examples are offered to validate the performance of the method.