Damage propagation in a masonry arch subjected to slow cyclic and dynamic loadings

In the present work, the damage propagation of a masonry arch induced by slow cyclic and dynamic loadings is studied. A two-dimensional model of the arch is proposed. A nonlocal damage-plastic constitutive law is adopted to reproduce the hysteretic characteristics of the masonry material, subjected to cyclic static loadings or to harmonic dynamic excitations. In particular, the adopted cohesive model is able to take into account different softening laws in tension and in compression, plastic strains, stiffness recovery and loss due to crack closure and reopening. The latter effect is an unavoidable feature for realistically reproducing hysteretic cycles. In the studied case, an inverse procedure is used to calibrate the model parameters. Then, nonlinear static and dynamic responses of the masonry arch are described together with damage propagation paths.


INTRODUCTION
ecent seismic events have led to an increased demand for the development of effective methods able to discern on the safety of existing masonry structures under dynamic loadings; in this context, arches and vaults have confirmed to be vulnerable with respect to the earthquakes showing the occurrence of extensive damage [1].In this respect, while the static analysis of arches is widely investigated, a better understanding of the response of masonry vaults and arches under dynamic loading conditions is necessary for the evaluation of the behavior of these structural elements when located in seismically active areas.The seismic capacity of masonry arches is usually studied considering the well-known Heyman's hypotheses [2], who assumed infinite compressive strength and no tensile resistance.This modeling approach is often adopted in the framework of limit analysis [3].The method provides analytical solutions and can lead to approximated estimates of the seismic capacity of arches.The results of the limit analysis have also been used to validate numerical simulations performed using the Distinct Element Method (DEM), which describes the structure as an assemblage of rigid blocks and nonlinear interfaces [4].Although, the DEM method allows an accurate description of the dynamic behavior of masonry arches composed by blocks, it could be however inadequate for studying three-dimensional structures made out of a large number of elements.R A further alternative is the use of the Finite Element Method (FEM), by adopting a micro-modeling approach [5,6] or a macro-modeling approach [7,8].The computational analysis of these structures subjected to dynamic loadings requires realistic stress-strain material models to satisfactorily reproduce their physical behavior.Research efforts on the static response of the masonry material under cyclic loadings aim at providing an efficient model capable of predicting the hysteretic characteristics of the entire element.Common models used to reproduce the cyclic behavior of cohesive material are based on damage mechanics, plasticity theory and coupling of both.As widely known, if the stress-strain laws with softening are used within the standard continuum theory, the obtained numerical results suffer from pathological sensitivity to the spatial discretization, e.g. to the size of the finite elements.Upon mesh refinement, the energy dissipated by the numerical model decreases and tends to extremely low values, sometimes even to zero.As remedy, regularized models, based on nonlocal continuum approaches, can be adopted [9][10][11][12].Recently, efforts have been made to develop dynamic tests which appear particularly promising for a deeper understanding of the response of masonry arches and vaults.Within this activity structural damage identification has been pursued to extract information on the health of these structures.Indeed, the presence of damage or of other inelastic phenomena modifies the overall structural dynamic response and the damage propagation potentially interacts dynamically with the element vibrations; in particular, changes in its behavior are associated to the decay of the mechanical properties of the system.Based on these considerations, several studies have been devoted to use the variations in the dynamic behavior to detect structural damage.Particular attention has been focused on the use of frequencies only, because of the simplicity of measuring them and, therefore, their experimental reliability.In this framework, dynamic analyses of damaged structures have been performed in [13] with the aim to detect the damage state of an arched masonry structure.In the present paper the damage propagation in a large-scale masonry arch induced by cyclic and dynamic loadings is analyzed.In particular, the masonry arch tested in [13] is considered as case study.A two-dimensional modeling of the masonry arch is proposed.The constitutive material model developed by Toti et al. [12] and extended to the dynamic analyses in [14] is used to reproduce the main hysteretic characteristics of the masonry material under cyclic or dynamic loadings, such as damage, inelastic strains, unilateral phenomena.An inverse procedure is used in order to calibrate the model parameters.The performance of the developed computational tool, in reproducing the experimental response of the arch and the damaging evolution, are investigated in the case of a static cyclic test.The damage propagation and the variations of the dynamic behavior (displacement, frequency contents) with respect to the undamaged condition are examined, when the studied mechanical system is excited by imposed base synchronous harmonic motions.

