Large displacement analysis of dry-jointed masonry structures subjected to settlements using rigid block modelling

Engineering Structures 148 (2017) 485–496

Contents lists available at ScienceDirect

Engineering Structures journal homepage: www.elsevier.com/locate/engstruct

Large displacement analysis of dry-jointed masonry structures subjected to settlements using rigid block modelling Francesco Portioli ⇑, Lucrezia Cascini Department of Structures for Engineering and Architecture, University of Naples Federico II, via Forno Vecchio 36, 80134 Naples, Italy

a r t i c l e

i n f o

Article history: Received 21 March 2017 Revised 26 June 2017 Accepted 29 June 2017

Keywords: Masonry block structures Settlement Support collapse displacements Rigid block modelling Incremental limit equilibrium analysis Linear programming

a b s t r a c t A discrete element model is proposed for large displacement static analysis of dry-jointed masonry structures subjected to foundation settlements. The numerical model is a two-dimensional assemblage of rigid blocks interacting at no-tension, frictional contact interfaces with infinite compressive strength. The static contact problem associated to the proposed model is formulated as a linear mathematical programming problem and an incremental procedure is implemented to take into account large displacements. Experimental tests were carried out on small scale panels with opening subjected to imposed settlements to validate the proposed model. Comparisons in terms of failure modes and maximum allowable support displacements just before the collapse are presented. An experimental case study from the literature related to a masonry arch on spreading support is also analyzed for validation. On the basis of the obtained results, the accuracy and the computational efficiency of the formulation adopted are discussed. Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction This paper is devoted to the development of a discrete rigid block model for the assessment of dry-jointed masonry structures subjected to foundation settlements or spreading support in the large displacement regime. The proposed model takes into account both opening and sliding failure at block interfaces. The numerical formulation was developed in the framework of incremental limit equilibrium analysis and linear programming was used to solve the associated contact problem. As such, few geometrical and mechanical parameters govern the response of the numerical model, namely rigid block configuration, load distribution and friction coefficient, leading to efficient and robust computational solutions. Crack patterns and distortions induced by foundation settlement or spreading supports represent one of the most frequent causes of damage in masonry structures. Evaluating consequences associated to foundation or support movements, both in terms of damage (i.e. cracks width) and collapse (i.e. amount of support displacements involving loss of stability), is one of the main questions that, in the past decades, has interested architects and engineers dealing with the assessment of historical and existing masonry constructions [1]. Emblematic examples that posed significant challenges for builders are the differential foundation settlements of Venetian ⇑ Corresponding author. E-mail address: [email protected] (F. Portioli). http://dx.doi.org/10.1016/j.engstruct.2017.06.073 0141-0296/Ó 2017 Elsevier Ltd. All rights reserved.

buildings caused by soft soils [2], the settlement mechanisms in the naves of the Cathedral of Milan due to subsidence [3] and in the Cathedral of Agrigento due to slope instability problems [4], just to cite a few [5]. Other notable examples of failure mechanisms caused by support movement are frequently concerned with vaults and domes on spreading supports [6,7], which in these cases can be also related to inherent structural deficiencies in addition to soil deformation problems. From the computational point of view, the prediction of the effects induced by imposed movements on masonry structures can be made by using a variety of modelling approaches. Those include finite element [8–16] and distinct element models [17–20] as well as rigid body spring models [21–23]. As an alternative to the above-mentioned modelling methods, limit analysis (LA) can be also conveniently used to investigate the behavior of masonry structures subjected to settlements. When compared to other computational modelling approaches, the use of LA formulations is mostly attractive for its ability to assess failure modes and associated loads directly, with no iterative calculation steps and simple assumptions on mechanical properties of constituent materials [24]. Several analytical and computational models based on the upper and lower bound formulations of LA theorems have been proposed in the literature. Analytical formulations based on thrust line analysis methods as well as numerical software tools for the analysis of masonry arches subjected to spreading supports can

486

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

be found in [25–29]. Strut-and-ties models and load path method have been proposed in [30–32] for the interpretation of the crack patterns in the case of landslides-induced settlements. Discrete element models based on fundamental theorems of limit analysis, such as the rigid block models developed in [33– 41], can be also used to evaluate the crack patterns and to investigate the likely causes of damage observed on dry-jointed masonry structures subjected to settlements [42–44]. It is well known that, in the case of limit analysis formulations, the solutions obtained are expressed in terms of infinitesimal displacements. As such, when the objective of the analysis is to predict the amount of support movement promoting the collapse, a key issue is concerned with numerical techniques adopted to take into account large displacements. Generally, these numerical procedures are used to magnify infinitesimal displacements which are obtained from limit analysis formulations in order to determine new position and contact configuration of the blocks. In the framework of limit analysis formulations, different incremental solution procedures have been proposed in the literature to take into account the effects of large displacements [45]. Most of these procedures have been developed for numerical models derived from analytical formulations with reference to special structural configurations (i.e. thrust line analysis method applied to circular arches) and under simplified assumptions for failure mechanism (i.e. hinging or rocking with sliding prevented). Examples of incremental solution procedures based on static and kinematic approaches for the large displacement analysis of masonry arches subjected to spreading supports can be found in [25,26,29]. In this context, it is worth noting that incremental solution procedures have been also proposed in [46] for non-linear static analysis of masonry buildings using strut and tie models. In this case evolutive models are used for in-plane loaded wall panels to take into account interlocking and friction dissipation as well as large displacement effects. To overcome limitations related to structural configuration and failure mechanisms, rigid block models adopting incremental solution procedures could be also used, considering that in this case sliding and crushing can be included in the formulation of failure conditions in addition to hinging/rocking failure. Large movements effects using a computational rigid block model were addressed by Gilbert in [47] to model soil-structure interactions in the case of masonry arch bridges and tunnel linings. The Author proposed a coupled analysis procedure which involves updating the geometry and soil pressures at successive iterations. In the framework of rigid block models, an attractive and promising method to deal with configuration changes in gross displacement analysis makes use of gap functions for the formulation of non-penetration conditions at contact points. The use of such functions is well consolidated in non-smooth contact analysis [48–51], and recent applications to the case of masonry structures under seismic loads can be also found in [52]. However, to the best of the authors’ knowledge, no application exists of these formulations to the case of large displacement static analysis of masonry structures subjected to support movement. With that in mind, the present study is aimed at developing a simple rigid block model for large displacement analysis of dryjointed masonry structures. The proposed model relies on an incremental limit equilibrium analysis formulation using gap functions for the prediction of failure mechanisms and for determining the amount of support movement promoting the collapse. A key feature of the proposed formulation is related to the simplified procedures which have been adopted for the specification of potential contact points and for the updating of the gap functions to take into account gross-displacements and to analyze the response of large numbers of rigid blocks with contact and friction. As such, the formulation adopted represents an extension to the

