The explicit version of the Material Point Method has been extended in order to model coupled hydromechanical saturated problems. MPM discretizes the continuum, which is considered as a saturated soilfluid mixture, by dividing it into particles or material points. The discrete movement equations are not solved at the material points. Instead a support mesh, built to cover the domain of the problem, is used. In this paper it is assumed that particles carry all the variables needed to represent the state of the continuum including the pore pressure as a variable associated with each particle. The particle pore pressure increment is calculated explicitly using the equation of fluid mass balance, from the particle volumetric deformation and the fluid velocity relative to the soil skeleton, at the particle location. The shape functions used for the mesh elements are usually the same bi-linear functions of the Finite Element Method and therefore the background mesh elements suffer the same drawbacks. These drawbacks include: volumetric locking for quasi-incompressible materials when four particles per cell are used, which is equivalent to four integration points in the finite element method, pressure instability for quasi-incompressible and low permeability materials and the generation of zero energy modes when one particle per cell is used, which corresponds to reduced integration in the finite element method. The MPM original version has also the disadvantage of generating "noise" in the solution when a particle pass from one cell to another. A simple procedure that can be used to reduce instabilities is to consider constant stress at each cell equal to the stress average of the particles which are in the cell at the instant k. In this case the internal forces are obtained in the same way as in the finite element method when one point of integration is used, using the gradient of the shape functions calculated in the cell center. In this work, to avoid volumetric locking and simultaneously achieve a stable behavior, internal forces and pressure increments at the nodes are calculated using the gradients calculated at the cell center. The procedure is completely explicit and has proved to be stable for the low permeability values used to model the foundation of Aznalcollar dam. The simulation of Aznalcollar dam progressive failure is presented as an example.