THE MASONRY ARCH CASE STUDY
he masonry arch studied during the experimental campaign presented in [13] is considered as case study.The masonry structure is built with standard clay bricks with size 100×50×25 mm3.The geometrical data of the tested circular arch, schematically reported in Fig. 1, are the following:  internal radius r = 0.77 m;  width w = 0.45 m;  thickness t = 0.05 m;  height h = 0.59 m;  span l = 1.50 m.Fig. 1 shows even the coordinate reference system adopted for the computations.Static tests to induce damage and dynamic identification tests are performed during the experimental campaign [13].During the static test, the damage of the arch is caused by the application of a point load located with a distance d from the abutment equal to 0.38 m (about a quarter span), where the lowest safety factor for arches is obtained.The point-wise load induces a damage evolution in the structure, which leads a four hinges collapse mechanism.Dynamic tests for modal structural identifications have been conducted on the masonry arch through vibration measurements [13].The frequencies and the modal shapes are estimated both for the undamaged arch and for damaged arch for different values of the increasing static force.As strains are measured at 11 points of the extrados and intrados of the arch, it was also possible to evaluate the centerline curvature mode shape.Here, the novel two-dimensional model is calibrated in order to reproduce the measured quantities on the experimental apparatus presented in [13], in particular the values of modal frequencies obtained by vibration measurements are compared in Tab. 1 with the same quantities obtained by the proposed numerical model.For sake of brevity, only the modal parameters related to first two in-plane modes of vibration of the arch computed in the undamaged condition are herein considered.

Nonlocal damage-plastic model
nonlocal damage-plastic model [12] is herein adopted to model the macroscopic cyclic behavior of the material point of the examined structure.The model is able to account for the tensile and compressive damage, the accumulation of irreversible strains and the unilateral phenomenon, i.e. the stiffness recovery and loss due to crack closure and reopening.The stress-strain constitutive relationship is given by the following expression: with where σ and σ are the stress tensor and effective stress tensor, respectively; ε , π and e are the total strain, the plastic strain and elastic strain, respectively; C is the fourth-order constitutive tensor; t D and c D are two damage variables which describe the stiffness degradation of the masonry in tension and compression, respectively; 1 ( ) e J tr  e is the first invariant of the elastic strain;   H x denotes the Heaviside function,   1 The evolution of the plastic strain is described by introducing a yield function representing a branch of a modified hyperbola and the loading-unloading conditions in the Kuhn-Tucker form: As a damage softening constitutive law is introduced, the localization of the strain and damage variables could occur.In order to overcome this pathological problem, to account for the correct size of the localization zone and, also, to avoid strong mesh sensitivity on the numerical results in finite element analyses, a nonlocal constitutive law is considered both for the compressive damage and the tensile damage.In particular, the evolution of the compressive damage c D is combined with the development of the plastic strain through the following cubic relationship: where u  is the final accumulated plastic strain associated with the compressive damage counterpart of the accumulated plastic strain given by the following expression: x y the compressive weight function which determines the influence of the point y on x ; the bracket symbol . the positive part of the A number; the parameter 2 c R is the compressive characteristic length, which contributes to the definition of nonlocal accumulated plastic strain.The evolution of the tensile damage parameter t D is governed through the following nonlocal exponential law: where 0  is a material parameter indicating damage threshold strain; eq  is the nonlocal equivalent strain defined as: x y defines the tensile weight function, with 2 t R the tensile characteristic length, which contributes to the definition of nonlocal equivalent strain.
Finally, the condition that the damage in compression must be lower than the one in tension is enforced t c D D  .The constitutive law defined in Eq. ( 1), for particular strain paths and damage conditions can be lead to discontinuities in the response of the material point.Indeed, for an isotropic cohesive material subjected to 1 J history and to a constant elastic shear strain xy e , the shear stress is J changes from positive to negative value or vice-versa.While it can be considered realistic to have a stiffer shear response when the material point is subjected to volumetric contraction with respect to the case of volumetric expansion, because of the positive effect of the friction in compression, this strong discontinuity of response of the material point can be considered undesirable from both a physical and a mathematical point of view.For this reason, a regularized form of the Heaviside function can adopted with the aim to avoid annoying jumps in the response; in particular, the following regularized form is considered for the Heaviside function: where h is a small parameter governing the regularization effect (h=0.0001  0.001).