case of large displacements of the computational model presented in [44], where infinitesimal displacements have been assumed to investigate the behavior of masonry block structures subjected to settlements. The objective of this study was to assess the ability of the proposed formulation in predicting the behavior of masonry structures subjected to foundation settlement by comparing numerical and experimental outcomes, and to evaluate the accuracy and computational efficiency of the adopted modelling approach when applied to large scale wall panels. The paper is organized as follows. The relationships governing the behavior of the rigid block model and the contacts formulation are presented in Sections 2 and 3. The adopted iterative solution procedure for large displacements is described in Section 4. In Section 5 a validation study is presented against ad-hoc experimental investigation involving small-scale masonry portal frames. Experimental results from the literature on a circular masonry arch subjected to spreading support are also used for validation. To assess the efficiency of the proposed model when applied to assemblages with a large number of blocks and contacts, in Section 6 the formulation is applied to numerical case studies of 2D wall panels with different load configurations. 2. The rigid block model for large displacement analysis To represent the masonry structure subjected to settlements or spreading supports, we assume an assemblage of rigid blocks i interacting at no-tension, frictional contact interfaces j with infinite compressive strength (Fig. 1). The rigid block behavior is governed by equilibrium equations, relating external and internal (contact) forces, and by kinematic equations, which ensure compatibility between contact displacement rates and block degrees of freedom. To model the interface behavior, failure criteria are defined on the basis of contact forces and flow rules are adopted to express the increment of displacement rates when failure occurs. In order to apply base movements to the masonry panels, an additional support block s is introduced in the numerical model (Fig. 1a). Although, in principle, different support displacement profiles could be modelled using the present formulation, in this study we assume that the settlement is governed by a single degree of freedom which is associated to this support block (i.e. vertical, horizontal or rotational movement). On the basis of the relations governing the rigid block assemblage, upper and lower bound problems are formulated and solved for the calculation of the failure mechanism and contact forces associated to the imposed movement. Large displacements are taken into account implementing a simple sequential solution procedure which is divided into increments to determine new positions of rigid blocks and contacts configuration, as detailed in the following sections. 2.1. Contact forces and kinematic variables A point contact model is adopted to represent interactions at interfaces j. The primary static variables are the internal forces acting at each contact point k, which are located at a vertex of interfaces j (Fig. 2a). These variables are collected in vector ck and include the shear force component tk and the normal force nk along the local coordinate axes:

ck ¼ ½ tk

nk T :

ð1Þ

Following usual limit analysis formulations, external loads f are expressed as the sum of known dead loads fD and live loads fL which act at the centroid of the blocks.

487

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

block ‘i’

block ‘i’

contact point ‘k’ k

j

k+1

j+1

k+3

interface ‘j’

fDs

k+2

block ‘i+1’

Fixed support

fLs=- fDs

Support block ‘s’ for selement

a)

b)

Fig. 1. (a) Rigid block assemblage; (b) Rigid block i, interface j and contact point k.

xi zi

fzi nk

z

f

i

unk

fxi

tk x

x0i

i

xi

utk a)

b)

Fig. 2. (a) Contact forces and (b) kinematic variables at block centroid i, and contact point k.

For rigid blocks i used to represent the brickwork subjected to settlement, the external loads are assumed to be equal to the known dead loads fDi:

f i ¼ f Di ¼ ½ f xi

f zi

T

f hi  :

ð2Þ

At the support block s the loads are work-conjugated to the degree of freedom of the moveable block and are expressed as the sum of the dead load fDs and live load fLs, increased by an unknown scalar multiplier a (Fig. 1a):

f s ¼ f Ds þ af Ls :

ð3Þ

The dead load fDs is assigned as an admissible base reaction (i.e., in the case of vertical movement, this is assumed to be equal to the resultant force calculated at the base along the foundation settlement considering a uniform distribution of vertical loads in the structure). The variable (live) load fLs is a force that is directed along the opposite direction and is expressed as the dead load fDs multiplied by the load factor a. As such, the following expression holds for the load applied at the support block:

Eq. (5) is the simple expression also adopted in [44] to introduce the foundation settlement in the model using dead and live loads at the support block as a function of the collapse load multiplier. The position at the centroid of block i is collected in the vector xi:

xi ¼ ½ xi

zi

hi T :

ð6Þ

As mentioned above, for the support block s a single degree of freedom is assumed, which is associated to the imposed movement (i.e. xs, zs or hs in the case of horizontal, vertical and rotational movement or a combination of them in the case of linearly varying displacement profiles). The kinematic variables associated in a virtual work sense to the contact forces are the displacement rates at the contact points, namely the tangential and normal displacement rates Dutk and Dunk (Fig. 2b). Those are collected in the vector Duk:

Duk ¼ ½ Dutk

Dunk T

ð7Þ

2.2. Equilibrium conditions

Assembling vectors fi and fs, external loads for the entire structure can be expressed as follows:

The equilibrium equations are expressed with reference to the configuration x0 related to a displacement increment no. inc. of the sequential solution procedure and can be formulated as follows:

f ¼ f D þ af L :

A0 c ¼ f

f s ¼ ð1  aÞf Ds :

ð4Þ

ð5Þ

ð8Þ

488

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

t k  lnk 6 0:

n0k

For separation condition at contact points, it is assumed that nk  0 (i.e. normal forces are assumed to be non negative in compression). In matrix notation, the limit conditions for sliding and opening failure can be written as:

R0ik Block ‘i’

t0k

g0k

ð12Þ

-t0k

y ¼ Y T c 6 0;

-n0k

ð13Þ

T

where Y is the matrix of failure conditions.

Block ‘i+1’ 2.5. Non penetration conditions Fig. 3. Initial gap at contact point k at position x0, position vector from centroid of block i and normal and tangent directions at contact point k.

