Influence of soil–structure interaction on seismic demands of historic masonry structure of Kashan Grand Bazaar

Heritage structures are an important part of the history of each community as well as an economic resource. This study aims to investigate the earthquake damage and failure mechanism of a historic Bazaar located in Kashan (center of Iran), which dates back to the eighteenth century. Numerical models were built by a macro method, considering two cases of fixed support and soil–structure interaction (SSI). To assess the structural behavior, modal, nonlinear static, and dynamic time history analyses were performed using a three dimensional finite element method. The results of SSI analysis show an intense weakness of the structure against the required demand based on the design spectrum. Also, foundation flexibility has a great impact on the vibration characteristics, displacement demand, and damage distribution. By comparing the results of time history analysis to the pushover procedure, a lower shear value and higher displacement are obtained. Furthermore, the damage process is as dome collapse alongside short arch, consequently long arch collapse followed by pier wall overturning.


Introduction
Historical buildings represent a set of specific heritage, economic and social values, which can be considered as parts of identity and continuity (Feilden 2007). Nowadays, the determination of international institutions in charge of recommendations for historical monuments is to protect the integrity and importance of such buildings Machat and Ziesemer 2017). Most heritage buildings are made of masonry that naturally has low capacity against lateral loads. In addition, seismic design concepts were usually not considered during the construction of historic structures. Thus, damage evaluation and vulnerability analysis of such structures are very crucial, especially in earthquake-prone areas. Due to the complications in geometry, boundary constraints and material characteristics of historic masonry constructions, their vulnerability assessments are often carried out using numerical tools (Peña et al. 2010;Valente and Milani 2016a;Roque et al. 2019;Aghabeigi et al. 2021;Aghabeigi and Farahmand-Tabar 2021). Page (1978) developed a finite element model for masonry including material nonlinearity and failure mechanisms of the joint. Lourenço et al. (1995) used contact element to micro-model the in-plane unreinforced masonry walls. In this method, the mortar joints and the masonry elements are defined separately as a continuum. Regarding troubles in understanding geometry, modeling, and the need for excessive computational power, the micromodel is appropriate for small-scale simulations (Lourenço 1997). Therefore, homogenization techniques are addressed, in which it is not possible to distinguish masonry units and mortar joints from each other, and masonry constructions can be modeled in an integrated manner (Marques and Lourenço 2014;Zizi et al. 2022). The simplified macro modeling homogenized strategy is very promising in analyzing large-volume models and is applicable for the safety assessment of masonry constructions in engineering practice due to the adequacy and computational efficiency.
Soil-structure interaction (SSI) is a procedure in which soil response influences structure movement and vice-versa (Kramer 1996). Depending on the structural characteristics, the features of the subsoil, and the frequency content of input motion, ignoring the SSI effects can be beneficial or detrimental (Tahghighi and Rabiee 2017;Tahghighi and Mohamadi 2020). In historical masonry buildings due to the limited height and stiff sections, the SSI could have a significant effect, especially for structures on soft soils (Piro et al. 2020). In engineering practice, two common approaches for subsoil modeling are numerical continuum methods and Winkler-based methods. Despite the simplicity and practicality of Winkler's hypothesis, it is basically unable to describe the actual characteristics of the soil medium due to the independent assumption of overlapping soil layers, the main focus on the adjacent soil touching the foundation, and the difficult estimation of spring stiffness values (Tahghighi and Konagai 2007;Far 2019). The numerical method, however, is the most prominent and rigorous solution to deal with complicated geometries and subsoil conditions, especially for 3D nonlinear problems (Tabatabaiefar et al. 2015;Sharari et al. 2022). Housner (1957) was one of the pioneers who studied the interaction of ground and buildings subjected to an earthquake. Maria et al. (2015) used experimental results from the Ghirlandina monitoring tower station to investigate the accuracy of simple methods to consider the SSI effects. Using the Winkler method, de Silva et al. (2015) studied the SSI impacts on the seismic behavior of the monumental heritage towers. The effects of SSI on the historical masonry bridge using the continuum finite element method was assessed by Güllü and Jaf (2016). Fathi et al. (2020a, B) studied the seismic behavior of the historical citadel of Tabriz and concluded that considering the SSI increases the probability of overturning, decreases the acceleration response and affects stress and deformation distribution. The influence of soil characteristics on the seismic damage of historical masonry minaret using the continuum model and the semi-infinite boundary element was investigated by Hökelekli and Al-Helwani (2020). de Silva (2020) evaluated the SSI effects on the input motion measure and the required demand of masonry towers against weak to strong earthquakes using the 3D finite difference method, where the foundation flexibility causes the structural demand enhancement and considerable tower tilt.
In this research, a part of Kashan Bazaar with minimal intervention and damage to the historic fabric of the structure was surveyed and the 3D geometric model was made based on the macro method (Systemes 2019). Finite element Abaqus code was applied to perform numerical simulations (Systemes 2014). Modal, nonlinear static, and time history analyses have been employed to evaluate the structural performance of Bazaar by considering two cases of fixed support and SSI. The results of the analysis show that SSI has a great influence on the dynamic properties, pushover capacity, failure mechanism, base shear, and displacement demands.

