NURBS-based kinematic limit analysis of FRP-reinforced masonry walls with out-of-plane loading

A three-dimensional (3D) general upper-bound limit analysis procedure for the determination of the collapse load of out-of-plane loaded masonry walls with Fiber Reinforced Polymer (FRP) reinforcement strips is presented. The geometry of a given FRP reinforced masonry wall of arbitrary shape is represented by its Non-Uniform Rational B-Spline (NURBS) description in the three-dimensional Euclidean space. The NURBS parameter space is partitioned by means of a number of possible fracture lines and the original reinforced wall geometry is subdivided into an initial set of rigid elements, accordingly. An upper-bound limit analysis formulation, accounting for the main characteristics of both masonry material and FRP reinforcement by means of homogenization techniques, is deduced. Internal dissipation is allowed along element edges only and the effect of vertical loads and membrane stresses is considered as well. Numerical experiments show that a good estimate of the load bearing capacity is obtained provided that the initial arrangement of yield lines is adjusted by means of a suitable Genetic Algorithm (GA).


INTRODUCTION
asonry walls represent the commonest structural element encountered in masonry constructions, largely used both as principal members [1] and non-structural internal partitions and infill panels [2]. A relevant issue in the field of structural analysis of existing buildings is the assessment of the load bearing capacity of out-of-plane loaded masonry walls. In fact, out-of-plane failures are mostly a consequence of horizontal actions, especially seismic loads, against which they generally exhibit poor resistance. Recent earthquakes in Italy (L'Aquila, 2009; Emilia 2012; Lazio 2016), demonstrated that out-of-plane failures stem from the lack of out-of-plane strength, especially for historical constructions [3][4][5]. For this reason, masonry walls in historical buildings often require external retrofitting interventions, in order to guarantee a sufficient level of out-of-plane strength. Since conventional retrofitting solutions -such as the use of steel plates and reinforced concrete plaster -are often impractical and may add undesirable mass to the existing structure, the use of Fiber Reinforced Polymer (FRP) strips is acquiring good popularity in the scientific and professional community. FRP retrofitting strips are characterized by durability, low invasiveness and good ultimate behavior. However, a limited body of work is available to date, addressing the numerical modeling of masonry walls reinforced with FRP strips [6,7]. On the one hand, finite element analyses have been proposed as a viable tool for the prediction of the load bearing capacity of FRP reinforced masonry walls [8]. However, finite element simulations are time expensive and require skilled users to correctly set the many material parameters required in order to provide reliable results. On the other hand, experimental tests carried out since the 70's on masonry walls with out-of-plane loading, have shown that collapse occurs upon formation of a welldefined pattern of linear cracks [9,10]. This evidence inspired approximated solutions relying on both the fracture line theory [11] and the yield line theory [12,13], which in fact are applications of the upper-bound theorem of limit analysis. Indeed, despite the well-established fact that masonry behaves very differently from a rigid-plastic material, limit analysis is among the most reliable tools for the assessment of the load bearing capacity of masonry walls [14]. After Heyman's assumption that masonry behaves as a no-tension material [15], several Italian scientists, among which we recall Como [14], Di Pasquale [16], Angelillo [17], Del Piero [18], have proved that the load bearing capacity of masonry structures can be estimated within the scheme of classic limit analysis under the principles of the theory of plasticity. Nevertheless, several mechanical features of masonry material are not accounted for by the no-tension model. First, masonry is a heterogeneous material with a nonisotropic behaviour, both in the elastic region and at collapse [19]. Second, masonry tensile strength, even if usually very low, is quite variable and uncertain; moreover, the assumption of infinite compression strength is often too strict. Third, experiments show that friction coefficient μ for masonry is relatively high [20]. Therefore, a non-associative flow rule must be enforced and, as a consequence, limit analysis theorems cannot be applied anymore [21][22][23]. For such reasons, a number of computational methods for the limit analysis of out-of-plane loaded masonry walls, even in the presence of FRP reinforcement, not based on the no-tension model, have been proposed in literature. For instance, we can recall several Finite Element methods (FEM) based on homogenized limit analysis [24,25]. In the present paper, we propose a new NURBS-based adaptive scheme for the homogenized kinematic limit analysis of masonry walls with out of plane loading, in the presence of FRP reinforcement. The approach, originally proposed by the Authors for both unreinforced and reinforced masonry arches and vaults (see [26][27][28][29][30][31][32][33][34]), allows to easily assess the out-ofplane failure mechanism for an FRP reinforced wall with openings of arbitrary geometry. Homogenization concepts [35] are employed to obtain out-of-plane homogenized failure mechanical parameters. NURBS (i.e. Non-Rational Uniform Bi-Spline) are the most common class of approximating basis functions employed in the field of 3D free form modeling [36]. An arbitrary FRP reinforced masonry wall can be represented by a NURBS description of its mid-surface, which can be obtained within any commercial free form modeler, and providing information about the structural thickness at each point of the surface. FRP reinforcement strips are modeled as NURBS surfaces as well. The NURBS parameter space is partitioned by means of a number of possible fracture lines and the original reinforced wall geometry is subdivided into an initial set of rigid elements, accordingly. A homogenized upper-bound limit analysis formulation, accounting for the main characteristics of both masonry material and FRP reinforcement, is deduced. Internal dissipation is allowed along element edges only and the effect of vertical loads and membrane stresses is considered as well. Since the discretization makes use of a very limited number of rigid elements, the actual failure mechanism can be found provided that the discretization is suitably adjusted by means of a Genetic Algorithm (GA) with nonstandard optimization tools [37]. In the process, possible delamination at the FRP/masonry interface is accounted for. Even though delamination is a typically brittle phenomenon, and thus, strictly speaking, it limits the applicability of limit analysis theorems, an equivalent ultimate shear strength for FRP/masonry interface is assumed in the framework of limit analysis, as suggested by the Italian norm regulating of retrofit interventions with FRP materials [38]. M With respect to existing methods, the proposed procedure exploits the properties of the NURBS representation of a masonry panel in order to build a very coarse even though exact rigid element discretization. In addition, the knowledge of the actual failure mechanism is not required in advance and the problem is governed by a relatively low number of variables thanks to the use of a very coarse discretization. Therefore, a high computational efficiency is attained. Moreover, since NURBS functions are widespread in the field of free form modeling, the proposed method could easily be integrated within existing commercial CAD software packages, which are popular among professional engineers and architects. The paper is organized as follows: Section 2 treats the NURBS discretization of a generic masonry wall with openings of arbitrary shape. In Section 3, the proposed upper bound limit analysis NURBS-based adaptive formulation is summarized, which allows to compute the collapse load for a set of given failure mechanisms also accounting for the presence of FRP strips. Finally, Section 4 illustrates the proposed approach by means of numerical simulations.

RIGID ELEMENTS
Differently from B-spline basis functions, NURBS basis functions allow the exact representation of the geometry of a wide class of curves such as circles, ellipses, and parabolas [36], as well as the surfaces that can be generated by these curves. In analogy to NURBS curves, it is possible to define a NURBS surface of degree p in the u-direction and q in the v-direction as the three-dimensional parametric surface defined as: is a lattice of control points in the Euclidean 3D space. Again, it is necessary to set a suitable set of weights , { } i j w and two knot vectors in both u and v directions. Most of the commercial three-dimensional surface modelers, such as Rhinoceros ® [39], make use of NURBS functions and their properties to generate and manipulate geometrical objects in the three-dimensional space. In the numerical simulations discussed in Section 4, both masonry walls and FRP reinforcement strips have been modeled within Rhinoceros ® as NURBS surfaces. The corresponding NURBS mathematical structure have been exported by means of the IGES (Initial Graphics Exchange Specification) standard [40], so that a suitable (rigid) discretization of the FRP reinforced masonry structure can be constructed in MATLAB ® , by exploitation of the properties of NURBS functions. In the obtained mesh, each element is a NURBS surface itself and is assumed to be rigid and a thickness is assigned to the elements by offsetting the mid-surface of the desired quantity inward and outward with respect to the normal direction. In fact, the obtained discretization is an assembly of rigid blocks. In the simplest cases, assigned a given masonry wall which can be described by a single NURBS surface, the counter-image w I of its NURBS representation in the 2D parameters space u-v is the square  [0,1] [0,1] . In general, the counter-image domain w I is bounded by NURBS two-dimensional curves, directly defined in the parameters space. The IGES format stores all the information needed to reconstruct the parameters space of a given NURBS surface. Subdividing the two-dimensional parameters space u-v, it is possible to construct a NURBS discretization of a given planar surface representing a masonry wall. More specifically, in the 2D parameter space it is possible to define a suitable lattice of nodes, which in fact are the counter-image of the nodes of the actual mesh. Suitable curves in the parameters space, A connecting the nodes define the edges of the elements in the three-dimensional Euclidean space. An appropriate connectivity matrix for the 2D counter-image of the mesh in the parameters space can be assembled.   Fig. 1(c) to the one represented in Fig. 1(d) by controlling the coordinates of the internal nodes. In the proposed approach this process will be suitably automated by means of a GA (see Section 3). Since masonry walls fail upon formation of a limited number of crack lines, the number of nodes belonging to the parameters space lattice can be relatively small [41]. In order to correctly treat openings, a simple algorithm is proposed to generate a robust mesh of the structure. This algorithm is illustrated and summarized in Fig. 2 for a façade wall with two openings containing curved edges. Fig. 2(a) shows the 3D NURBS representation of the wall, while its counter-image w I in the parameters space Fig. 2 The surface area of i E can be determined by integration as: where u S and v S are the partial derivatives of surface ( , ) u v S along u and v directions. On the other hand, the center of mass of i E may be computed as:  Integrals (3) and (4) on a generic rectangular domain can be approximated by Gauss numerical cubature coupled with a suitable isoparametric approach, see Fig. 3 and [26] for more details.

Rigid blocks kinematic limit analysis model
iven the NURBS rigid blocks discretization model of the FRP reinforced wall, an upper bound limit analysis formulation can be provided, where internal dissipation is allowed only at the element edges. This assumption has proven to be adequate for the general case of FRP reinforced masonry shells with out-of-plane loading [25]. This Section summarizes the proposed limit analysis formulation.
 is the multiplier of live loads. The Linear Programming (LP) problem associated to the upper-bound formulation of limit analysis requires the minimization of  under suitable constraints. The unknown of the problem are the set of elemental generalized velocities and plastic multipliers at the interfaces. Geometrical constraints are imposed by prescribing the values of the generalized velocities at nodes belonging to the element free edges. Therefore, the geometric constraints can be expressed in terms of the generalized velocities at the barycenter of the element containing those nodes, in the following equality form (see [26] for more details): where ,

eq geom
A is the matrix of geometric constraints and , t is the transversal tangential unit vector. On each point k P of a given interface, between two elements  E and  E , a compatibility relation must hold in the form: nn ns nt σ is the stress vector in the local reference system at k P , ( ) f σ is a suitable yield function and  λ is the (unknown) vector of plastic multipliers. In Eqn. (7),   u represents the relative velocities at points k P in the local reference system. The same quantity u in the global reference system is defined as: being  i P u and  i P u are the velocity vectors at k P on the two sides of the interface. Furthermore, the following relation must hold: where  R is a rotation matrix. The yield surface ( ) f σ can be obtained by means of a homogenization technique based on the so-called Method of Cells (MoC), first introduced in [42] for unidirectional composites and recently applied to inplane loaded masonry walls [43]. The reader is addressed to [26] FRP NURBS elements are assumed rigid, as well. Therefore, dissipation is allowed along interfaces between adjacent elements only, and is related to longitudinal stresses in the fibers direction. Again, with the aim of imposing plastic compatibility along FRP-FRP interfaces and correctly evaluating the internal dissipation rate, the mid-line of each interface has discretized into an assigned number  ( 1) F sd N of collocation points k P (see Fig. 5). On each point k P , the local reference system ( , , ) F F F n s t is introduced, F n being the outward unit vector normal to the interface, F s is the longitudinal tangential unit vector and F t the transversal tangential unit vector. Velocity jumps are 3 3  allowed along direction F n only. Typically, due to their negligible thickness, FRP strips buckle when subjected to the lowest compression stresses. Therefore, different tension and compression limit stresses are assigned, denoted as  according to [38]) and   0 FRP f respectively. The internal dissipation rate at the j-th FRP-FRP interface is computed as: where s is the FRP strip thickness, i L the FRP-FRP interface mid-line and   

