FDEM Simulation on the Failure Behavior of Historic Masonry Heritages Subjected to Differential Settlement

Abstract

Historic masonry heritages, such as cathedrals, colonnades, and arch bridges, were constructed with individual components (e.g., stones, bricks, other materials) bound together with, e.g., mortar, and they are very vulnerable to foundation settlement, especially differential settlement which occurs frequently in engineering practice. These masonry structures are discontinuous, and therefore, their behavior under differential settlement is highly nonlinear and complex. In this study, the combined finite-discrete element method (FDEM) is employed to simulate the failure behavior of historic masonry heritages subjected to support differential settlement. In the FDEM models, structures are discretized into elements where FE formulation is incorporated, resulting in an accurate estimate of structural deformation and interaction forces. In addition, a fracture model is employed for masonry blocks. Numerical examples are given and compared with results from the literature, showing that the FDEM is applicable and reliable in simulating the failure behavior of historic masonry heritages. Further analyses including block fracture reveal that fracturing can decrease the capacity against settlement significantly.

1. Introduction

Many historic heritages in the world are dry-joint masonry structures, e.g., Hadrian’s Wall in the UK and Temple of Olympian Zeus in Greece. They were made of individual bricks or stones, which were bound together by mortar. However, it is widely found that the mortar degraded severely or lost completely over ages, and even more, it was not used at all. In civil engineering, foundation settlement is common and has detrimental effects, especially for such dry-joint historic structures due to a lack of bonding strength between individual bricks or stones. Thus, it is of theoretical and practical significance to study the behavior, especially the potential damage or failure mode, for such precious heritage masonry structures under settlement.

In Heyman [1], leaning towers caused by uneven foundation settlement were investigated, and a rule was developed for the maximum inclination for the safety of masonry towers. Atamturktur et al. [2] detected structural damages caused by the settlement of buttresses in Beverly Cathedral, whose arch crowns and walls separated and deformed seriously. Milani et al. [3] performed case studies on three inclined masonry clock towers, and results showed that tilting significantly reduces the bearing capacity of masonry structures and increases seismic vulnerability. Drougkas et al. [4] investigated the cracking of the nave wall of St. Jacob’s Church by foundation settlement. Some further research on the influence of foundation differential settlement on structural integrity can be referred to [5,6,7,8,9,10]. The respective authors have demonstrated that differential settlement increases the probability of structural failure, especially when there is no mortar or the strength of mortar is low, and thus these dry joint assemblies are suitable for considering a discrete modeling approach to masonry.

Experimentation is a direct way of studying the responses of masonry structures by differential settlement. Giardina et al. [11] performed a 1/20 scaled model test on a stone masonry façade to evaluate its failure mechanism subjected to differential settlement. Portioli and Cascini [12] investigated the collapse of rectangular masonry walls under foundation differential settlement experimentally and identified the locking failure mode. Romano and Ochsendorf [13] studied and compared the mechanical behavior of various Gothic masonry arches under differential settlement through experiments, and the obtained results showed that pointed arches could withstand larger support displacement than circular arches.

Due to high costs in conducting physical tests, numerical simulation has become increasingly popular in recent decades. The finite element method (FEM) has been employed in evaluating the behavior of masonry structures caused by differential settlement. Alessandri et al. [14] analyzed the cracking of a masonry façade under foundation differential settlement based on 2D homogenized nonlinear finite element model. Landolfo et al. [15] employed similar homogenized FE model and predicted the failure modes of two-story masonry wall façades induced by differential settlement. Truong-Hong and Laefer [16] devised FE models and investigated the influence of window shape and size, block orientation, and lintels on the failure of masonry walls caused by excavation subsidence. Malena et al. [17] combined the piecewise rigid displacement (PRD) method and the FE approach and studied the failure of masonry arch bridge caused by pier displacement.

On the other hand, the discrete element method (DEM) is a discontinuous computational technique for analysing the responses of masonry structures. Bui et al. [18] employed the DEM in their research and simulated the in-plane and out-of-plane behavior of dry-jointed masonry walls under support differential settlement. Baraldi et al. [19] developed a DEM model and evaluated the nonlinear behavior of masonry panels with regular textures subjected to in-plane forces. Sarhosis et al. [20] studied the structural behavior of a two-story colonnade under static and dynamic loads by using the commercial DEM code-UDEC and identified the major factors affecting the stability of colonnades. Foti et al. [21] simulated the collapse of masonry cross vaults induced by support displacement by using a commercial DEM software-3DEC, and compared their results with those from tests.

