Earthquake faults are the result of severe strain localization in rocks deep down in the earth’s crust. This localized deformation is controlled by the size of the microstructure and various Thermo-Hydro-Mechanical (THM) couplings, whose modeling is central for understanding earthquake nucleation and seismic energy release. First, we model this challenging system using the Cosserat theory and by considering large shear deformations during seismic slip. In particular, we start our analysis by using the normalized coupled system of partial differential equations that include the THM couplings for the case of a Cosserat continuum. We then perform a bifurcation analysis, which indicates that traveling shear bands are possible inside the fault gouge. Next, we derive our non linear mesh independent numerical results accounting for the influence of large displacements by using and Adaptive Lagrangian Eulerian (ALE) procedure. We introduce viscoplasticity in our numerical analyses, which leads naturally to a "rate and state" frictional phenomenology. Depending on the boundary conditions, our numerical results show (a) frictional restrengthening and (b) the emergence of traveling shear bands along the thickness of the fault, leading to oscillations in the fault’s frictional response. Existing numerical analyses and the established models of uniform shear if Lachenbruch and shear on a mathematical plane of Rice, don’t capture this behavior. Therefore, at a second step, we extend the classical model of thermal pressurization to incorporate different strain localization modes, temperature and pore fluid pressure boundary conditions. This extension leads to a Volterra integral equation, which is solved semi-analytically. Our numerical findings (Stathas & Stefanou, 2023) show a good agreement with recent experimental results, that insulate thermal pressurization from other weakening mechanisms. This shows the relevance of the chosen THM mechanisms and of our models for studying fault friction and for improving the current understanding of this complex phenomenon based on computational mechanics.