A consistent formulation for unilateral contact problems including frictional work hardening or softening is proposed. The approach is based on an augmented Lagrangian approach coupled to an implicit quasi-static Finite Element Method. Analogous to classical work hardening theory in elasto-plasticity, the frictional work is chosen as the internal variable for formulating the evolution of the friction convex. In order to facilitate the implementation of a wide range of phenomenological models, the friction coefficient is defined in a parametrised form in terms of Bernstein polynomials. Numerical simulation of a 3D deep-drawing operation demonstrates the performance of the methods for predicting frictional contact phenomena in the case of large sliding paths including high curvatures.