In this paper, we present a solution to the contact problem for two elastic bodies with friction loaded by an arbitrarily varying force in 2D. The general statement of the problem concerns the following two important challenges. On one hand, the consideration of an arbitrary contact geometry implies that a multitude of individual arbitrarily-shaped contact zones appears, which can split and merge as a result of changes in the normal loading. On the other hand, the existence of friction invokes the presence of a tangential force and, consequently of a complicated traction distribution inside each of the contact zones. The complexity of these traction distributions results from the fact that, generally, each contact spot contains stick and slip zones with various boundary conditions defined by the previous tangential contact interaction forces. The presence of two classes of boundary conditions makes the problem memory-dependent and hysteretic. In spite of this, by using both conventional and original methods of contact mechanics, it is possible to drastically simplify the problem. First of all, the geometric difficulty residing in the description of arbitrary multiple contacts can be successfully dealt with by applying the reduced friction principle that expresses the tangential solution through the normal loading curve, for both normal and tangential forces kept constant. The principle is valid for a broad range of geometries and makes it possible to replace an arbitrary contact geometry by an equivalent system of two axisymmetric bodies. Secondly, the hysteretic character of the equivalent axisymmetric problem can be addressed through the method of memory diagrams which replaces the complex traction distribution in the contact zone by a simple internal functional dependency comprising the same memory information. Following any change of loading, the internal “memory function” evolves in accordance to a definite procedure that is based on the force balance equation written for small increments of the normal and tangential forces. In this paper, we formulate the method, prove its validity, and illustrate its numerical implementation with some simulated examples of hysteretic load–displacement curves. These obtained results might be of use for applications related to solids with cracks and contacts showing acoustical nonlinearity, and for any modeling requiring physics-based boundary conditions to describe and predict the effect of mechanical contacts with friction.