## 2 Numerical Formulation  
The model used in this work is an extension of the plastic damage model developed by [PERSON] et al. [25] and presented in [PERSON] et al. [26].  
The model is based on continuum damage mechanics and uses a mixed implicit integration scheme.  
The adopted model is an orthotropic model of pure tension/compression damage (d\({}^{+}\)/d\({}^{-}\)) based on the continuous model to accurately reproduce the nonlinear shear response of masonry walls. Failure surfaces are defined through two scalar quantities calculated using Equations (1) and (2):  
\[\tau^{-}=H(-\overline{\sigma}_{\text{min}})\bigg{[}\frac{1}{1-\alpha}\bigg{(} \alpha\overline{I}_{1}+\sqrt{3\overline{J}_{2}}+k_{1}\beta\langle\overline{ \sigma}_{\text{max}}\rangle\bigg{)}\bigg{]} \tag{1}\]\[\tau^{+}=H(\overline{\sigma}_{\rm max})\bigg{[}\frac{1}{1-\alpha}\bigg{(}\alpha \overline{I}_{1}+\sqrt{3\overline{I}_{2}}+\beta\langle\overline{\sigma}_{\rm max }\rangle\bigg{)}\frac{\sigma_{t}}{\sigma_{p}}\bigg{]} \tag{2}\]  
where \(\tau^{-}\) represents the equivalent compression stress, \(\tau^{+}\) represents the equivalent tensile stress, \(\overline{\sigma}_{\rm max}\) represents the primary effective stress, \(\sigma_{t}\) represents the tensile strength of the units or mortar joints, \(\sigma_{p}\) represents the peak compressive strength of the units or mortar joints, \(I_{1}\) is the first invariant of the effective stress tensor, \(I_{2}\) is the second invariant of the effective deviatoric stress tensor, and \(k_{1}\) is the ratio between the biaxial and uniaxial compressive strengths.  
Equivalent stresses and damage evolution in tension and compression are defined by uniaxial stress-strain laws. The uniaxial stress laws are defined as shown in Figure 1 and are comprised of two parts: a linear part [(0; 0)--(\(\epsilon_{0}\); \(\sigma_{t}\))] with \(\sigma_{t}\) equal to the stress resistance of the unit or mortar joints and \(\epsilon_{0}\) equal to the corresponding deformation at \(\sigma_{t}\) and of a softening branch that depends on the tensile fracture energy \(G_{\rm t}\) and the characteristic length of the finite element \(\rm l_{dis}\).  
The uniaxial compression law is defined as in Figure 2 and is characterized by five parts: a linear part [(0; 0)--(\(\epsilon_{0}\); \(\sigma_{0}\))], a hardening section [(\(\epsilon_{0}\); \(\sigma_{0}\))--(\(\epsilon_{p}\); \(\sigma_{p}\))], two softening parts [(\(\epsilon_{\rm P}\); \(\sigma_{\rm P}\))--(\(\epsilon_{k}\); \(\sigma_{k}\))] and [(\(\epsilon_{k}\); \(\sigma_{k}\))--(\(\epsilon_{u}\); \(\sigma_{u}\))], and a final residual [(\(\epsilon_{u}\); \(\sigma_{u}\))--(+\(\infty\); \(\sigma_{u}\))].  
With:  
\(\sigma_{0}\) equal to the compression strength at the beginning of hardening;  
\(\epsilon_{0}\) equal to the deformation at the beginning of hardening;  
\(\sigma_{\rm P}\) equal to the peak compression strength of the units or mortar joints;  
Figure 1: Uniaxial stress laws.  
Figure 2: Uniaxial compression laws.  
\(\sigma_{\rm p}\) peak deformation;  
\(\sigma_{\rm r}=\sigma_{\rm u}\) equal to residual strength;  
\(\varepsilon_{\rm r}\) equal to the residual deformation;  
\(\sigma_{\rm u}\) equal to the ultimate deformation;  
\(\sigma_{\rm k}\) e \(\varepsilon_{\rm k}\) evaluated as intermediate control points.  
All the numerical simulations shown in this work are completed using the OpenSees [27] solver, where the authors have implemented the proposed constitutive model. Pre- and post-processing are carried out with STKO software [22].