Granular materials are omnipresent in many fields ranging from civil engineering to food, mining and pharmaceutical industries. Often considered a fourth state of matter, they exhibit specific phenomena such as segregation, arching effects, pattern formation, etc. Due to its potential capability of realistically rendering these behaviors, the Distinct Element Method (DEM) is a very enticing simulation technique. Indeed it makes it possible to analyze and observe phenomena that are barely if at all accessible experimentally. DEM works by tracking every particle in the system individually, maintaining for each a trajectory influenced by external factors such as gravitation or contacts with boundary objects and by the interactions with other grains. The mathematical problem of identifying pairs of grains that interact and locating precisely where the contact occurs is highly dependent on the shape of the grains. We focus in this thesis on 3D spherical grains and use dynamic weighted Delaunay triangulations to track the collisions. The triangulation is built on the centers of the grains and evolves to follow their motion. We prove that all potentially colliding pairs of spheres are adjacent in the triangulation. As there are 6n to 8n edges for n spheres in most practical cases, the complexity of the collision detection becomes linear instead of quadratic in the number of particles, with a small overhead in maintaining the triangulation with efficient local operations. For the physical problem of realistically rendering the collision in a numerical contact model suitable for computer simulation, we have used widely accepted theories such as the viscoelastic model of Cundall, but have also tested some recent, more sophisticated developments in the field. The collision detection and contact models have been implemented in a modular DEM simulation code with advanced features in data structures storing the triangulation, in numerical robustness of the geometric computations, and in parallel processing on shared memory computers. Optimal packing of powders is important in many industrial processes, yet no theoretical result exists when dealing with grains of different sizes. We have performed simulations of such cases and could compare our results with experimental data. Preliminary results have been obtained regarding the relation between the size and proportion of grains and the density of the packing. Other simulations have also been performed, such as the granular flow through an hourglass. As no efficient simulation method is currently known for non-spherical 3D grains, we propose an intermediate approach of gluing spheres together into arbitrary shaped clusters and show some examples based on this approach.