  Working within the framework of the minimal supergravity model with gauge
coupling unification and radiative electroweak symmetry breaking, we evaluate
the cosmological relic density from lightest neutralinos produced in the early
universe. Our numerical calculation is distinct in that it involves direct
evaluation of neutralino annihilation cross sections using helicity amplitude
techniques, and thus avoids the usual expansion as a power series in terms of
neutralino velocity. Thus, our calculation includes relativistic Boltzmann
averaging, neutralino annihilation threshold effects, and proper treatment of
integration over Breit-Wigner poles. We map out regions of parameter space that
give rise to interesting cosmological dark matter relic densities. We compare
these regions with recent calculations of the reach for supersymmetry by LEP2
and Tevatron Main Injector era experiments. The cosmologically favored regions
overlap considerably with the regions where large trilepton signals are
expected to occur at the Tevatron. The CERN LHC $pp$ collider can make a
thorough exploration of the cosmologically favored region via gluino and squark
searches. In addition, over most of the favored region, sleptons ought to be
light enough to be detectable at both LHC and at a $\sqrt{s}=500$ GeV $e^+e^-$
collider.