Kashan Grand Bazaar
In ancient Persian language, the word Bazaar comes in the form of "Vazar" and "Vakar", meaning community and place for trade (Karimi et al. 2015). First, Bazaar was formed seasonally and then weekly near large villages, and over time it changed from temporary to permanent constructions. Kashan Bazaar is one of the most famous and glorious historic Bazaars and the largest extensive monument in the center of Iran (Fig. 1). It is noteworthy that Kashan Grand Bazaar, as a national heritage with the registration number 1284 in August 1976, is in preparation for the UNESCO world heritage list. The history of Kashan Bazaar predates the tenth century. The destructive earthquakes of Kashan in 1755 and 1778 caused serious damages to the Bazaar which, was reconstructed over the following 5-6 years (Ambraseys and Melville 1982).
As a masonry construction, Kashan Bazaar with a linear direction consists of a series of interrelated buildings and consecutive domes and arches with a width of 3-5 m as illustrated in Fig. 2. The complex structure of the Bazaar includes covered passages with abundant markets on both sides; the primary path is named "Rasteh". Domes of the main Rasteh are from arch types, which is one type of dome and arch. Dome and arch are the outstanding Iranian methods of vault construction (Golabchi and Javani Dizaji 2018). Arches that are perpendicular to the Bazaar's direction are called long arches and arches which are in the Bazaar direction are called short arches. Figure 3 gives the typical architectural parts of the Bazaar.