Finite element mesh
A two-dimensional plane stress finite element model of the arch has been constructed.The opportune choice of the element size has been investigated demonstrating that a finite element model characterized by a total of 320 four node isoparametric quadrilateral elements (Fig. 3a) is accurate.The chosen mesh assures a satisfactory numerical result, particularly with reference to the robustness in the frequency evaluation obtained by modal analysis (see the following section).Indeed, considering refined meshes characterized by a higher number of finite elements, the first fundamental modal frequencies do not change significantly.Concerning the constitutive laws, the above introduced cohesive continuum formulation is considered to model the macroscopic behavior of the masonry material.

CALIBRATION
n inverse procedure, based on the results of the dynamic and static tests carried out in [13], is herein developed in order to derive some material properties of the arch.
The Young's modulus E , the Poisson's ratio  and the mass density  are set through a finite element updating procedure, using the modal structural properties estimated with the dynamic tests.In particular, the following objective function  to be minimized is considered: ,exp -1 , , with 2 where N denotes the number of the autofrequencies, ,i r  is the residual, .Fig. 3b-c show the modal shapes and the corresponding autofrequencies obtained by the developed updating technique.The comparison between the experimental modal parameters and the numerical one, provided in Fig. 2 and in Tab. 1, confirms a good accordance in term of curvature mode shapes and autofrequencies.Moreover, in Tab. 2 the comparison between the considered mesh with 320 elements and two different meshes of 160 and 640 elements, respectively, shows that the chosen discretization assures the converge of the results. and k are depicted in Fig. 4 and Fig. 5, respectively, where the numerical response, plotted in term of nodal reaction versus the prescribed displacement, is compared with the experimental data.The set parameters governing the nonlinear behavior of the material are reported in Tab. 3, with

Natural frequencies
It is also possible to observe from Fig. 4 and Fig. 5 that the elastic properties (i.e.E and  ) defined through the updating procedure provide a good estimation of the structural initial stiffness computed through the static analyses.

Static-cyclic test
he mechanical behavior of a masonry arch under cyclic loading is analyzed.The static-cyclic test performed in [13] is reproduced.The loading history reported in Tab. 4 is applied to the displacement v in order to simulate the cycles of the experimental test.The comparison between the numerical result and experimental data of the global mechanical behavior is illustrated in Fig. 6.Observing this figure, it is possible to remark that the proposed modeling approach is able to catch with good approximation the actual behavior of the structural system occurred during the loading and unloading phases.Fig. 7 shows the tensile damage map obtained by the simulation in correspondence of the two time steps indicated in Fig. 6 with the letter A and B, respectively.From Fig. 7, it appears evident that the collapse mechanism developed by the arch is characterized by the formation of four hinges, two on extrados and two on intrados.With reference to [13], the numerical model well reproduces the failure mechanism developed by the arch during the static test, indeed, the tensile damage assumes close to 1 in the zones where the fracture opening occurred.[mm] 0 0.38 0.04 0.42 0.07 0.45 0.1 0.9 0.24 v Table 4: Loading history for the displacement v