where A0 is the equilibrium matrix corresponding to contact forces with coefficients determined by the position of contact points and geometry of rigid blocks, according to [35,40]. With reference to contact point k and block i in Fig. 3, the contact equilibrium matrix can be expressed as:

" A0ik ¼

^t 0k ^ 0k RT0ik n

^ 0k n RT ^t 0k

# ð9Þ

0ik

Non penetration conditions at potential contact points are formulated following the approach proposed in Krabbenhoft et al. [48] for dynamic analysis of granular materials and also adopted in [52] for the analysis of rocking behavior of masonry block structures. The conditions are formulated in a linearized form, assuming the same geometry for two subsequent increments. In this case, the condition can be expressed imposing that the normal component of the relative displacement at contact point k has not to be greater than the initial gap, as follows:

^ T0k Duk 6 g0k n

ð14Þ

^ 0k are, respectively, the column vectors with compowhere ^t 0k and n nents of tangent and normal to contact interface and R0ik is the position vector of contact point k with respect to the block centroid.

^ 0k is the initial normal associated with the surface j and g0k where n is the initial gap. For contact interactions a complementarity condition is included as follows:

2.3. Geometric compatibility conditions

 T  ^ 0k Duk  g0k ¼ 0 nk n

Geometric compatibility is another constraint to be included in the mathematical programming problem, relating displacement vector Dx to contact displacement rates Du. In analogy to Ferris and Tin-Loi [35] and under the assumption of small displacement rates, the vector Du which collects the relative displacements at potential contact points can be expressed by the following relationship:

Du ¼ AT0 Dx;

ð10Þ

with:

Dx ¼ x  x0 ;

ð11Þ

being x0 and x the initial and the new position vectors, at increment no. inc. and inc. + 1, respectively. 2.4. Contact failure conditions In the present concave contact formulation, behavior at contact points undergoing sliding failure is governed by failure conditions which are expressed according to the Coulomb friction law (Fig. 4b):

nk

tk µnk

Provided that nk  0, this condition ensures that contact forces occur only if the gap is closed otherwise are zero. The complementarity condition is represented in Fig. 4a, where it shown that positive (compression) values of normal forces are allowed when the gap is closed. Considering that an infinite compressive strength is assumed, the values of normal forces are inherently unbounded and are determined on the basis of equilibrium and sliding failure conditions. For gaps larger than zero, the contact is open and the normal force is equal to zero. Eqs. (14) and (15) express in matrix form the so-called Signorini unilateral contact conditions in terms of displacement rates. 2.6. Flow rule of displacement rates The flow rule describes the kinematic behavior at a contact point when failure occurs. This behavior is expressed by relationships between the vector of displacement rates Du and the vector k of non-negative flow multipliers associated to each failure condition. The form of the mathematical program arising from the set of governing equations is strictly related to the assumptions on flow rule. In this study, an associative behavior is assumed for flow rules. This implies that the vectors of resultant displacement rates associated to the flow multipliers k are perpendicular to the relevant failure surfaces (normality rule). As such, and also considering Eqs. (12) and (14), the following conditions hold for the associative behavior at a contact point k:

Dutk ¼ kskþ  ksk utk

g0k a)

b)

Fig. 4. Contact laws for normal reaction force (a) and shear force (b) at point k.

ð15Þ

ð16Þ

Dunk ¼ lkskþ  lksk  kok þ g 0k ; being ksk+, ksk- and kok the flow multipliers corresponding to positive, negative sliding and opening.

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

For the whole structure, it can be posed:

Du ¼ Yk þ g 0

ð17Þ

where g0 is the vector collecting initial gaps [0 g0k] at each contact point and k  0. 2.7. Complementarity conditions Activation of positive values of flow rates in case of sliding failure and contact opening relationships (15) can be expressed using the following complementarity condition:

yT k ¼ 0

ð18Þ

This condition, which also holds component wise considering the signs of y and k, imposes that flow rates occur only when stress resultants lie on the relevant failure surfaces otherwise are zero. 2.8. Energy dissipation Considering that positive energy has to be dissipated by live loads and to simplify the formulation of the mathematical program, it is assumed [53]:

f TL Dx ¼ 1

ð19Þ

3. Linear programming formulation The Eqs. (8), (10), (13), (17), (18) and (19) represent the Karush– Kuhn–Tucker (KKT) optimality conditions associated to the following linear programming problems [54]: T

min  f D Dx T

s:t: f L Dx ¼ 1

ð20Þ

Du ¼ AT0 Dx Du ¼ Yk þ g 0 ;

kP0

and:

max a  c T g 0 s:t: A0 c ¼ f D þ af L

ð21Þ T

 y P 0; y ¼ Y c; corresponding to the kinematic-based problem and its dual forcebased problem, respectively. It is interesting to note that, when g0 is omitted, previous optimization problems reduce to the standard forms of upper and lower bound formulations of limit analysis theorems presented in [35]. It is also worth noting that the above optimization problems are similar to those which arise for contact dynamics when the effect of inertial and dynamics forces are neglected [48,52]. Static variables and collapse load multipliers were obtained from the solution of the optimization problem (21). The kinematic variables in (20) were derived as Lagrange multipliers associated to the solution of the static problem (21), according to the procedure illustrated in the following section. 4. Solution procedure and implementation To take into account the effects of large displacements, an incremental solution procedure based on optimization problem (20) and (21) was implemented. The implemented algorithm can be summarized as follows:

489