Besides FEM and DEM, other numerical methods, such as limit analysis and the PRD method, have also been used to study the failure behavior of masonry structures under differential settlement. Gagliardo et al. [22] investigated the failure mechanism of a masonry church façade under support differential settlement based on the rigid block limit analysis. Iannuzzo et al. [23] addressed the stability of 2D masonry structures under large support displacement with the PRD method.

Lately, the combined finite-discrete element method (FDEM) is an advanced numerical approach developed by Munjiza in 1990s [24]. In FDEM models, structures are fully discretized into number of elements, and FE formulation is incorporated within these elements, enabling accurate estimate on structural deformation and interaction forces. The details about the FDEM can be found in the literatures [25,26,27]. Recently, the FDEM has been used to simulate the failure of brittle/quasi-brittle solids under static/quasi-static and dynamic loads [28,29,30,31,32,33]. Regarding the failure of masonry structures, Chen et al. [34] employed the FDEM to investigate the collapse of dry-jointed masonry arches subjected to support movement, and investigated the effects of geometry and friction coefficient. Chen et al. [35] simulated the behavior of masonry walls subjected to support differential settlement with the FDEM, considering floor load and block fracture. Pepe et al. [36] investigated the effect of geometry, opening and region of differential settlement on the failure behavior of masonry structures with the FDEM. Smoljanović et al. [37,38,39,40] also used the FDEM to model and analyze masonry structures.

In this study, the behavior of historic masonry structures due to foundation differential settlement was simulated with the 2D FDEM program “Y” developed by Munjiza [41]. It was designed to demonstrate some concepts explained in [25]. In order to use the FDEM program “Y”, an input file describing the investigated structure needs to be prepared firstly. A flowchart is shown in Figure 1. Relevant data is updated in the database after each time step, and the database is accessible for the computation of next time step.

The aim of this paper is to show the capability of FDEM on modeling historic masonry heritages subjected to differential settlement. This paper intends to provide a new and effective tool to predict the potential damage or failure mode of historic masonry heritages under various differential settlement scenarios, which is highly beneficial to protect these precious historic heritage structures against settlement risks and also to enrich the literature on FDEM masonry applications. The present study did not require Heyman’s hypothesis, which assumes zero tensile strength, infinite compressive strength, and no sliding of blocks [42], and furthermore, the influence of block fracture was taken into account. The simulation results were compared and validated with the data from the literature, providing insights into the protection of these historic masonry structures. Layout for the rest of the paper is as follows. Section 2 introduces the methodology, particularly the fundamentals of the FDEM. Section 3 presents numerical examples on three types of historic masonry structures (i.e., Natività della Beata Vergine Maria church in Bondeno, Pompeii colonnade, and Spanish Deba arch bridge) subjected to foundation differential settlement. In Section 4, the fracture behavior of masonry units, which was ignored in Section 3, is included, and the simulation results are compared with the counterpart results in Section 3. Finally, concluding remarks are summarized in Section 5.

2. Methodology

In the 2D FDEM program “Y”, simple constant strain triangular elements are used to reduce the complexity in contact detection and interaction. These elements are linear elastic and cannot be further divided. Separation between any pair of adjacent elements is achieved through the element interface, characterizing the failure criterion of masonry structures. In the following, the methodology, including motion of element, contact and fracture model, is introduced with details.

2.1. Equation of Motion

According to Newton’s second law of motion, the governing equation of an element i is given by Equation (1):

$$m_i{ (\mathbf{a}-\mathbf{g})=\mathbf{F}}_i$$

(1)

where m_i is the element mass; g is the gravitational acceleration; a is the acceleration except g; F_i is the net external force. Based on Equation (1), the position and velocity of element i can be obtained for each time step.

2.2. Contact

In the FDEM, contact is categorized as contact detection and contact interaction. The Munjiza-NBS contact detection algorithm [25,43] is employed. Contact interaction follows contact detection and determines contact forces between detected contact couples.

Contact forces are determined using the penalty function method [44]. Figure 2 shows two entities in contact in 2D with an overlapping area S. Penetration area dA yields the infinitesimal contact force df, as shown in Equation (2):

$$d\mathbf{f}=-{d\mathbf{f}}_t{+d\mathbf{f}}_c$$