I FRP k
and   

I FRP k
are plastic multipliers at k P on the interface. Delamination phenomenon depends on a large number of parameters related to materials with very different mechanical behavior (glue, brickwork, fibers). For the sake of simplicity, we adopt the provisions suggested by the Italian technical norm [38]. More precisely, we account for delamination by means of a conventional approach, consisting in suitably limiting the longitudinal tensile stress in the FRP strip -i.e. if the bond length lb is larger than optimal bond length le, the FRP strip design traction strength ffdd is computed as: whereas if l b l e f fdd is computed as In Eqn. (13) and (14), ffdd and ffdd,rid denote the design bond strength and the reduced design bond strength respectively, while E FRP is the elasticity modulus of FRP, t FRP is the thickness of FRP strip, γ fd is a safety factor equal to 1.20, γ M is the masonry partial safety factor, assumed equal to 1.0, lb is the FRP strips bond length and is the optimal bond length.
The quantity Γ Fk in Eqn. (13) represents the specific fracture energy of the FRP reinforced masonry undergoing delamination. A sound estimation of the specific fracture energy is still an open research topic, since mechanical properties of masonry are widely viariable. Nevertheless, as depicted in Fig. 6, the Italian norm prescribes the use of a bilinear τ b -slip constitutive law, allowing the assessment of the limit shear stress b f , provided that ultimate slip is known (tipically assumed equal to 0.2 mm). As clearly pointed out in many studies, a damaging material model would be the optimal choice when studying failure of FRP strips on masonry [44,45]. However, such model cannot be employed in the framework of limit analysis, which, as is well known, relies on the assumption of perfect plasticity. Nevertheless, in agreement with suggestions contained in [38], limit analysis still remains the tool of choice in order to easily estimate the bearing capacity of a given masonry structure. This fact is also reflected in the provisions provided by both the Italian Building Code [46] and the related explicative instructions [47]. Moreover, a number of studies in the literature proved that limit analysis is still capable of reliably assess FRP reinforced masonry structures, see e.g. [48,49]. Finally, FRP-masonry interfaces are NURBS surface belonging to the faces of the masonry wall. To account for dissipation along such third type of interface, a given number along the local axes is well defined. As before, with the aim of applying the limit analysis theorems, an associate flow rule on the FRP-masonry interface must hold. Again, the failure surface to be used on FRP-masonry interfaces can be linearized as being the number of linearizing planes employed). In this case, flow rule (7) is specialized to: Therefore, the total number of unknown plastic multipliers for FRP-masonry interfaces is equal to On each FRP-masonry interface i , associated to the surface i S , the internal dissipation rate is computed in the local reference system as: Moreover, the non-negativity of each plastic multipliers must be enforced by means of the additional constraint: Finally, we must impose a normality condition, requiring that the power dissipated by a unitary live load is equal to one, i.e.: Therefore, the LP problem associated to the proposed upper-bound formulation reads: under geometric constraints (6), compatibility constraints (10), (12) and (15), non-negativity of plastic multipliers constraints (17) and the normality condition (18).

