The Monte Carlo (MC) method is a numerically efficient way to solve the Boltzmann Transport Equation (BTE), which governs the motion of electrons and holes in semiconductors. The method itself, now ubiquitous in computational science, was originally devised by Stanislaw Ulam to solve the BTE for transport of neutrons in the fissile material of the A-bomb. Since these pioneering times in the mid 1940's, the popularity of the MC method has grown with the increasing availability of faster and cheaper digital computers. Its application to our specific problem of electronic transport in semiconductors is first due to Kurosawa in 1966 (J. Phys. Soc. Jpn., Suppl. 21, p. 424 (1966)). Shortly afterwards the Malvern, UK, group (Boardman, Fawcett, Hilsum, Swain, among others) provided the first wide application of the method to the problem of the "Gunn effect" in GaAs. Applications to Si and Ge boomed in the 1970s, thanks to work performed at the University of Modena, Italy. The review articles by C. Jacoboni and L. Reggiani, Rev. Mod. Phys., vol. 55, No. 3, pp. 645-705, July 1983, and by Peter J. Price, Semiconductors and Semimetals, vol. 14, pp. 249-308 (1979) provide a deeper historical and technical perspective. The popularity of this method is testified to by the fact that, as of 1995, some 20 groups worldwide are active in MC simulation of electron transport in Si.

The main feature of the MC method is that, contrary to "drift-diffusion" or "hydrodynamic" models which treat electrons and holes as a fluid, the carriers are treated as point-like semiclassical particles (an assumption which eventually breaks down in yet smaller devices) moving within the device under the action of the external field and of the collision processes. Their individual trajectories are computed by following three simple steps:

1. By moving each carrier in "free flights" between collisions under the action of the electric field. This is done by solving Newton's equations of motion in the crystal. The band structure of the solid enters heavily in this kinematic step.
2. By letting the carriers scatter (with impurities, lattice vibrations, other carriers...) conditionally at the end of free flights. The probability of a collision is given by the total scattering rate (i.e., the number of collision per unit time).
3. Finally, by selecting a new "state" (i.e., the new energy, speed, and direction of motion) after a collision.

The name Monte Carlo derives from the algorithm used to determine when a carrier will collide at the end of a free flight, and to select its "final state" after the collision: The scattering probability and the probability distribution of the final states are computed using quantum mechanics. A particular event (collision or no collision, which type of collision, which final state) is selected randomly, by comparing the probability of occurrence of that event to a random number, as when throwing dice or playing roulette.

The figure below shows a typical electron trajectory, with collisions labeled by the type of scattering responsible for each collision. The asterisks denote scattering with the Si-Si$O$₂ interface, and labels "ta/te", "la/le", and "oa/oe" indicate absorption/emission of transverse, longitudinal acoustic or optical phonons.

As the trajectories are computed, estimators of interesting physical quantities (e.g., carrier velocity, energy, current density) are obtained.

A major advantage of the MC method is its ability to provide a solution of the BTE without making additional approximations on the basic physics ( band structure and scattering processes). Thus, it provides a perfect laboratory for computer experimentation on the basic properties of electronic transport. The accuracy in the physical description comes at a price of hunger for computing time and of a loss of numerical accuracy in the solution, in the form of stochastic noise. Standard algorithms are available for reducing this noise. DAMOCLES makes extensive use of these "variance reduction" techniques which allow the study of "rare events" (high-energy tails of the distribution function, low-density regions in the device, etc.)

DAMOCLES uses an "ensemble" MC algorithm, consisting of the calculation of trajectories of many (10,000 to 100,000) particles simultaneously. Their positions are recorded at the end of each free flight, in order to gather information also about the charge density within the device. This is used to solve numerically the Poisson equation, yielding the electric field used to move particles in free flights...and so on in a self-consistent loop.