(2)

where subscripts ‘t’ and ‘c’ relate to the target and the contactor, respectively. df_t and df_c are given in Equations (3) and (4), as:

$${d\mathbf{f}}_t{=-E}_p{\mathit{grad}\phi }_c{(\mathrm{P}}_c)dA$$

(3)

$${d\mathbf{f}}_c{=-E}_p{\mathit{grad}\phi }_t{(\mathrm{P}}_t)dA$$

(4)

where P_c and P_t are the points sharing the same coordinate on S, and they belong to the contactor and the target, respectively; φ_c and φ_t are pre-defined potentials; grad is the abbreviation of gradient; E_p is the penalty term. The total contact force f is given by Equation (5).

$${\mathbf{f}=E}_p{\mathbf{\int }}_{\mathrm{S}}{[\mathit{grad}\phi }_c{(\mathrm{P}}_c{)-\mathit{grad}\phi }_t{(\mathrm{P}}_t\left)\right]dA$$

(5)

2.3. Fracture Model

In the 2D FDEM, cracks are assumed to commence and develop only along element boundaries [45]. Interfaces are defined between each pair of unseparated elements. The breakage of two unseparated elements is determined by the interface deformation δ (Figure 3). When δ reaches the elastic limit δ_p, the damage initiates. When δ equals to the ultimate value δ_c, two elements separate completely. Although cracks can only initiate and develop along the element boundaries, mesh dependency can be reduced when unstructured fine mesh is used [46].

The fracture model employed in this study is similar to Hillerborg’s [47] and in a cohesive fashion [48]. In the FDEM, two unseparated elements initiate separating when σ = f_t, and at this time, δ = δ_p. Afterwards, σ decreases with increasing δ. When δ = δ_c, σ = 0 and the two elements are completely separated. The process that σ decreases with increasing δ is named with ‘strain softening’, which is characterized using a descending σ–δ curve. Crack onset and crack development are determined by actual stress and structural deformation since no pre-existing flaws or notches are required. The σ–δ relation is given in Equation (6) [25]:

