This paper presents a systematic method by which closure relations for the ensemble-averaged equations of disperse two-phase flows of solid spheres can be derived. The method relies on the direct numerical simulation of three flow situations: equal forces or couples applied to the spheres, and an imposed macroscopic shear flow. A crucial aspect of the work is that it focuses on systems that are spatially non-uniform on average. It is shown that, due to this feature, several new terms arise in the constitutive relations that would vanish for a uniform system. For example, while the usual effective viscosity is recovered in the closure of the stress tensor, it is found that other terms are also present, which confer a markedly non-Newtonian nature to the rheological constitutive equation. (C) 2001 Elsevier Science Ltd. All rights reserved.