  We investigate the effect of three-body correlations on the phase behavior of
hard rectangle two-dimensional fluids. The third virial coefficient, $B_3$, is
incorporated via an equation of state that recovers scaled particle theory for
parallel hard rectangles. This coefficient, a functional of the orientational
distribution function, is calculated by Monte Carlo integration, using an
accurate parameterized distribution function, for various particle aspect
ratios in the range 1-25. A bifurcation analysis of the free energy calculated
from the obtained equation of state is applied to find the isotropic
(I)-uniaxial nematic (N$_u$) and isotropic-tetratic nematic (N$_t$) spinodals
and to study the order of these phase transitions. We find that the relative
stability of the N$_t$ phase with respect to the isotropic phase is enhanced by
the introduction of $B_3$. Finally, we have calculated the complete phase
diagram using a variational procedure and compared the results with those
obtained from scaled particle theory and with Monte Carlo simulations carried
out for hard rectangles with various aspect ratios. The predictions of our
proposed equation of state as regards the transition densities between the
isotropic and orientationally ordered phases for small aspect ratios are in
fair agreement with simulations. Also, the critical aspect ratio below which
the N$_t$ phase becomes stable is predicted to increase due to three-body
correlations, although the corresponding value is underestimated with respect
to simulation.