Dynamic test
Numerical analyses are conducted to investigate the damage propagation induced by a sinusoidal synchronous motion applied to all nodes of base , where 0 u and 0 f indicate the displacement amplitude and the frequency of the harmonic imposed displacement, respectively.The propagation of the damage in the arch is analyzed for frequency ratio 0 1 f f   equal to 0.7 and 1.3.For each frequency ratio, two values of the displacement amplitude are considered: 0 0.10 mm u  and 0 0.11 mm u  .The equations of motion are integrated until the time 5 s is reached, with a time step equal to 0.001 s (approximately 1/14 of the second natural period of vibration) in order to correctly take into account the dynamics of the first two modes.The comparison between the damage model and the elastic model is provided in term of time histories of the displacement v evaluated at the extrados of the key-section of the arch.The effect of damage on the dynamic response is analyzed through the Fourier spectral analyses.Fig. 8 and Fig. 9 show the displacement response and Fourier spectra for all the considered forcing functions, adopting both elastic and the proposed damage model.Observing these figures, it can be remarked that: the arch collapses before 5 s for 0.7 ; the development of the damage in the nonlinear model leads, away from the failure mechanism, only to a decrease in time of the natural frequencies, while the frequency contents become richer for all cases close to the collapse condition; the reduction of the frequency contents became faster with increasing the forcing displacement amplitude.In particular, for 0.7   with 0 0.11 mm u  at the time equal to 1.5 s the first natural frequency of the damaged system is approaching the lower forcing one, the arch goes in resonance following its first modal shape and, as consequence, collapses.On the contrary, for the other case, characterized by lower forcing displacement amplitudes, the structural system does not reaches after 5 seconds the critical condition; the damage of the arch leads only to energy dissipation.In fact, the peak amplitudes associated to the frequency contents appear lower than the ones obtained with the elastic model, where the energy dissipation is always null.Instead, the collapse of the arch for the case 1.3

 
is not due to the resonance phenomenon, but it occurs because the higher forcing frequency leads to a larger base acceleration and, as consequence, higher inertia forces.In particular, it can be remarked that for these three analysis cases, before 5 s, the arch tends to develop a four hinges collapse mechanism compatible with the first mode of vibration.Indeed, the tensile damage variable assumes higher values exactly in the four zones where the maximum curvature values of first mode occur (Fig. 2).The damage variable reaches first the unitary value at the two support areas because in these regions the bending curvatures are larger than the ones of other two zones (Fig. 2).

CONCLUSIONS
he present work focuses on the development of a computational strategy for studying the damage propagation in a masonry arch induced by slow cyclic and dynamic loadings.A nonlocal plastic damage model is adopted in order to reproduce the cyclic macroscopic behavior of the masonry material.An inverse procedure is used to define the material model parameters.In particular, the mass density and the elastic properties are set through an updating procedure, using the first two modal natural frequencies identified by dynamic tests.The found solution assures a good accordance of the modal parameters in term of curvature mode shapes and natural frequencies of the first two modes of vibrations between experimental and numerical quantities.The calibration of post peak parameters of the proposed formulation are set on the basis of the experimental nonlinear response determined during a static test.A numerical example concerning a slow cyclic test shows the effectiveness of the developed arch models in reproducing both the experimental response and the collapse mechanism.Finally, the masonry arch is excited by sinusoidal imposed base motions to study the effect of both the amplitude and the frequency content of the base motion on damage propagation.The numerical computations conducted in the time and analyzed in the frequency content, show that the damage model is able to catch the dissipation of energy and the reduction of the natural frequencies.In particular, the base motions at forcing frequency higher than the first natural one is accompanied by large accelerations which induce the collapse of the structure.In the case of forcing frequency lower than the first natural one due to the stiffness degradation, even in the presence of small imposed acceleration, the collapse occurs due to an on-going resonance mechanism between the decreasing first natural frequency and the lower forcing frequency.

Figure 1 :
Figure 1: Geometry of the masonry arch.

Fig. 2 Figure 2 :
Fig.2shows the direct comparison between modal curvatures ''  , experimentally measured and numerically evaluated, plotted along the ratio s/s A , in which s and s A being the curvilinear abscissa and the total length along the arch centerline, respectively.

where 1  and 2  1 B
are the principal values of the effective stress tensor σ ; the bracket symbol . denotes the negative part of the number; Y  is the uniaxial compressive strength of the concrete related to the hyperbola asymptotes Y  and B are parameters governing the shape of the yield function; in the following analyses, it is always assumed  ;   is the plastic multiplier.
with G the shear modulus.If t c D D  the value of shear stress xy  undergoes a sudden jump, when 1 e of the autofrequencies deriving by the numerical computation and experimental test, respectively.In particular, in the expression(10) of the objective function only the frequency of the first mode is taken into account.Specifically, the solution is found using an optimization function in MATLAB and it is characterized by  4000 Mpa E

Figure 3 :
Figure 3: Finite element model of the masonry arch: a) mesh; b) mode shape and frequency for mode 1; mode shape and frequency for mode 2.

Table 2 :
Natural frequencies obtained with three different discretizations of the arch model.

Figure 4 :
Figure 4: Response of the masonry arch: comparison between the experimental result data and numerical results obtained for different values of the 0  .

Figure 5 :
Figure 5: Response of the masonry arch: comparison between the experimental result data and numerical results obtained for different values of the k .

Figure 8 : 7 
Figure 8: Time history of the displacement:

Table 1 :
Comparisons between the first two experimental and numerical natural frequencies.