GA mesh adaptation scheme
The mesh is adaptively adjusted by means of a GA, allowing the estimation of the minimum collapse multiplier, therefore determining the actual collapse mechanism. A GA is a metaheuristic algorithm for solving both constrained and unconstrained optimization problems, which mimics the process of natural selection and biological evolution [50]. The GA iteratively modifies a population of individual solutions. At each step, individuals are randomly selected from the current population and used as parents to produce next generation offspring, so that the population evolves toward an optimal solution. Here, individuals are a set of parameters defining a mesh and the objective function is the collapse load multiplier.

NUMERICAL EXAMPLES
n this Section we apply the proposed GA-NURBS limit analysis procedure to the study of the ultimate behavior of out-of-plane loaded masonry walls reinforced with FRP strips. The numerical simulations focuses on two FRP reinforced masonry wall specimens subjected to two-way bending with and without openings, fixed at the bottom edge and simply supported along vertical edges. The geometry of the walls and the disposition of FRP reinforcement strips is shown in Fig. 7. In particular, the models reproduce specimens SB01 and SB02 of a set of experimental tests on five unreinforced solid clay bricks masonry panels (labeled from SB01 and SB05) carried out in [51]. All the panels of the experimental campaign have dimensions 5615x2475x102 mm. The panels were loaded in the out-of-plane direction by airbags inflated until failure with increasing pressure p. Masonry parameters adopted in the simulations are taken from [51] and [52] (see Tab. 1); no experimental data are available from the literature concerning such panels in presence of FRP reinforcement. As for the reinforced case, two 100mm wide horizontal strips are disposed on the walls extrados, at the top and at the bottom, with the aim of preventing the opening of vertical and inclined cracks observed experimentally in the unreinforced panels. The initial NURBS discretization of the unreinforced panel is made of ten quadrangular elements, obtained by subdividing the parameters space starting from a 4x6 nodes lattice. While the four vertexes are kept fixed, the GA allows to evaluate the optimal position of the remaining fourteen free nodes, minimizing the collapse load multiplier and, thus, obtaining the actual failure mechanism.  Every node position is controlled by two parameters, with the exception of one-parameter edge-nodes and the four fixed vertex-nodes. Thus, also relying on the symmetry of the problem, the optimization problem is governed by fourteen parameters, which can be reduced to nine for to symmetry. A collapse load  2.69 p kN/m 2 has been coputed for the unreinforced case. Fig. 8(a) and Fig. 8(b) respectively depict the computed failure mechanism obtained by means of the proposed GA-NURBS approach and the homogenized FE limit analysis technique proposed in [25]. Moreover, a good agreement can be found when comparing the obtained results from the proposed GA-NURBS approach with the outcomes of both original experiments and different numerical procedures found in the literature [27]. In particular, it can be seen that the proposed GA-NURBS approach slightly overestimates the actual collapse load.  [51]. (a-c) Unreinforced and reinforced case, respectively: collapse mechanism with the GA-NURBS approach; (b-d) Unreinforced and reinforced case, respectively: collapse mechanism with the homogenized FE limit analysis proposed in [25].   [51]. (a-c) Unreinforced and reinforced case, respectively: collapse mechanism with the GA-NURBS approach; (b-d) Unreinforced and reinforced case, respectively: collapse mechanism with the homogenized FE limit analysis proposed in [25].
On the other hand, Fig. 8(c) and Fig. 8(d) respectively represent the computed failure mechanism for the reinforced case, obtained by means of the proposed GA-NURBS approach and the homogenized FE limit analysis approach proposed in [25]. A collapse load p = 4.74 kN/m 2 has been obtained for the reinforced case, which is in agreement with the collapse load found in [25], i.e. 4.58 kN/m 2 , as shown in Fig. 9 (for  0.30 b f MPa). Fig. 9 also shows the variation in the collapse load multiplier due to a variation in the FRP-masonry bond strength f b , compared with the values computed in [25]. Then, the unreinforced panel SB02 with a rectangular opening experimentally tested in [51] is analyzed. The initial NURBS discretization of the wall is made of sixteen quadrangular elements, after subdividing the parameters space from a 4x6 lattice of nodes. Again the four vertex nodes are fixed. The GA allows to estimate the optimal position of the remaining twentyone nodes, by minimization of the collapse load multiplier, thus obtaining the actual failure mechanism. Every node position is controlled by two parameters, except for the (one parameter) edge-nodes. Finally, appealing to symmetry, the number of governing parameters is reduced to fourteen. A collapse load  2. 19 p kN/m 2 has been computed with the proposed GA-NURBS scheme. Fig. 10(a) and Fig. 10(b) respectively represent the computed failure mechanism obtained through the proposed GA-NURBS approach and the homogenized FE limit analysis approach proposed in [25]. A good agreement is also found between the results from the proposed GA-NURBS approach with the outcomes from both original experiments and different numerical procedures found in the literature [27]. On the other hand, Fig. 10(c) and Fig. 10(d) represent the computed failure mechanism for the reinforced case, computed by means of the proposed GA-NURBS approach and the homogenized FE limit analysis approach proposed in [25], respectively. A collapse load p = 4.55 kN/m 2 has been computed for the reinforced case, which is in agreement with the collapse load found in [25], i.e. 4.30 kN/m 2 , as shown in Fig. 11 (for  0.30 b f MPa). Fig. 11 also shows the dependence of the collapse load multiplier on the FRP-masonry bond strength f b , compared with the values computed in [25].

CONCLUSIONS
e presented a new GA-NURBS based approach for the upper-bound limit analysis of FRP reinforced masonry walls with out-of-plane loading and arbitrary openings, in which the properties of NURBS functions are exploited to provide an efficient adaptive limit analysis scheme, allowing to easily evaluate both the collapse load multiplier and the failure mechanism of any given FRP reinforced out-of-plane loaded masonry wall, building on its 3D NURBS model. The use of NURBS functions to construct a rigid blocks discretization, allows to easily port the proposed approach into any commercial modeling environment. Moreover, differently from existing procedures implemented in commercial software packages, the GA-NURBS approach does not require an a-priori knowledge of the actual failure mechanism. By means of numerical simulations and comparisons with both experiments and numerical results from the literature, the approach has proved to be capable to accurately predict the ultimate capacity of any FRP reinforced masonry wall with out of plane loading, by using very few elements, and, therefore, maximizing computational speed. W