The behaviour of timber connections is very complex because of the interaction between brittle and ductile failure modes that develop within the contact areas between the timber and the dowel. Simulating this behaviour numerically requires the use of multidimensional failure criteria as well as the use of coupled constitutive equations accounting for both non-linear isotropic hardening and isotropic ductile damage. Such a formulation is proposed in the present work and implemented within a Vumat user-defined subroutine in Abaqus/Explicit. The model is used to simulate the behaviour of a double shear timber steel connection. Good agreement was found between the FE results and experimental ones, which showed the good capability of the model to predict the onset of ductile damage and growth. It was found that failure results from ductile defects initiation, growth and propagation inside narrow shear bands wherein the plastic strain is highly localised.