(1) Solve the LP problem defined in Eq. (21) at increment no. inc. and configuration x0 and obtain base reaction load multiplier a and contact forces c; (2) Obtain kinematic variables Dx from Lagrange multipliers associated to the solution of the dual problem (20); (3) Determine new positions of the blocks x = x0 + Dx and contact gaps g0; (4) Set up a new optimization problem on the basis of the new configuration at increment no. inc.=inc. + 1 and repeat from step (1). The procedure is repeated until the maximum allowable displacement at support is attained, i.e. when stability loss of structural parts occurs with sudden changes in the base reaction values or control point displacements. Simple assumptions were introduced in the model in order to simplify the implementation of the formulation. The assumptions take advantage of the bond pattern, i.e. rectangular blocks with regular staggering, and are concerned with the specification of potential contact points and with the updating of gap functions during time history analysis. For the specification of potential contact points, four contact points k, . . ., k + 3 located at the vertexes of adjacent blocks per each pair of interfaces j and j + 1 are used (Fig. 5a). For the gap functions, the distances of contact points from potential contact surfaces have been measured as indicated in Fig. 5b. In order to simplify the implementation of the formulation, it is assumed that directions of normal and tangent vectors to contact interfaces are parallel to x and z axis and do not change with the new configurations of contact surfaces. Under the adopted assumptions, contact interactions are extended into the large displacement regime and, at the same time, the search of new contact points at each increment of the sequential solution procedure is avoided. Considering that the number of potential contact points which is monitored per each interface is fixed, it should be noted that the number of static and kinematic variables associated to the formulations of problems (20) and (21) does not change during the analysis. This is a key feature which simplifies the implementation of formulation proposed significantly, because the size of matrices governing equilibrium, kinematic and contact conditions does not change neither with the deformed configuration of the rigid block assemblage nor with the contact points state, i.e. open (Dunk < g 0k ) or closed gap (Dunk ¼ g 0k ). As such, the proposed algorithm for contact detection represents an extension to the case of gross displacement of the modelling approach adopted in [52] for the analysis of the rocking behavior of dry jointed masonry structures, where just two contact points (rather than four as it is in the present study) are monitored per each interface pair j and j + 1 and block i. Considering that only interfaces which overlap in the starting configuration are monitored (i.e. interfaces j and j + 1 in a running bond texture), shortcomings of the implemented procedure are related to the possibility that contact interfaces which do not overlap in the original configuration come into contact during the analysis in the large displacement regime (i.e. interfaces j and j + 2 in Fig. 5, when the amount of sliding displacement is larger than the length of the overlapping joint). In this case, contact interactions are no more taken into account. On the basis of the proposed formulation, a computer code has been developed in MATLAB which provides as outputs the contact forces and kinematic variables as well as the plots of the failure mechanisms at different displacement increments. The linear programming problems (20) and (21) were solved using the primal–dual interior-point solver in MOSEK [55,56].

490

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

Block ‘i’

R0ik node ‘k’

j+2

k+2

node ‘k’

R0ik+1 interface ‘j’

g0k

k+1

t0k

g0k+3

k+1 k+3

k+2

k+3

j+1

Block ‘i’

n0k

Block ‘i+1’

Block ‘i+1’

a)

b)

Fig. 5. (a) Specification of contact points for rigid block i used in the present study for the implementation of the adopted formulation; (b) initial gaps g0 for contact points k and k + 3.

A PC containing a 3.3 GHz Intel Xeon E3-1245 processor with 8 GB of RAM was used to perform numerical analyses.

5. Validation against experimental tests Two experimental case studies were analyzed to validate the proposed formulation. The first case study is concerned with a set of small-scale tuff masonry portal frames. In this case, ad-hoc experiments were carried out to investigate failure modes and maximum allowable displacements at moving supports just before the collapse. The second example is taken from the literature and is concerned with a masonry arch subjected to spreading supports.

5.1. Portal frames Three configurations of small-scale dry jointed masonry panels with opening were tested (Fig. 6). The panels were made of tuff blocks with dimensions of 100  100  50 mm. The wall specimen investigated in the first configuration is eight bricks width (800 mm), one brick thick (100 mm) and twelve courses high (600 mm). The dimension of the opening was two bricks width (200 mm) and seven courses high (350 mm). In the first configuration, the horizontal displacement at the movable support (the timber table supporting the left pier in Fig. 6) was free. For the second and third configuration, the width of the wall panels was nine bricks (900 mm) and the dimensions of the openings were three

8-9 bricks

Wall panel

Liing table

2-3 bricks

12 courses

Plaorm scale

Rollers Fixed support

Hydraulic Jack Mechanical dial gauge

a)

b)

c)

Fig. 6. (a) Experimental test set-up; (b) front view of the test set up; (c) rear view.

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

a)

491

b)

Fig. 7. Wall sample 8B_FH_2: (a) failure mode of the specimen just before the collapse at vertical displacement = 100.0 mm; (b) plot of the computed failure mode.

a)

b)

Fig. 8. Wall sample 9B_FH_1: (a) failure mode of the specimen just before the collapse at vertical displacement = 129.0 mm; (b) plot of the computed failure mode.

a)

b)

Fig. 9. Wall sample 9B_RH_1: (a) failure mode of the specimen just before the collapse at vertical displacement = 159.0 mm; (b) plot of the computed failure mode.

bricks width (300 mm). In these other configurations, free and restrained horizontal displacement conditions at the movable support were investigated, respectively. The brickwork was arranged in a stretcher bond pattern with all bricks laid as stretchers and half-bats at the beginning or at the end of alternate courses. A timber lintel was used to span the openings. For each wall panel configuration, three tests were carried out. Specimen labels are defined according to the number of bricks width (8B and 9B), to restraint conditions on horizontal displacements at the support (FH and RH in the case of free and restrained horizontal displacements) and to sample number (i.e. the ‘8B_FH_1’ label refers to the sample number 1 with eight brick width, with free horizontal displacements at the movable support, so does the ‘9B_RH_1’ in the case of restrained horizontal displacements).

The test set-up adopted for foundation settlements on the investigated wall panels was arranged according to [44]. The bearing support system was designed in order to impose vertical movements to the left masonry pier in Fig. 6a and b. The pier is supported by a lifting table placed on a hydraulic jack which was used to impose vertical movements under displacement control (Fig. 6a). The lifting table is made of a pair of timber tables with three cylindrical steel rollers positioned in between in order to allow displacements along the horizontal direction as well in the case of specimens 8B_FH and 9B_FH. As for specimen 9B_RH, the three cylindrical steel rollers positioned in between the pair of timber tables were removed in order to prevent the horizontal displacements.

492

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

Table 1 Experimental wall panels: Base reactions and maximum allowable vertical displacements at collapse. Wall sample

Model size (b  c)

8B_FH 9B_FH 9B_RH

91  944 95  992 95  992

fDs (N)

fs (N)

d (mm)

Proposed

Mean exp.

Proposed

Mean exp.

Proposed

Mean exp.

246.1 259.5 259.5

