An efficient numerical algorithm for discrete dislocation dynamics simulations in two-dimensional, finite polygonal domains is presented. The algorithm is based on a complex boundary integral equation method. By use of the fast multipole method, linear complexity and storage requirement are achieved. This method has not, to the present author?s knowledge, previously been used in such simulations. Convergence studies show that the algorithm is accurate and numerically stable. Results from uniaxial load and bending moment load simulations at different loading rates are presented. The effect of finite size is studied. The results show that higher loading rate gives less yielding, and that a smaller specimen is harder than a larger one. This is in agreement with well-known results, and demonstrates that the dislocation dynamics model can describe important features of the physical problem. The cut-off velocity, that is the maximum velocity of the dislocations, is an important model parameter. In the present paper, it is shown that a four times higher cut-off velocity than was previously deemed sufficient is needed to obtain results independent of the cut-off velocity for the bending moment load simulations. ? 2002 Elsevier Science B.V. All rights reserved.