We study a classical fluid of nonspherical molecules. The components of the fluid are the ellipsoidal molecules interacting through the Gay–Berne potential model. A method is described, which allows the Percus–Yevick (PY) and hypernetted-chain (HNC) integral equation theories to be solved numerically for this fluid. Explicit results are given and comparisons are made with recent Monte Carlo (MC) simulations. It is found that, at lower cutoff l max , the HNC and the PY closures give significantly different results. The HNC and PY (approximately) theories, at higher cutoff l max , are superior in predicting the existence of the phase transition in a qualitative agreement with computer simulation.