$$\sigma =\{\begin{array}{cc}\left[2\frac{\delta }{{\delta }_p}-{\left(\frac{\delta }{{\delta }_p}\right)}^2\right]f_t & 0 \le \delta \le { \delta }_P\\ f_{t }z & {\delta }_P< \delta \le { \delta }_c\end{array}$$

(6)

The expression of z is given by Equation (7) [49]:

$$Z=\left[1+{\left(c_1{D}_t\right)}^3\right]e^{-c_2{D}_t}-D_t\left(1+c_1^3\right)e^{-c_2}$$

(7)

where D_t is the fracture index in tension within [0, 1], and has the form of

$$D_t=\{\begin{array}{ll}\left(\delta -{\delta }_p\right)/\left({\delta }_c-{\delta }_p\right)& {\delta }_p< \delta \le {\delta }_c\\ 1& \delta > {\delta }_c\end{array}$$

(8)

In Equation (7), c₁ and c₂ are constant parameters. For masonry structures, c₁ = 3.0 and c₂ = 6.93 were suggested by Zivaljic et al. [50]. Thus, the z–D curve is plotted in Figure 4, and the area beneath the curve is the Mode I strain energy release rate G_f, i.e.,

$$G_f={\int }_{{\delta }_p}^{{\delta }_c}\sigma d\delta$$

(9)

Sliding (denoted by s) may occur when two adjacent elements are subjected to shear deformation. When s reaches the elastic limit s_p, correspondingly, the shear stress τ reaches the ultimate strength f_s. Similar to the tension, τ decreases with the increasing s, and τ drops to 0 when s reaches the ultimate sliding s_t. Thus, the definition of τ is given in Equation (10),

$${\tau =(1 - D}_s{) f}_s$$

(10)

where D_s is the damage index about shear failure within [0,1], and it is given as:

$$D_s=\{\begin{array}{ll}\left(s-s_p\right)/\left(s_t-s_p\right)& s_p<\left|s\right|\le s_t\\ 1& \left|s\right| > s_t\end{array}$$

(11)

To include both tensile and shear damages, a composite damage index D is defined and given in Equation (12):

$$D=\sqrt{D_t^2+D_s^2}$$

(12)

The above is the failure criterion used in the FDEM analysis of dry-joint masonry structures.

3. Numerical Examples

Historic masonry structures under differential settlement are simulated with the FDEM. The church façade in Bondeno, Pompeii colonnade and the Spanish Deba arch bridge are selected in this study, representing various types of historic heritage structures. The church façade is composed of masonry bricks with openings for doors and windows of circular or pointed shapes. Such façades are frequently found in the world, especially in Europe. The Pompeii colonnade is a representative of stone-stack structures, with small-number but large-size stone blocks forming columns and epistyles. Similar structures include Parthenon Temple in Greece, Stonehenge in the UK, etc., and finally, the Spanish Deba arch bridge represents wide-existing masonry arch bridges, and many of them are still in service today. These masonry structures were built hundreds or even thousands of years ago, and have deteriorated over the ages. In this regard, selected cases are highly worthy of being investigated and provide a good overlook towards existing historic masonry heritages under support differential settlement. Fracture is not considered in this section, and the failure of these historic heritages is attributed to the excessive deformation of the structure or instability induced between masonry blocks.

3.1. Church Façade

Gagliardo et al. [22] studied the response of a masonry façade in Nativita della Beata Vergine Maria church (see Figure 5a) in Bondeno, Italy under differential settlements using rigid block limit analysis. The church façade is 22.0 m wide and 19.0 m high, consisting of 1873 masonry blocks. The FDEM mesh is shown in Figure 5b. The right support is assumed to settle downwards at a velocity of 0.05 m/s, while the left support is fixed. In total, there are 3890 elements in the model, and the material parameters in FDEM analysis are given in

Table 1. Material properties of the church façade.

Young’s Modulus (GPa)	Poisson’s Ratio	Unit Weight (kN/m3)
10.0	0.2	18.0

based on [22]. The friction coefficient is 0.6, and the time step is 5 × 10⁻⁷ s.

Figure 6 shows the failure mode from the FDEM simulation for the differential settlement d = 0.2 m, with the result from limit analysis [22] for comparison. It is observed that cracks are mainly found in the middle and the right part of the church façade where the foundation settles. Further, cracks are mostly oblique, connecting the middle/right circular window or the door and the far field. The FDEM result is in good agreement with that from limit analysis [22]. The failure process of the church façade is presented in Figure 7. When d = 0.15 m, the main cracks in the middle and the right part of the façade appeared; when d = 0.245 m, the top circular window collapsed partially due to excessive deformation; when d = 0.48 m, the windows on the right collapsed, and the right part of the façade tended to rotate clockwise; and finally, when d = 0.685 m, the right part of the façade collapsed significantly while the left part was still largely stable.

By simulating the failure behavior of masonry church façade under differential settlement, the FDEM not only reproduced the crack mode of the church under foundation settlement, but also predicted the collapse process and the ultimate displacement, providing guidance on maintenance and reinforcement of such historic heritages.

3.2. Pompeii Colonnade

The remain structure of Forum in Pompeii (see Figure 8a), which is a group of two-story colonnades, is investigated. The geometry and the FDEM mesh are shown in Figure 8b. There are seven columns, denoted by UC1, UC2, UC3 in the top story and LC1, LC2, LC3, LC4 in the bottom story. Four supports are denoted by S1, S2, S3, and S4. Two reference points, i.e., RP1 and RP2, are selected for discussions. There are 94 elements in total, and the material properties are given in

Table 2. Material properties of the colonnade.

Young’s Modulus (GPa)	Poisson’s Ratio	Density (kg/m3)
40.0	0.25	2680.0

according to [20]. Each lower column is placed above a movable support so that different differential settlement combinations can be simulated. The settling rate of supports is assumed to be 0.05 m/s. The friction coefficient is 0.6, and the time step is 5 × 10⁻⁷ s.

Several differential settlement scenarios are investigated: (1) one single support settles, (2) two supports settle simultaneously, and (3) three supports settle simultaneously. In total, there are 14 cases (see

Table 3. Cases of support settlements.

Case No.	Case 1	Case 2	Case 3	Case 4	Case 5	Case 6	Case 7
Settled support	S1	S2	S3	S4	S1 & S2	S1 & S3	S1 & S4
Case No.	Case 8	Case 9	Case 10	Case 11	Case 12	Case 13	Case 14
Settled support	S2 & S3	S2 & S4	S3 & S4	S1 & S2 & S3	S1 & S2 & S4	S1 & S3 & S4	S2 & S3 & S4

). Failure mode and the corresponding ultimate support displacement for each case are presented in

Table 4. Failure of the colonnade due to differential settlement without considering fracture.

Case 1	Case 2	Case 3	Case 4
Case 5	Case 6	Case 7	Case 8
Case 9	Case 10	Case 11	Case 12
Case 13	Case 14

and Figure 9, respectively. The collapse of upper columns and the fall of blocks in epistyle are observed in most cases. Cases 1, 7, and 13 are the most critical cases in their own category, and the corresponding ultimate support settlements are 170.0 mm, 175.0 mm, and 140.0 mm, respectively. Among the 14 cases, Case 13 exhibits the lowest ultimate support displacement, i.e., the colonnade structure has the lowest capability against differential settlement. For Cases 5 and 10, the colonnade structure did not collapse at the end of the simulation, i.e., the ultimate settlement is larger than 700 mm.

The failure of the three most critical cases (i.e., Cases 1, 7, and 13) were all attributed to the falling of block B1 in the epistyle. When S1 settled in Case 1, the left columns (UC1 and LC1) tilted to the left and B1 fell, as shown in Table 4. Similarly, in Case 7, UC1 and LC1 inclined to the left and B1 fell down too. On the other hand, another group of columns (UC3, LC3, and LC4) leaned to the right, and B2 and B3 fell, due to the differential settlement of S4. The second columns on the left (UC2 and LC2) were almost straight. In Case 13, UC1, and LC1 still leaned to the left, while the right part of the colonnade was nearly in a vertical translation motion, and the moving inconsistence between S2, S3, and S4 leaded to the falling of B2 between them. The differential settlement of S3 and S4 accelerated the falling of B1, reducing the ultimate support displacement, and hence, the capability of Case 13 against differential settlement is smaller than that of Cases 1 and 7. The differential settlement of S1 is involved in all the three critical cases. Thus, S1 is the most critical member for the structure subjected to differential settlement. Additionally, the falling of B1 is determinant in the three most critical cases, and the differential settlement of other supports in Cases 7 and 13 accelerate this process.

3.3. Arch Bridge

In Malena et al. [17], the damage of the Deba arch bridge in Spain (Figure 10a) due to a pier settlement was studied. The geometry and the FDEM mesh of the four-span masonry arch bridge are shown in Figure 10b, and an uneven settlement is applied to pier P-2 with δ₁ = 0.90 m and δ₂ = 0.86 m, while other piers and the abutments are kept intact. There are 1528 elements in total. Material properties are given in

Table 5. Material properties of Deba arch bridge.

Young’s Modulus (GPa)	Poisson’s Ratio	Density (kg/m3)
3.18	0.23	2500.0

based on [17], and the friction coefficient is 0.6. The time step is 2.5 × 10⁻⁷ s.

Figure 11 presents the FDEM simulation result against that from the geometrical survey [17]. After pier P-2 settled with δ₁ = 0.90 m and δ₂ = 0.86 m, five hinges (i.e., A, B, C, D, and E) were found in the masonry arch bridge. Some differences are still observed, since the actual structure and the differential settlement are complicated, and the FDEM simulation is performed in 2D, not in 3D. However, in spite of these differences, the FDEM simulation results are highly consistent with results from the geometrical survey [17], especially on hinges and their locations and the crack development. Thus, the FDEM is reliable in simulating the damage of masonry arch bridges under pier differential settlements. If P-2 settled further proportionally, the deformation of the arch bridge increased and the hinge positions changed. The vertical displacements of the bridge at different settlement stages are presented in Figure 12. It can be observed that the bricks separated clearly (Figure 12b,c), hinge A and B moved up towards the crown, similarly for hinge D. Hinge E moved downwards to pier P-3, and hinge C disappeared.

4. Further Analyses

Units in large-size masonry structures may fracture after differential settlements. In this section, the fracture of masonry units is also included using the model addressed in Section 2.3, and the effect on the failure behavior of masonry structures is studied.

4.1. Fracture of Colonnade

The fracture of the Pompeii colonnade due to differential settlements is investigated in this section, and a refined mesh with a characteristic size of 0.1 m is used (see Figure 13). Additional material parameters for fracture are: tensile strength f_t = 2.0 MPa, shear strength f_s = 1.0 MPa and the strain energy release rate G_f = 100.0 J/m². Other material parameters and the time step are kept as same as in Section 3.2. In this case, δ_p = 4.0 × 10⁻⁶ m, δ_c = 2.0 × 10⁻⁴ m, s_p = 2.0 × 10⁻⁶ m, and s_t = 4.0 × 10⁻³ m.

Similar to the case definitions in Table 3, the 14 cases are named as Case 1′, Case 2′, …, and Case 14′, respectively.

Table 6. Failure of the colonnade due to differential settlement considering fracture.

Case 1′	Case 2′	Case 3′	Case 4′
Case 5′	Case 6′	Case 7′	Case 8′
Case 9′	Case 10′	Case 11′	Case 12′
Case 13′	Case 14′

shows the failure of the colonnade structure subjected to the combinations of support differential settlement considering fracture, and the corresponding ultimate support displacement for each case are presented in Figure 14. In addition to the collapse of upper columns and the fall of blocks in epistyle observed in non-fracture cases in Section 3.2, fracture of blocks was identified. The most critical cases were still Case 1′, 7′, and 13′ in their respective categories. However, the corresponding ultimate support displacements were much smaller (i.e., 146.3 mm, 155.0 mm, and 136.3 mm, respectively), suggesting that fracture makes the colonnade more vulnerable to differential settlements. It may be explained that fracture makes multiple micro-cracks possible in masonry blocks, which reduce the stability of the entire colonnade structure to some extent. Among the 14 cases considering fracture, Case 13′ is still the most critical one with the lowest capacity against differential settlement.

Comparing Table 4 and Table 6, it is worth noting that the failure mechanisms are largely the same, despite some upper-story columns or blocks in the epistyle fractured. Figure 15 shows the development of the vertical displacement of reference points RP1 and RP2. Obviously, block fracture results in earlier failures of the colonnade structure. In some cases, (Case 4′, 6′, 9′, 12′), the ultimate support settlement when considering fracture is slightly larger than that when fracture is not considered. One possible reason is that the movement of masonry blocks were stuck by micro-cracks and small fragments, leading to the sense that a larger ultimate support settlement was obtained. Fracture of blocks makes the failure of colonnade more complex and different failure mechanisms may occur. Fragments may prevent the structure from moving, i.e., increasing the stability to some extent, as was observed in Chen et al. [34]. FDEM modeling considering fracture of blocks revealed the influence of fracturing on the failure of colonnade-type structures and provided beneficial supplements to traditional rigid DEM modelling.

4.2. Fracture of Arch Bridge

In this section, the Deba bridge is investigated as a monolithic bridge, and a new fine mesh with a characteristic size of 0.4 m is created (see Figure 16a). Material parameters about fracture are: tensile strength f_t = 1.0 MPa, shear strength f_s = 1.0 MPa, and the strain energy release rate G_f = 25.0 J/m². Other material properties and the time step are kept as the same as in Section 3.3. In this case, δ_p = 8.0 × 10⁻⁶ m, δ_c = 1.0 × 10⁻⁴ m, s_p = 8.0 × 10⁻⁶ m, and s_t = 1.0 × 10⁻³ m. Figure 16b shows the damage of the monolithic arch bridge under uneven settlement of P-2. Comparing Figure 11b and Figure 16b, the cracking paths and cracking locations are similar to each other. The ultimate pier settlement of the monolithic arch bridge is quite small (i.e., δ₁ = 0.04 m, δ₂ = 0.038 m), suggesting that non-fractured masonry arch bridge can withstand a much larger pier settlement since it can undergo a much larger deformation than the monolithic one. The differences between monolithic and masonry modeling strategies of a same structure have been captured by the FDEM.

5. Concluding Remarks

The FDEM was employed to analyze the failure behavior of historic masonry heritages subjected to differential settlement. Simulation results were compared with those from the literature. In Section 3.1, crack paths and the failure process of an Italian church façade under support differential settlement were obtained. In Section 3.2, the failure and ultimate differential settlement of Pompeii colonnade were presented, and the most critical cases were identified. In Section 3.3, the hinge positions and failure mode of a four-span Spanish arch bridge subjected to an internal pier differential settlement were captured. The influence of fracture on structural responses was studied, showing that fracture will reduce the capacity against differential settlement significantly. Additionally, deformation capacity of masonry structures is better than that of the monolithic counterpart built with the same configuration. In general, the FDEM has proven to be a useful tool in analyzing the failure behavior of historic masonry heritages by differential settlement.