A simple means of computing the rate of conformational space sampling and energy transfer in computer simulations of biomolecules using replica molecular dynamics is described. The method is based on the idea that in an ergodic system trajectories should be self-averaging-properties measured over two independent trajectories must average to the same result. Replica molecular dynamics simulation is used to calculate the generalized ergodic measure and the rate of self-averaging for the force and potential energy for the S-peptide and RNase A enzyme over a range of temperatures from 40 to 400 K. The results clearly demonstrate that even on a short time scale on the order of 10 ps, several distinct conformational states are sampled. The ergodic measures are used to obtain quantitative estimates of the rate at which conformational substates separated by relatively small barriers (on the order of a few kcal/mol) are sampled. Examination of the ergodic measure for nonbonded and dihedral angle forces proves that the time required for effective conformational space sampling is long (especially motions involving long length scales) compared to realizable computational times at all temperatures. The atomic force ergodic measure is evaluated for a harmonic system of normal modes and shown to provide a direct means of calculating the second moment of the vibrational density of states for the protein using a short dynamics trajectory. Finally, the instantaneous normal mode spectrum is calculated for the S-peptide as a function of temperature. A simple model of the potential energy hypersurface is developed and used to interpret the fraction of unstable modes in terms of the distribution of energy barriers separating the various peptide conformational substates. The distribution of energy barriers has a constant density of low-energy barriers and a Poisson distribution of high-energy barriers. The resulting energy barrier distribution is used to calculate the number of dihedral angle transitions expected in a dynamic trajectory, and the results are in good agreement with those found in the simulations. This study contains the first semianalytic method for extracting the distribution of barrier heights in systems with complex energy landscapes. The implications of our study for biomolecular simulations are discussed.