A Monte Carlo approach is presented for the treatment of the pairing force in nuclear systems. This method is computationally efficient and very simple to implement. Numerical results are given and compared to an exact calculation and to the predictions of Bardeen-Cooper-Schrieffer (BCS) theory. (author) 21 refs.; 1 fig.; 2 tabs.