244.3 253.5 254.1

221.9 240.9 252.1

224.4 236.0 247.2

104.1 129.0 168.0

95.3 132.3 153.7

a)

c)

Diff. d proposed vs. exp. (%)

CPU time (s)

9.2 2.5 9.3

37.2 51.6 59.7

b)

d)

Fig. 10. Portal frames: (a–c) Base reaction vs. vertical support displacement; (d) Base reaction vs. horizontal support displacement for portal frames 8B_FH and 9B_FH.

The stroke of the hydraulic cylinder is 200 mm and the size of the timber tables is 600 mm  200 mm. The hydraulic jack was positioned on a platform scale with a maximum capacity of 600 N in order to measure variations of dead load distribution when the settlement is imposed. In order to measure the vertical displacements at the movable support, a mechanical dial gauge is held rigidly in a fixture which is undisturbed by the operating force of the hydraulic jack. Indeed, the gauge is supported by a rigid bearing that is separated from the test apparatus and the plunger is fixed to a steel plate that is connected to the lower timber table of the rolling bearing. The adopted dial gauge reads 0.1 mm over a range of 200 mm. A displacement rate of 0.1 mm per second was applied to the test apparatus and corresponding weight variation were registered during the test. The tests were monitored by means of two video cameras positioned in front of the panel and on the right hand side of it. The failure modes observed from testing just before the collapse for the three wall panel configurations are shown in Figs. 7a–9a. The vertical movement imposed to the left pier involves the rotation of the timber lintel and of a part of the supported spandrel around the right end support the timber lintel. This rotation is associated to the formation of stepped diagonal cracks at the top of the right pier and to vertical cracks in the middle part of the left pier due to sliding failure. In the case of wall panels 8B_FH and 9B_FH, this rotation is also associated to an horizontal rigid body movement of the left pier, which is free to slide.

The collapse of the wall panels is attained when the relative sliding displacement between the left bottom corner of the lintel and the upper right corner of the supporting tuff block below is greater than the staggering length (i.e. half a block length in this case, given the running bond pattern). For the numerical simulations, unit weights of tuff masonry blocks and timber lintel were taken as 12.3 and 4.3 kN/m3. The friction coefficients at tuff-to-tuff and tuff-to-timber joints were taken as 0.72 and 0.66, respectively, in accordance to values measured experimentally. In the case of wall panels 8B_FH and 9B_FH, a zero value of the friction coefficient was adopted at the interfaces between the movable supports so that horizontal displacements were allowed, according to experimental set-up. The starting values of dead loads applied at the movable support in the numerical model were equal to 246.1 and 259.5 N. Those values correspond to the weight of half the wall panel 8B and 9B which was assumed to act on the support block in the initial configuration. The results of the numerical analysis are reported in Figs. 7b– 9b, 10 and in Table 1. The comparison of failure modes with testing outcomes in Figs. 7–9 shows that the collapse mechanisms predicted by the numerical models for the two wall panels are in good agreement with the experimental tests. In all cases, the predicted failure modes involve the rotation of the spandrel over the right pier and the vertical and horizontal movement of the left support. According to experimental tests, the maximum displacement capacity is attained when the relative sliding displacement of the left bottom corner of the timber lintel overcomes half the length of the supporting tuff block.

493

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

The experimental and calculated base reaction-displacement curves are plotted in Fig. 10. On the basis of the experimental curves, the behavior of the wall panels can be schematized in three phases. In the first phase, the response is characterized by high-rate decreasing values of the reaction-displacement curve. In this phase, the failure mechanism starts to develop progressively within the wall panel, with random activation of failure (i.e. toppling and sliding) at dry joints due to the variability of mechanical properties and geometric imperfections. In the second phase, the failure mechanism is completely developed in the whole wall panel and the magnitude of the base reaction is approximately constant with increasing displacements at the movable support. During the third phase, the base reaction increases slightly with increasing movements at base support due to the effects of large displacements. As such, it is worth noting the similarity of base reactionsdisplacement curves with the curvilinear response of out-ofplane bended masonry panels under seismic loads, where analogous relationships are obtained for lateral forces in the large displacement regime and similar idealizations (i.e. tri-linear models) are adopted [57]. The starting values of dead loads fDs applied at the movable support in the numerical models are also reported for comparison in Fig. 10a–c. It can be noted that these values are slightly different from the base reactions corresponding to the starting configuration (i.e. zero movement at the support) measured from experimental tests. Clearly, this is due to inherent imperfections of the tuff masonry block assemblages which involve a slight nonsymmetric distribution of dead loads also in the starting configuration of the wall panels. Slight differences between the predicted and experimental value of the base reactions can be observed for small values of support displacement. This can be ascribed to the assumption of indeterminate (and infinitesimal) displacement rates for the activation of the failure mode which is inherent to each increment of the proposed sequential numerical procedure. Conversely, the results of experimental tests show that large displacements are necessary for the mobilization of the complete failure mechanism. Comparisons with the results obtained from the limit analysis formulation proposed in [44] which is associated to indefinite displacement are also reported to point out differences with testing outcomes and the new model.

Most importantly for the scope of the present study, the comparison with experimental results shows that the values of the computed base reactions and the maximum vertical movement imposed at the base block just before the collapse closely match those from testing. Differences between predicted and experimental results are lower than 10% in terms of maximum displacements and even lower in the case of base reaction at collapse, which are reasonable values considering that the effects of geometric imperfections also influence the behavior of tested specimens. Fig. 10d shows the curves of the base reaction versus horizontal displacement at the movable support computed from the numerical model in the case of specimens 8B_FH and 9B_FH. It can be noted that in these cases the predicted horizontal displacements are about 28.0 and 31.0 mm just before the collapse. It is clear that the free component of the horizontal displacement at the support plays an important role on the maximum allowable vertical displacement. The comparison between specimens 9B_FH and 9B_RH shows that the vertical displacement just before the collapse increases of about 16% in the case of the fixed horizontal displacement. 5.2. Masonry arch on spreading supports The second experimental case study is a small scale masonry arch on spreading supports presented in [25]. The same case study was also analyzed by other authors to test the accuracy of other numerical and analytical formulations based on limit analysis theorems [29]. The application of the proposed approach to this case study shows the potentialities of the adopted point-contact model, which can be applied to non-rectangular or even irregular block shapes as well, by simply using contact points located at the corners of each interface along the edges of the rigid body. The model of the arch was made of voussoirs cast as concrete blocks and had a mean radius of 220 mm. The radial thickness was 50 mm. The arch was tested by spreading the support till collapse. The rigid block model developed for numerical simulation is shown in Fig. 11a. The unit weight was taken as 25.0 kN/m3, according to [25], and friction coefficient was set equal to 0.7. The numerical results in terms of hinges position and displacements just before the collapse and relevant experimental outcomes are presented in Fig. 11b and Table 2. According to experimental

