Abstract Shear zones and competent layers and boudins represent viscosity heterogeneities within the rock mass. Differences in viscosity impel differences in strain rates between the background material and the heterogeneities. In this work, we represent the viscosity heterogeneities as elliptical inclusions. We use the Kolosov‐Muskhelishvili equations for the incompressible viscous flow problem in and around elliptical inclusions. Systematic investigation of the stress equilibrium at the matrix–inclusion interface shows that the mean stress, equivalent to the total pressure, is not continuous across viscosity boundaries. The results predict that pressure and stress perturbations depend strongly on the orientation of the elliptical heterogeneity with respect to the far-field stresses. A viscosity ratio of 10 between the inclusion and the surrounding material is sufficient to produce pressure discontinuities approximately equal to the magnitude of the effective shear stress of the strongest rock under the considered physical conditions. Comparison of the analytical solutions with thermo-mechanical models confirms pressure incongruity and suggests that dynamic parameters such as pressure and effective shear stress vary spatially and temporally within deforming, two-viscosity rock systems. As a corollary, the dependence of metamorphic phase equilibria on thermodynamic pressure implies that shear zones, taken as weak inclusions, and boudins, taken as hard inclusions, may record non-lithostatic pressure during deformation.