A thermodynamic method is proposed to generate sequences of optimal non-ideal thermal energy recovery systems (TERS). The two-fold objective of using the thermal energy recovered from hot process streams primarily as heating power and then as shaft power is considered. Shaft power generation might be a technical goal by itself. However, it always gives a way of measuring the possibilities for the original set of streams to be efficiently integrated with additional process streams. If T v is the absolute temperature of the steam utility and T o is that of the cold utility, then the total power recovered is assessed as G = Q E( T v - T o)/( T v) + τ, where Q E is the recovered heating power and τ is the shaft power that the system is able to generate. Each solution in a sequence of optimal non-ideal TERS will show: (i) the maximum value of G for the total exchange area A required and (ii) the maximum value for Q E, provided that condition (i) has already been met. Therefore, generation of shaft power is subsidiary to thermal energy recovery as heating power. It is shown that for optimal networks the rate of internal generation of entropy, σ, attains its minimum value compatible with constraints describing: (i) the transformations to be operated on the process streams and (ii) the subsidiary character of shaft power generation to thermal energy recovery as heating power. A systematic procedure is deviced to incorporate these constraints to the functional form for σ making use of the operating line concept. Optimal non-ideal networks are derived from a family of operating lines that are extremals (minimals) for the functional σ subject to the prescribed constraints. A procedure to evolve from this type of optimal networks to new ones showing a maximum for the ratio β = Q E/ A for each given size of the exchange area, is also outlined. Applications to a test problem available in the literature and to an azeotropic distillation are also made.