Fig. 11. Masonry arch on spreading supports [25]: (a) Dimensions (mm); (b) predicted failure mode just before the collapse.

Table 2 Masonry arch on spreading supports: comparison of numerical and experimental results [25]. Sample

Circular arch

Model size (b  c)

17  68

Hinge position

% Span increase at collapse

CPU time (s)

Exp.

Proposed

Exp.

Proposed

Diff

56°

56°

15.4%

16.5%

7.1%

8.9

494

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

outcomes, small outward movement of the abutments in the numerical model causes the formation of three hinges, which are associated to a stable rigid body mechanism. The predicted position of these three hinges is in full accordance with experimental tests. The maximum span increase is attained when the stability conditions of the bottom arch segments are exceeded, i.e. when the location of the thrust line is associated to a five hinge failure mechanism. The predicted support displacement at collapse for the investigated arch is 64.3 mm. From the comparison with experimental outcomes reported in Table 2 and also considering the effects of imperfections [25], it can be concluded that also in this case the predicted results are in good agreement with those obtained from testing. 6. Numerical case studies

1.00

In this section, the numerical case study of the portal frame shown in Fig. 12a is analyzed in order to show capabilities of the proposed formulation and to evaluate its computational efficiency when full scale structures are investigated. The model is composed of 568 blocks with dimensions 260  100  250 mm and 6364 contact points. Unit weight of 16.0 kN/m3 and friction coefficient of 0.6 were adopted for numerical analysis. The portal frame is subjected to a linearly varying vertical movement at the block supporting the left pier. The linearly varying vertical movement is associated to a rotation of the supporting block around a point which is positioned at a distance of 2.60 m on the right side from the center of the left pier. A friction coefficient equal to 0.6 was also assumed at the support block interfaces to prevent free horizontal displacements as in the case of the wall panel 9B_RH.

The effects of different loading conditions and the efficiency of strengthening ties on the maximum displacement capacity of the portal frames subjected to settlement are evaluated with respect to three different structural schemes. In the first scheme, the structure is subjected to dead loads associated to unit weight of masonry blocks. In the second case the effects on the maximum displacement capacity of a uniformly distributed floor load of 50 kN/m are examined. In the third scheme the effects of a strengthening tie positioned at the lintel course are taken into account. To model tie effects, two opposite forces of 70.0 kN where applied at the center of the blocks where the strengthening element is assumed to be positioned. To take into account the effects of large displacements, the direction of these forces is updated at each increment of the sequential procedure. The failure mechanisms predicted by the numerical analyses are shown in Figs. 12b–d. It can be seen that introducing distributed loads and tie forces significantly influences the failure mechanism and the maximum allowable displacement capacity. For the case of the wall panel with ties, a gap can be noted in the numerical displaced configuration between the left pier and the spandrel (Fig. 12d). This gap is a consequence of the associative behavior assumed for sliding failure, which involves dilatancy between contact interfaces, as well as of the simple assumptions adopted for contact point specification [52]. Despite this geometric separation, contact forces still exist between the left pier and the spandrel, as it can be also noted from the displaced configuration of the left pier and from the diagonal crack at its base. In Fig. 13a the vertical component of the displacement at the left bottom corner of the lintel (which is the control point in this case) is plotted against the displacement at the support block for the three wall panel configurations. It can be seen that the collapse displacement predicted for the wall panel subjected to unit weight

3.30 2.30

Control point

1.30

2.60

1.30 5.20

a)

c)

b)

d)

Fig. 12. Portal frame subjected to linearly varying movement at the bottom of the left pier: (a) geometry and loading conditions; (b) Failure mode at incipient collapse under self weight loads; (c) under self weight and floor loads; (d) under self weight and floor loads and with tie forces.

495

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

b)

a)

Fig. 13. (a) Displacement curves of the numerical portal frame under different loading conditions; (b) Base reaction vs. vertical support displacement.

Fig. 14. Sensitivity of the response of the numerical wall panel to friction coefficient.

is 206.1 mm and that the failure mechanisms of the other configurations are stable under support displacements up to 300.0 mm. The curves of the base reaction versus vertical displacement at the support are shown in Fig. 13b. In the case of the wall panel under self-weight loads, a gap in the value of the base reaction at the collapse displacement can be noted, which corresponds to the weight of the spandrel supported by the left pier just before the collapse. The difference between the values of the base reactions in the cases of wall panels with floor load is clearly due to the effect of ties, which deflect the flow of forces to the fixed right pier. Considering that in the proposed formulation the sole mechanical parameter influencing the structural response is the friction coefficient, a sensitivity analysis was also carried out on the response of the numerical case study varying m in the range 0.5– 0.7. The results of the sensitivity analysis are plotted in Fig. 14. It is shown that varying the friction coefficient within the range slightly influences the collapse displacement (by up to 10 percent). For computational efficiency, the average CPU time required to solve each displacement increment and to complete the analysis were 4.4 and 222.6 sec, respectively, which are reasonable values considering the size of the rigid block model and the number of the variables involved.