The finite element model
The historical masonry structure under study and its foundation is composed of clay bricks and lime-cement mortar. Due to the lack of a reliable experimental result, the mechanical properties of the masonry material were evaluated considering the poor level of information in technical standards (NTC 2008) and the test results on similar historical buildings in Iran (Aghabeigi and Farahmand-Tabar 2021;Fathi et al. 2020b), see Table 1. Masonry has complex behavior and the implementation of an effective and simple model for this material is under investigation (Degli Abbati et al. 2019). The specificity of each masonry unit in geometry and type of material is one of the principal limitations in this way. Accordingly, the use of isotropic models based on damage plasticity is generally accepted for masonry (Valente and Milani 2016b;Lee and Fenves 1998;Lubliner et al. 1989;Hillerborg et al. 1976;Karimi et al. 2016Karimi et al. , 2017. Figure 4 shows the masonry stress-strain relationship and the damage evolution. In this model, when the masonry is unloaded at any point in the softening branch of the stress-strain curve, its unloaded response and modulus of elasticity of the material are weakened. The decrease in modulus of elasticity is indicated by two damage parameters, d c and d t . The amount of the damage parameters varies from zero, which represents the material without damage, to one, which represents a complete failure, and it is calculated according to the following equations:   where σ c = compressive stress, σ t = tensile stress, E 0 = preliminary elasticity modulus, ε t = tensile strain, ε t pl = plastic tensile strain, ε c = compressive strain and ε c pl = plastic compressive strain. The dilation angle controls the quantity of plastic volumetric strain made in the course of plastic cutting. It is accepted that the value of dilation angle is steady during plastic rupture and its value for masonry is adopted equal to 10° (Valente and Milani 2016b). Eccentricity is a parameter associated with the slope of the flow potential function. If zero is taken, the noted function is converted to a line; it also makes difficulties in the convergence of the problem with very small amounts. σ b0/ σ c0 is biaxial to uniaxial compressive strength. A value between 1 and 1.27 is recommended for this parameter (Karimi et al. 2017). Selecting less than 1 caused a divergence, but by selecting more than 1.27 no change was found. Viscosity is a parameter that enables the preservation of the solving problem after breaking and decreasing stiffness. K c is the second stress invariant ratio and a value between 0.5 and 1 is proposed for this parameter (Karimi et al. 2017). Table 2 summarizes the parameters adopted for the damage plasticity model in the numerical simulations. Table 3 lists the mechanical properties of the subsoil layers (Municipality of Kashan 2010). The Mohr-Coulomb model is taken into consideration as elasto-plastic behavior for soil (Coulomb 1776). This model is widely used in geotechnical studies, and its accuracy for sandy soils (Bazaar site soil) has been confirmed in comparison with laboratory results (Jin et al. 2019;Yang et al. 2012;Silva et al. 2018;Ghadimi Chermahini and Tahghighi 2019). Within this material model, the fracture is controlled by the maximum shear stress that depends on normal stress. Geometry plays a key role in numerical modeling. In order to be more precise, two survey methods have been used in this study. First, images are taken from different distances of the selected arches, then the images are analyzed utilizing the ImageJ computer program (Schneider et al. 2012) and the size of various segments is gotten. To make certain the results are accurate, the size of various segments was then collected using point transfer and triangulation methods (Salih et al. 2012). It is noted that the collected values in both methods had a good overlap with each other. From geotechnical investigations, the subsoil is taken into consideration as five distinct layers and bedrock at a depth of 20 m (Table 3). Based on previous studies (Ghosh and Wilson 1969), the length and width of the soil profile have been taken 81.5 m and 56 m, respectively. Also, 10 m was considered for both the length and width of infinite boundary elements of CIN3D8 (Systemes 2014). Infinite elements were used to model the semi-infinite environment of the subsoil.

Fig. 4 Constitutive laws for masonry
In order to reduce the size of the complex three dimensional model and consequently the computational efforts, only a four-span segment of the Rasteh Bazaar was selected for numerical simulation (Fig. 5). As shown, the 3D finite element models were constructed under two different cases of SSI and fixed support in the ABAQUS software. The soil-foundation interaction was described with the aid of contact elements utilizing normal and tangential behavior. The friction coefficient was set to 0.577 (Municipality of Kashan 2010). Stress continuity and displacement compatibility on the boundary among soil layers were taken using 1161 Tie elements (Systemes 2014). 2415 and 28,259 four-node tetrahedral element (C3D4) has been employed to mesh the foundation and structure above. Subsoil meshing was performed using 9364 eight-node brick elements (C3D8) as well. The finite and infinite elements used in this study are shown in Fig. 6. For the sake of implementing proper boundary conditions, the movement of the foot of the piers was closed in translation directions for the fixed support model. For the SSI case, the horizontal sides of the subsoil have been blocked perpendicular to the plane for pushover analysis. The freedom degrees of the bedrock were blocked in all three translational directions for both models, with and without SSI.

Validation experiment
The primary step in the finite element method is to guarantee that the numerical model works correctly. In order to validate the numerical method, the finite element results are compared with laboratory tests on masonry walls by Karimi et al. (2016). The masonry wall includes solid brick and mortar via one over one volumetric ratio. The damage plasticity theory was employed for the nonlinear modeling of masonry. Configuration, boundary constraints, material characteristics, and applied forces have been considered in line with the experiment. Table 4 gives the characteristics of the masonry wall. Figure 7 depicts the masonry wall numerical model. A ten-node brick element (C3D10) was used for the wall  meshing and a four-node 3D rigid element (R3D4) was used for meshing U-shaped steel profiles. Figure 8 compares the numerical and laboratory hysteresis curves. Furthermore, Fig. 9 compares the envelope of numerical and laboratory analyses with the capacity curve from the pushover procedure. It is observed that the simulation results agree well with the laboratory measurements.

Analyses performed
Eigenvalue frequency analysis is frequently used for controlling accuracy and obtaining the vibrating characteristics of the model. For an initial assessment of the structure, modal analysis was performed on two models; fixed support and SSI cases. Because of the high capability of the Lanczos procedure in analyzing complex problems (Systemes 2014), this technique is utilized for eigenproblem solution in the present study to calculate natural frequencies and mode shapes. Pushover is a suitable and compatible procedure for the analysis of masonry structures. The analysis result is the capacity curve wherein the vertical axis describes the base shear in terms of top lateral displacement on the horizontal axis. Control points have been taken based on Fig. 10 To perform this analysis, first, the gravitational forces were taken into consideration, and then the horizontal forces were applied independently on the horizontal X and Y axes of Envelope curves of experiment and numerical analyses results the structure. According to NTC regulations (NTC 2008), the lateral load distribution is described using two configurations of G1 (linear distribution of acceleration at height) and G2 (steady distribution of acceleration at height). G2 and G1 configurations are usually critical in low-rise and tall structures, respectively (Castellazzi et al. 2017;Acito et al. 2014Acito et al. , 2016. Thus, for the Bazaar model, horizontal load distribution was assumed based on the G2 configuration. The capacity curve stops on the correspondent movement of the base shear, identical to the peak value of 85%. The implicit integration solver has been used to overcome the numerical convergence complexities. The vulnerability of the Bazaar structure is investigated by method N2 (Fajfar 1999). In this method, respectively, the multiple-degree-of-freedom system is transformed into an equivalent system with one-degree-of-freedom, the capacity curves are converted to bilinear load-displacement diagrams, and the elastic demand spectra are transformed into the inelastic demand spectra. Eventually, the inelastic demand spectra and the capacity curves are united in Acceleration Displacement Response Spectra (ADRS). The required elastic demand spectra are extracted according to Standard 2800 (2015) for two earthquake hazard levels; i.e., Design Basis Earthquake (DBE) and Maximum Credible Earthquake (MCE) with return periods of 475 and 2475 years, respectively.
Following the previously described analysis procedures, nonlinear dynamic time history analyses were carried out. Due to the lack of strong earthquake records in Kashan, three consistent record sets were carefully selected using the PEER-NGA database (PEER 2020) representing possible seismic scenarios on the bedrock. These records were chosen by considering the features of historical earthquakes in the studied region, such as fault rupture mechanism, magnitude and distance (Tehran Padir Consultant Engineering 2012). The details of selected records are given in Table 5. The vertical component of records was neglected. Earthquake records have been scaled to the reference peak ground acceleration (PGA) of 0.27 g for the Bazaar site from seismic hazard assessment. Figure 11 presents the seismic hazard map of Kashan in terms of PGA on the seismic bedrock. It must be noted that earthquake waves are affected by local site conditions as they pass through the soil layers from the bedrock toward the free-field surface. Hence, to perform time history analysis for the fixed support model, site response analysis was performed using the scaled bedrock records and the soil layer properties in Table 3 (Hashash et al. 2020). Figure 12 shows the response spectra for 5% damping of the longitudinal (Z-direction) and transverse (X-direction) components of the scaled acceleration records on the free-field and the bedrock levels. Two orthogonal components of horizontal ground motions were applied to the model once, one on the X-axis and the other on the Z-axis, and once again vice versa. The results for the critical cases are then considered from each record set. To reduce the computation time for the time history analyses, the significant duration of earthquake records has been defined as between 5 and 95% of the Arias intensity (Trifunac and Brady 1975). 1 3 6 Results and discussion Table 6 gives the vibration period of the structure for the ten primary modes in both SSI and fixed support cases. Regarding Table 6, the mass participation factors differ in the same modes with the effects of SSI. It is observed that the third, fourth, and fifth mode shapes along the X-axis carry the largest participation factors and have the most influence on fixed support responses, while the most effective modes for SSI model are the third, fourth, and eighth modes. In addition, the most important modes for the fixed support model along the Z-axis are the first and seventh modes; while the most effective modes in the SSI model are the first and ninth modes. To verify the results, the relation of the fundamental period according to the NTC regulation, which is described as T 1 = 0.0187 h, has been used. By placing the height of the structure, h, which is equal to 5.12 m, in the mentioned equation; 0.096 s is obtained, which is very close to the fixed support first period from numerical analysis, which is 0.093 s. In order to perform a preliminary study on the effect of finite element mesh size on response accuracy, frequency analysis was carried out on a large number of fixed support models with different element sizes. Figure 13 illustrates a sensitivity analysis concerning the first three vibration periods and the element number. The investigation shows that convergence is achieved once the number of finite elements is increased. Figure 14 shows the first five shape modes for the given finite element models. It can be seen that the important modes in terms of displacement and rotation considering the SSI are different from the fixed support case. By comparing the first 60 modes of both models, only two common modes were found. Thus, it is concluded that SSI affects the shapes of the modes. Table 7 shows a comparison between common modes. For the Fig. 12 Acceleration response spectra of scaled records (bedrock and free field): a X-direction, and b Y-direction  X-direction, it is observed that the fifth mode of the fixed support is similar to the third mode of the SSI model, and in Z-direction, the first mode of both models (fixed support and SSI) is similar. According to Table 7, the vibration period of the modes related to the model considering SSI is much longer than the fixed support model. A period ratio is usually provided by seismic codes to evaluate the response of soil-structure systems. According to Table 7, the period ratio obtained from numerical simulations for the first two common modes in the directions X and Z are respectively, 3.26 and 2.32. On the other hand, the ASCE 7 standard (2016) provides the period ratio as follows: where T = period of the SSI model, T = fixed support period, K * fixed = stiffness of the fixed support model, h = structure height, k x = foundation transition stiffness and k θ = foundation rotation stiffness. The transitional and rotational stiffness of the foundation was calculated using guidelines given in NIST (2012). The soil properties are according to the geotechnical report in Table 3. For this purpose, the equivalent modulus of elasticity of layered subsoil was estimated as 34.32 MPa using the thickness weight averaging of the elasticity modulus of each soil layer to a depth of 20 m (Bowles 1996;Som and Das 2003). Through field survey and personal communication with relevant experts in the field, the dimensions of the foundation have been taken 14 m × 9 m × 2 m (see Fig. 5a). Finally, the ratio of the period of the SSI model to the fixed support model in the directions X and Z equals 3.27 and 2.86. It is observed that the period ratio derived by ASCE 7 is in good accordance with the results of numerical analysis.
A comparison between the capacity curves of the structure was performed in both situations, considering the SSI and the fixed support models. Figure 15 shows the capacity curves from the pushover analysis. On each response curve, three different points A, B, and C are marked; structure yielding occurrence (A); the appearance of maximum tensile stress at the heel of the masonry wall, ε tu = 0.15%, which is accompanied by an obvious change in the diagram slope (B), and the peak displacement value (C). As shown in Fig. 15, the SSI effects produce a reduction of the base shear and initial stiffness and increase the displacement capacity. The effects of SSI in the X-direction reduce 78.0% and 94.4% of the maximum base shear and stiffness, respectively, as well as increase 3.66 and 5.17 times the displacement at the yield and collapse points. In addition, for direction Z, SSI reduces 49.7% and 72.3% of the maximum base shear and stiffness, respectively, as well as increases 1.49 and 2.11 times the corresponding displacements. As expected, the soil compliance has a higher impact within the X-direction due to the greater stiffness provided by the piers' presence. It is observed that by increasing the rigidity of structure in the X-direction, the difference between pushover curves in the fixed support and the SSI models will be increased. The disparity in pushover curves with increasing structure rigidity indicates that the structure-foundation system is getting softer due to the capacity mobilization of a larger portion of supporting soil. Figure 16 shows the progress of damage by increasing load in three stages corresponding to points A, B, and C on the pushover curves. According to Fig. 16, the damage onset is accompanied by tension cracks in the footwall joint. Then, in the X-direction, by creating cracks on the right and left sides of the short arch and dome, the dome and the short arch collapse and finally the whole structure is destroyed by the overturning of the pier walls. In the Z-direction, with the growth of cracks in the area of the bearing walls-long arch-short arch-dome and their connection to each other, as well as the cracks in the connection of the short arch into the wall, the entire roof of the Bazaar collapsed and the entire structure is destroyed by the overturning of the masonry walls.
In the fixed support model, it is observed that the damage is more distributed in the entire structure because of the rigid base condition. The results of analyses indicate that plastic points and damage distribution are fundamentally modified by the SSI effects. As a consequence, the flexible foundation substantially affects the tension strain arrangement and afterward the local failure mechanisms. The numerical results agree generally with the observations arising from the field surveys which have recently been carried out for seismic assessment of Kashan historical Bazaar (Lazizi and Tahghighi 2019).  Figure 17 compares the capacity spectrum and the demand spectrum within the ADRS style. As shown, in design-basis earthquake demand, the fixed support model continues linearly in the X-direction whereas it hardly withstands collapsing in the Z-direction. However, the SSI model collapses completely under the DBE demand. Moreover, the structure does not have the required capacity for MCE demand at all. Comparing the results of the fixed support model with the SSI case shows that soil deformability leads to increment displacement, reduce base shear and factor of safety. Thus, the above-mentioned numerical simulations illustrate that the SSI effects are damaging to the structure seismic behavior, and the SSI omission overestimates the structure capacity and induces irrational responses.
For the given ground motion records, Fig. 18 shows the distribution of tensile damage in two models; fixed support and SSI cases. According to this figure, the state of the structure in the fixed support model is better than in the SSI model. Due to the Duzce and Tottori records, the structure maintains its stability and it is destroyed under the Manjil record. As indicated by response spectra in Fig. 12, the extent of damage can be attributed to the amplification and significant frequency content of the Manjil free field record, particularly in the low period range (less than 0.5 s), which corresponds to the fundamental period range of the Bazaar structure. It is observed that the state of the structure in the SSI model is critical so that the structure goes under complete collapse under all ground motion records. It is noteworthy to mention that due to the significant difference between fixed support and SSI fundamental frequencies, a different trend is observed for the SSI response in the same earthquakes. This time history results can be considered closely matched with the results of the N2 method, which confirm each other. Examining the form of damage, it can be concluded that the destruction progress of the structure under all three records is as follows: The collapse of the dome with a short arch and then the collapse of the long arch, and finally, the overthrow of the upright walls. The agreement between the time history and pushover procedures is verified by the tensile damage distributions obtained at the end of the numerical simulations, reported in Fig. 19. The damage pattern corresponding to the nonlinear dynamic simulations is similar to the one observed in the pushover approach, but severe damage is mainly concentrated on the upper parts of the structure along with the masonry dome and arches. In addition, localized damage at the base of the piers is observed during the nonlinear static analysis. Note that, damage in time history analysis is a collection of damage in both directions X and Z. Since the structure is subjected to earthquake record components in both directions, this result is rational. Generally, it can be expressed that the structure collapse is not only due to the growth of cracks in each area but also due to the joining of the cracks, which causes the failure of the whole structure. As mentioned before, this is because the continuous and interconnected performance of the force transition system has been a key part of the stability of the Bazaar.
The base shear using time history analysis for the selected seismic record set is shown in Fig. 20. According to this Figure, the foundation flexibility reduces the base shear compared to the corresponding fixed support model. The value of this decrease for the Duzce record in the directions X and Z is respectively, − 4.4% and 44.8%; for the Manjil record in the directions X and Z are respectively, 48.8% and 65.2%; and for the Tottori record in Fig. 17 Demand versus capacity spectra in directions X and Z for models with and without SSI under DBE and MCE hazard levels the directions X and Z are respectively, 79.5% and 67.0%. The maximum base shear from the pushover analysis (PA) for each direction is also shown as a horizontal line in the diagrams. It is observed that the maximum base shear time history is less than the maximum base shear from the pushover analyses, except in the SSI model where, depending on the ground motion record, the base shear time history may go more than the maximum base shear in the pushover analysis (Fig. 20a). Figure 21 compares the displacement time histories due to the considered ground motion records. As shown, the displacement demand in the SSI model has increased compared to the fixed support model, in contrast to the base shear demand (Fig. 20).  Fig. 18 Tensile damage distribution due to earthquake records for fixed support and SSI models (Top view and Bottom view) This increases for the Duzce ground motion record in the directions X and Z is respectively, 9.43 and 5.28 times. For the Manjil record in the directions X and Z is respectively, 8.19 and 6.57 times. For the Tottori record in the directions X and Z is respectively, 3.87 and 3.05 times. According to Fig. 21, the displacement time histories in the fixed support model are less than the PA procedure. However, for the SSI model, the displacement values have increased to such an extent that in some places it reached the displacements from the pushover analysis. In addition, the displacement in the case  Fig. 19 Comparison between the damage distributions for fixed support and SSI models: a pushover procedure, and b nonlinear time history analysis of fixed support along X-axis is much less than the Z-axis; but, this difference is not observed in the model considering SSI. It is noted that such pattern was also observed in the pushover capacity curves as discussed earlier.
In the case of a flexible foundation, a portion of the superstructure's vibrating energy is transmitted to the subsoil region, and it can be dissipated due to radiation and hysteretic damping. Evaluation of dissipated energy by the SSI model will be helpful to be referred to as a major reason for reducing base shear and increasing displacement demands of the structure. Figure 22 shows the time history diagrams of total energy dissipation by fixed support and SSI models. It is observed that foundation flexibility significantly intensifies the energy dissipated within the soil-structure system subjected to the given earthquake ground motion records.

Conclusions
In the present study, the seismic behavior and failure mechanism of the historical structure of Kashan Grand Bazaar was assessed by considering the effects of SSI. Modal, nonlinear time history and pushover analyses were carried out using detailed 3D finite element models. The masonry and soil medium have been considered by the damage plasticity and the Mohr-Coulomb material models, respectively. The validity of the numerical method was evaluated by comparison with experimental results. The main conclusions are as follows: • Soil compliance significantly affects the vibration periods and mode shapes, where the fundamental period of the SSI model was amplified up to 2.32 times compared to the fixed support case. • For DBE demand, the fixed support model sustains stability in the X-direction whereas it hardly withstands against collapse in the Z-direction. The SSI model collapses completely under the DBE demand. Also, both structure models do not have the required capacity for MCE demand. • The failure mechanism and damage distribution were essentially influenced by the foundation flexibility. The collapse of the structure is not only due to the growth of cracks in each component, such as dome, short arches, long arches, and piers, but also due to the joining of the cracks in the force transition system. • From pushover analysis, considering SSI caused a decrease of 70% in the ultimate base shear and an increase of up to 5.17 times in the displacement demands. • From time history analyses, considering SSI reduced the base shear by as much as 80% and increased the displacement up to 9.43 times. Thus, in time history analysis, less base shear and higher displacement values were obtained rather than through pushover analysis. • The results obtained in this study showed that SSI was detrimental to the structure behavior against earthquake loading, and its omission overestimates the structure capacity and induces irrational responses.
As described, the obtained results of this study are useful to evaluate the SSI effects on the seismic behavior of heritage structures. However, the employed numerical model still needs verification against additional experimental studies such as the mechanical properties of masonry and field pushover tests to generalize the findings for retrofitting recommendations.