and uses linear programming to solve the associated incremental equilibrium problem. A simple procedure for specification of potential contact points and for the updating of the gap functions was proposed to avoid the search for new contact between blocks. The procedure assumes that contact forces are applied to the vertexes of potential contact interfaces between blocks and is able to take into account large displacement and rotations. The model was validated against experimental data from tests carried out within the present study and also from the literature. A good agreement in terms of failure mechanisms, base reactions and collapse displacements was observed, with percentage differences less than 10 percent. However, a remarkable difference between experimental and numerical curves of base reaction versus vertical displacement was observed for small values of vertical support movement. The difference is due to the assumption of indeterminate and infinitesimal displacement rates for the activation of the failure mode which is inherent to each increment of the proposed sequential numerical procedure. Applications to numerical case studies showed that the computational efficiency and the stability of the implemented procedure are suitable to model static gross displacement behavior of masonry block assemblages subjected to settlement with a large number of rigid bodies. The computational efficiency of the developed model is mainly related to the use of a static formulation in the large displacement regime. This allowed for the use of linear programming to cast the corresponding contact problem, for which fast and robust solver exists in the literature. However, as a drawback, considering that inertial and dynamic effects were not taken into account, it is clear that the model adopted is not able to analyze the progressive collapse of the structure once the maximum allowable displacement has been attained. Moreover, earlier local failure modes (i.e. involving the incipient collapse of few blocks rather than global collapse) might result in numerical convergence problems which could prevent the analysis of the whole behavior of the structure. To overcome previous issues, further development of the present study could be devoted to the analysis of settlements using quasi-static modelling in the framework of contact dynamics [52]. More refined algorithms could be also implemented for contact point specification on the basis of different approaches available in the literature, such as those used in the case of the distinct element method.

7. Conclusions Acknowledgements A novel 2D formulation for large displacement analysis of dryjointed masonry block structures subjected to foundation settlements was presented. The formulation is based on the contact point model for interfaces (for which a no-tension and associative frictional behavior with infinite compressive strength is assumed)

The financial support of PRIN 2015 Programme by the Italian Ministry of Education, Universities and Research (MIUR) is gratefully acknowledged for funding the research project ‘‘Protecting the Cultural Heritage from water-soil interaction related threats

496

F. Portioli, L. Cascini / Engineering Structures 148 (2017) 485–496

– PERICLES” (Prot. No. 2015EAM9S5), which is the main framework of the study presented in this article. The authors wish also to express their gratitude to Mario Torricella from the Laboratory of the Department of Structures for Architecture for his assistance and support in the preparation of the specimens, test set-up and throughout the execution of experimental investigations. Finally, the authors are grateful to the three anonymous reviewers for their constructive comments and suggestions which helped to significantly improve the quality and readability of this article. References [1] DeJong MJ. Settlement effects on masonry structures. In: Proc of 10th Int. Conf. on Structural Analysis of Historical Constructions (SAHC 2016), 13–15 September 2016, Leuven, Belgium. CRC Press. [2] Foraboschi P. Specific structural mechanics that underpinned the construction of Venice and dictated Venetian architecture. Eng Fail Anal 2017;78:169–95. [3] Cardani G, Coronelli D, Angjeliu G. Damage observation and settlement mechanisms in the naves of the Cathedral of Milan. In: Proc. of 10th Int. Conf. on Structural Analysis of Historical Constructions (SAHC 2016), 13–15 September 2016, Leuven, Belgium. CRC Press. p. 623–30. [4] Valore C, Ziccarelli M. The preservation of Agrigento Cathedral. In: Proc of the 18th Int. Conf. on Soil Mechanics and Geotechnical Engineering, 2–6 2013, Paris. p. 3141–44. [5] Calabresi G. The role of geotechnical engineers in saving monuments and historic sites. In: Kerisel lecture, Proc of the 18th Int. Conf. on Soil Mechanics and Geotechnical Engineering, September 2–6 2013, Paris. p. 71–83. [6] Foraboschi P. Resisting system and failure modes of masonry domes. Eng Fail Anal 2014;44:315–37. [7] Foraboschi P. The central role played by structural design in enabling the construction of buildings that advanced and revolutionized architecture. Constr Build Mater 2016;114:956–76. [8] Amorosi A, Boldini D, de Felice G, Malena M, Sebastianelli M. Tunnellinginduced deformation and damage on historical masonry structures. Geotechnique 2014;64(2):118–30. [9] Alessandri C, Garutti M, Mallardo V, Milani G. Crack patterns induced by foundation settlements: Integrated analysis on a renaissance masonry palace in Italy. Int J Archit Heritage 2015;9(2):111–29. [10] Reccia E, Milani G, Cecchi A, Tralli A. Full 3D homogenization approach to investigate the behavior of masonry arch bridges: the Venice trans-lagoon railway bridge. Constr Build Mater 2014;66:567–86. [11] Giardina G, van de Graaf AV, Hendriks MAN, Rots JG, Marini A. Numerical analysis of a masonry façade subject to tunnelling-induced settlements. Eng Struct 2013;54:234–47. [12] Giardina G, Hendriks MAN, Rots JG. Sensitivity study on tunnelling induced damage to a masonry façade. Eng Struct 2015;89:111–29. [13] Caliò I, Marletta M, Pantò B. A new discrete element model for the evaluation of the seismic behaviour of unreinforced masonry buildings. Eng Struct 2012;40:327–38. [14] Caliò I, Cannizzaro F, Pantò B, Oliveto F. Evaluation of foundation settlements in masonry buildings using non-linear static analysis carried out with 3DMacro Software. In: Proc. of XV Convegno ANIDIS – L’Ingegneria Sismica in Italia. Padova; 2013. [in Italian]. [15] Rossi M, Calderini C, Lagomarsino S. Experimental testing of the seismic inplane displacement capacity of masonry cross vaults through a scale model. Bull Earthq Eng 2016;14(1):261–81. [16] Milani G, Rossi M, Calderini C, Lagomarsino S. Tilting plane tests on a smallscale masonry cross vault: experimental results and numerical simulations through a heterogeneous approach. Eng Struct 2016;123:300–12. [17] Al-Heib M. Distinct element method applied on old masonry structures. In: Peep Miidla, editor. Numerical modelling. InTech; 2012. [18] Bui TT, Limam A. Masonry walls under membrane or bending loading cases: experiments and discrete element analysis. In: Civil-comp proceedings; 2012. [19] McInerney J, DeJong MJ. Discrete element modeling of groin vault displacement capacity. Int J Archit Heritage 2015;9(8):1037–49. [20] Bui TT, Limam A, Sarhosis V, Hjiaj M. Discrete element modelling of the inplane and out-of-plane behaviour of dry-joint masonry wall constructions. Eng Struct 2017;136:277–94. [21] Galassi S, Paradiso M, Tempesta G. Non-linear analysis of masonry structures subjected to external settlements. Open J Civ Eng 2013;3:18–26. [22] Baraldi D, Cecchi A. A full 3D rigid block model for the collapse behaviour of masonry walls. Eur J Mech A Solids 2017;64:11–28. [23] Bertolesi E, Milani G, Casolo S. Homogenization towards a mechanistic Rigid Body and Spring Model (HRBSM) for the non-linear dynamic analysis of 3D masonry structures. Meccanica 2017. http://dx.doi.org/10.1007/s11012-0170665-6. [24] Gilbert M, Smith CC, Pritchard TJ. Masonry arch analysis using discontinuity layout optimisation. Proc Inst Civ Eng: Eng Comput Mech 2010;163 (3):155–66.

[25] Ochsendorf JA. The masonry arch on spreading supports. Struct Eng 2006;84 (2):29–35. [26] Block P, Ciblac T, Ochsendorf J. Real-time limit analysis of vaulted masonry buildings. Comput Struct 2006;84(29–30):1841–52. [27] Como M. Statics of historic masonry constructions. New York: Springer; 2013. [28] Pugi F, Galassi S. Seismic analysis of masonry voussoir arches according to the Italian building code. Ingegneria Sismica (Int J Earthq Eng) 2013;30(3):33–55. [29] Coccia S, Di Carlo F, Rinaldi Z. Collapse displacements for a mechanism of spreading-induced supports in a masonry arch. Int J Adv Struct Eng 2015;7:307–20. [30] Palmisano F, Elia A. Assessment of masonry buildings subjected to landslide by using the load path method. Int J Civ Eng 2014;12(2 A):312–30. [31] Palmisano F, Elia A. Behaviour of masonry buildings subjected to landslideinduced settlements. Int J Struct Eng 2014;5(2):93–114. [32] Palmisano F, Elia A. Shape optimization of strut-and-tie models in masonry buildings subjected to landslide-induced settlements. Eng Struct 2015;84:223–32. [33] Livesley RK. A computational model for the limit analysis of three-dimensional masonry structures. Meccanica 1992;27(3):161–72. [34] Gilbert M, Melbourne C. Rigid-block analysis of masonry structures. Struct Eng 1994;72(21):356–61. [35] Ferris M, Tin-Loi F. Limit analysis of frictional block assemblies as a mathematical program with complementarity constraints. Int J Mech Sci 2001;43:209–24. [36] Orduña A, Lourenço PB. Three-dimensional limit analysis of rigid blocks assemblages. Part I: Torsion failure on frictional joints and limit analysis formulation. Int J Solids Struct 2005;42(18–19):5140–60. [37] Whiting E, Ochsendorf J, Durand F. Procedural modeling of structurally-sound masonry buildings. ACM Trans. Graph. 2009;28(5):1–9Article 112. [38] Portioli F, Cascini L, D’Aniello M, Landolfo R. A rigid block model with cracking units for limit analysis of masonry walls subject to in-plane loads. In: Civilcomp proceedings; 2012. p. 99. [39] Portioli F, Cascini L, Casapulla C, D’Aniello M. Limit analysis of masonry walls by rigid block modelling with cracking units and cohesive joints using linear programming. Eng Struct 2013;57:232–47. [40] Portioli F, Casapulla C, Gilbert M, Cascini L. Limit analysis of 3D masonry block structures with non-associative frictional joints using cone programming. Comput Struct 2014;143:108–21. [41] Portioli F, Casapulla C, Cascini L. An efficient solution procedure for crushing failure in 3D limit analysis of masonry block structures with non-associative frictional joints. Int J Solids Struct 2015;69–70:252–66. [42] LimitState Ltd. RING Manual. Sheffield UK; 2014. [43] Acikgoz S, Pelecanos L, Giardina G, Aitken J, Soga K. Distributed sensing of a masonry vault during nearby piling. Struct Control Health Mon 2017;24(3):art. no. e1872. [44] Portioli F, Cascini L. Assessment of masonry structures subjected to foundation settlements using rigid block limit analysis. Eng Struct 2016;113:347–61. [45] Gilbert M. Limit analysis applied to masonry arch bridges: state-of-the-art and recent developments. In: Proc of the 5th Int. Conf. on Arch Bridges (ARCH’07), 12–14 September 2007, Madeira, Portugal. [46] Foraboschi P, Vanin A. Non-linear static analysis of masonry buildings based on a strut-and-tie modelling. Soil Dyn Earthq Eng 2013;55:44–58. [47] Gilbert M. Gross displacement mechanism analysis of masonry bridges and tunnel linings. In: Proc. of the 11th Int. Brick/Block Masonry Conference, 14– 16 October 1997, Tongji University, Shanghai, China. [48] Krabbenhoft K, Lyamin AV, Huang J, Vicente da Silva M. Granular contact dynamics using mathematical programming methods. Comput Geotech 2012;43:165–76. [49] Krabbenhoft K, Huang J, Da Silva MV, Lyamin AV. Granular contact dynamics with particle elasticity. Granular Matter 2012;14(5):607–19. [50] Huang J, da Silva MV, Krabbenhoft K. Three-dimensional granular contact dynamics with rolling resistance. Comput Geotech 2013;49:289–98. [51] Meng J, Huang J, Sheng D, Sloan SW. Granular contact dynamics with elastic bond model. Acta Geotech 2016. http://dx.doi.org/10.1007/s11440016-0481-5. [52] Portioli F, Cascini L. Contact dynamics of masonry block structures using mathematical programming. J Earthquake Eng, http://dx.doi.org/10.1080/ 13632469.2016.1217801. [53] Maier G, Nappi A. A theory of no-tension discretized structural systems. Eng Struct 1990;12:227–34. [54] Lloyd Smith D. Linear programming. In: Lloyd Smith D, editor. Mathematical programming methods in structural plasticity. New York: Springer-Verlag; 1990. [55] Andersen ED, Roos C, Terlaky T. On implementing a primal-dual interior-point method for conic quadratic optimization. Math Program 2003;95(2). [56] The MOSEK optimization tools manual. Version 8.0. MOSEK ApS, Denmark; 2017. . [57] Doherty K, Griffith MC, Lam N, Wilson J. Displacement-based seismic analysis for out-of-plane bending of unreinforced masonry walls. Earthquake Eng Struct Dynam 2002;31(4):833–50.