Thomas E. Boothby¹ and H. Sezer Atamturktur²

¹ Professor of Architectural Engineering, The Pennsylvania State University – University Park, PA, USA 16802
² Graduate Student, The Pennsylvania State University – University Park, PA USA 16802

Citation

Boothby, T., Atamturktur, S., (2007), “A Guide for the Finite Element Analysis of Historic Load Bearing Masonry Structures,” submitted to 10th North American Masonry Conference, St. Louis, Missouri, USA.

Abstract

The assessment of unreinforced masonry structures, especially in arched or vaulted forms, is difficult to undertake in practice. A structural engineer is often confronted with a standing, apparently competent, structure that seems to defy most of the rules of structural behavior as incorporated in modern building codes. The engineer must then choose between reinforcing the structure according to a modern understanding of material strength and structural behavior or trying to make sense of the behavior and anticipated strength of the structure on a more fundamental level. A manual has been developed to introduce a sympathetic structural engineer to some of the principles of unreinforced masonry and to provide some basic instructions in preparing a finite element model of complex vaulted masonry structures using widely available modern tools of structural analysis (Boothby et al. 2006). Based on the manual, this paper gives guidance on the application of finite element tools for linear and non-linear assessment of arches. Detailed instructions are given for development of geometric models, solid models and meshing and entering material properties, boundary conditions and loads for models of complex three-dimensional structures, such as domes and vaults.

Introduction

The repair and maintenance of complex masonry monuments is becoming a significant and frequent engineering problem. Engineers, when called upon to assess the load resistance under gravity, occupancy or snow loads or due to multi-hazard actions (earthquake, wind, blast), are investigating a structure built based on rules of geometric ratios or obsolete methods of structural analysis, instead of modern building codes. Often, engineers are not able to justify why and how certain heritage structures, stand based on modern structural engineering principles. A recent example is reported for the elliptical arch during the Catoctin Creek Aqueduct Restoration. The numerical model developed by the consulting engineers predicted that from inception the structure was unstable, even though the structure was intact (Biemiller 2006).

The analysis of historic load-bearing masonry structures must proceed from a fundamental understanding of the premises of the original design and a fundamental understanding of the mechanical behavior of arched and vaulted masonry systems. The analysis may be based on appropriate limit states. In this case, allowance must be made for non-linear behavior, such as the strength of an arch beyond the appearance of the first crack. Account must be taken of the particular properties of the material used in construction: of the inherent compressive strength of stone masonry, of the limited but finite tensile capacity of masonry assemblies, and of the ductility inherent in high lime content mortars used in historic construction.

Analysis of a two-dimensional arch can be completed conveniently by hand methods using the principles of limit states analysis (Heyman 1984, Heyman 1995). In this approach, an arch collapses only when it develops four hinges, that is, points where the eccentricity of the internal force is approximately half the depth of the arch. The safety of an arch can be assured by finding a statically admissible distribution of internal forces for which the thrust line is contained within the arch. The location of the thrust line can be estimated using simple calculations of plane forces under arch loading due to selfweight, superimposed dead load and live load. Alternatively, a collapse mechanism can be imposed on the arch, and the load associated with this collapse mechanism can be calculated, recognizing that this load is an upper bound on the collapse load. In the past, behavior of masonry structure was appropriately modeled by classical methods like graphic analysis or limit theory.

With the development of recent computational scheme, numerical methods have gained acceptance in the structural engineering community. The analysis of an arch can be implemented in frame analysis software by developing frame elements endowed with the geometric and material properties of the arch ring and applying appropriate load and support conditions. The output is examined for the relationship between axial force and bending moment. Simple results of the eccentricity of internal forces can be inferred, or more complex determinations of the required compressive capacity of the masonry can be made, based on elasticity or plasticity assumptions for the behavior of the masonry.

Finite element (FE) method, a widely established computerized numerical method, has been one of the tools employed extensively for the analysis of large scale masonry buildings. Until recently, however, their accuracy has not been explored. In the process of FE analysis, there are two main sources for uncertainty: (a) the determination of true physical characteristics of the actual structure, (b) the utilization of the FE software options. The first group of uncertainties challenges the determination of the actual physical properties of the structure–the assessment of the geometry, boundary conditions and material properties—which are particularly challenging when dealing with historic masonry vaults or domes. The latter group includes the user preferences during the analytical model development: level of geometric simplification, selection of element types, decisions on mesh refinement, and choice of analysis type. The present paper intends to address these two sources of uncertainty and to guide the process of accurate analytical model development of complex vaulted masonry structures.

Analytical Modeling of Masonry Vaulted Structures

The process of developing an analytical model of any three-dimensional structure begins with the reproduction of its geometry and followed by the creation of the solid model in a computer platform. Once the appropriate element types are selected within the finite element software, the solid model is broken down into finite pieces through a systematic procedure of meshing. In this phase, the material input is also entered to the components of the model. With the application of loading and boundary conditions, the initial model is obtained. The solution of initial model must be then compared against known information about the structure. If the correlation is satisfactory the analytical model and its solutions are assumed to be in agreement; otherwise the uncertain parameters in the model are calibrated, so that the analytical model reproduces the reference data. When a finite element model with an acceptable predictive ability is obtained, the sensitivity of the model to the variations of the parameters is assessed through a sensitivity analysis. The following paragraphs discuss the analytical model development procedure with particular attention to large-scale masonry vaulted structures.

Geometric Model

To develop an accurate structural model, it is necessary to determine the physical dimensions of the system. Building survey measurements or available construction drawings can be used to locate reference points that define the curved geometry of a three-dimensional vaulted structure. The location of the reference points in a three-dimensional coordinate system must be entered manually to the program. Connections between these keypoints, depending on the software, can take the form of straight lines, arcs or splined curves.

In general, the analysis of load-bearing masonry structures does not require an extraordinary level of detail in the determination of the physical dimensions of the structure. The thickness of piers or walls supporting the structure should be determined at the base and at any important (>10-20 %) changes in thickness. On a load-bearing arch, a significant proportion of required analytical work can be completed with the dimensions of span, rise and thickness (Figure 1), along with a knowledge of the general arch profile, whether segmental (circular), three-centered, or four centered. The profile of load-bearing arches, unless seriously distorted due to loading or settlement effects, can be adequately characterized by span/rise and 1/4 point rise (Figure 2). Where filling or diaphragm walls are present, rudimentary measurements of the height and thickness of the walls are also necessary.

The dimensions of rib vaults can be adequately determined on the basis of the bordering arches of the vault (Figure 3). The vault ribs generally have a constant thickness; so any means of determining the dimensions of a rib, whether by remote measurement, consulting drawings, or direct measurement, is sufficient. When these have been determined, the actual shape of the vault between ribs can be generated with a few additional measurements taken from the curved web surface. Groin vaults require the measurement of arch parameters in the direction of the two intersecting barrel vaults. Domes are almost always sections of a sphere, and a satisfactory model of a dome can be created using a set of span, rise and thickness measurements. Vaults and domes are often provided with filling from the supports up to the haunches, and the height of the filling also needs to be measured. The thickness of webbing may be difficult to determine unless there are holes in the vault. The use of impact-echo or other non-invasive tools for the determination of thickness are available (Sansalone 1997, Schuller 1997).

In advance of the modeling effort, it is important to determine the level of simplification that will preserve realistic results. Complex molding profiles or multiple colonnettes on the ribs are practically impossible to model accurately, and the excessive effort spent on attempting to model these shapes is not justified by the increased accuracy of the results obtained. In general, it is preferable to reduce complex cross-sectional shapes to a simpler rectangular shape of approximately the same cross-sectional area and moment of inertia. An example of such reduction for the cross section of the piers of the Cathedral of Saint John the Divine is illustrated in Figure 4.

Solid Model

Depending on the preferred FE software, a solid model can be created with modeling utilities in the software interface (Figure 5a) or with alternative computer platforms, such as CAD programs (Figure 5b). There are significant interactions between the geometric model creation, solid model development and the choice of elements for use in structural analysis. For instance, the preference of frame, solid or shell elements for a component will determine whether it is modeled as a line, an area or a volume object, and thus the dimension and geometry of the solid model. Deciding on the element types for the structural components early in the process helps to avoid unnecessary modifications. In the next section, guidelines are given for the element type selections for typical structural components of a masonry vaulted system.

Element Types

Generally speaking, there are three types of elements available in most analysis programs: line, area and solid elements. The preference of one over another is dependent on the physical characteristics and scale of the members.

Line elements represent a beam or bar by a single line, usually along the centroid of the member. Depending on the abilities of the analyses software, the element may be straight or curved. A spatial frame element has up to six degrees of freedom–three displacement and three rotation–at each end. Shell elements are three- or four-cornered elements, with nodes at each corner and are well suited to model singly curved or doubly curved surfaces of constant or varying thickness. More refined shell elements with mid-side nodes are also available in advanced computer programs. The nodes of shell elements each have six degrees of freedom: three translational and three rotational. Alternatively, solid elements have three dimensions and hence a greater number of nodes than frame or shell elements. Tetrahedral solid elements with four nodes and so-called brick elements with a minimum of eight nodes are available in most FE software. Refined versions of these elements, including one mid-side node, can also be found in advanced programs. For most solid elements, nodal degrees of freedom include translations only.

Vault webbing can successfully be modeled by shell elements, as they have reduced degrees of freedom compared to solid elements and include all of the load effects that characterize vault webbing, including the interactions between membrane forces and bending moments (Figure 6). Modeling a thin vault web with solid elements is impractical, because the requirement of at least three solid elements through the thickness to preserve a reasonably accurate stress field and the limitations on the aspect ratio (the ratio between key dimensions, such as length, width, thickness) of elements result in an excessively fine mesh. Vault ribs, on the other hand, can be effectively modeled as solid elements, but require at least a nine-element cross section. For simplified models, beam elements can also be used to model vault ribs. Walls may be modeled using either solid elements or shell elements, depending on the scale and the situation. Because of the irregular shape of surcharge volumes, tetrahedral solid elements are virtually necessary for these features. Although, tetrahedral elements generally are inaccurate for stress computations, the stresses in the surcharge do not significantly influence the assessment of the structure.

Meshing

On completion of the solid model and selection of the element types, the solid model is broken down into individual elements according to a systematic procedure known as meshing the model. The shape and size of the elements impacts the solution. A mesh that is too coarse can produce inaccurate solutions, while a mesh that is too fine will result in problems with program limits on the number of nodes or elements or will result in excessive run times. The aspect ratio of the elements must also remain within reasonable limits: an aspect ratio of less than two is desirable, while a ratio of greater than four is unacceptable. Similarly, the angles of corners can be neither too acute nor too obtuse (ANSYS 2003). Some programs also give warning prompts when aspect ratio or corner angle limits are exceeded.

When solid elements are used, the requirements of the aspect ratio determine the required fineness of the mesh. Beyond this, experience indicates that, at the scale of most monumental masonry structures, an element dimension of approximately 10-30 cm results in an appropriate level of mesh refinement (Atamturktur 2006).

Particularly difficult problems in meshing are present for rib vaults at the intersections of vault webbing and ribs and at the springing of vaults. Closely spaced multiple ribs result in a very narrow vault web. At these locations, it is impractical to respect aspect ratio limitations, given that the shell elements need to share nodes with the solid arch rib elements. It has been found that omitting the vault webbing in the very narrow parts is a convenient way to manage this modeling problem (Erdogmus 2004).

Material Properties

The real behavior of masonry is non-linear, in the sense that the stress-strain law, even for elastic behavior in compression, may show softening. It is by far simplest to work with linearly elastic, homogeneous, isotropic material properties and significant progress in the assessment of a masonry structure can be made using this simple form of constitutive law (Atamturktur 2006). In that case, the principal material properties required for structural modeling and assessment are the compressive strength, tensile strength, density and modulus of elasticity of the masonry assembly. Poisson's ratio input is also required for analysis programs, but is of secondary importance. A generic value of 0.10 to 0.20 for Poisson's ratio is satisfactory.

The principal material properties of a masonry assembly can be obtained by consulting reference material, by conducting tests on specimens of material extracted from the structure, or by methods of non-destructive evaluation. Although it is tempting to think that the most reliable procedure for obtaining the mechanical properties of the masonry units is to extract a specimen from the structure, especially when tension-controlled properties are considered, obtaining reliable values is practically impossible. However, compressive strength may be relatively unaffected by the breaking of mortar joints during specimen extraction, and it is possible that compressive strength and compressive modulus of elasticity tests may be able to be conducted on such specimens.

Loading

Loads can be applied to nodes as thermal, inertial, concentrated or distributed forces in the form of nodal displacement or nodal force. The arrangement of loading conditions depends on analysis type. Static analysis, typically used to examine the structural behavior under gravity or service loading conditions, requires definition of the gravitational acceleration as well as existing external loads. In modal analysis, the results are independent from force input, and the only possible loading condition is a zero displacement constraint, in other words, boundary condition definitions. However, transient and harmonic analysis options, which may be used to simulate the seismic, wind or blast events, require a user-defined time dependent loading in addition to the provided boundary conditions.

Boundary Conditions

Boundary conditions or support conditions have significant influence on the computed results. In general, it is expedient to model only a portion of the structure, and it is very difficult to assign correct boundary conditions as that modeled portion of the structure is supported elastically. The boundary conditions, for a masonry building, are dependent on physical properties and configuration of the material, rather than intentionally designed 'pins' or 'points of fixity,' as in a steel structure. In many cases the choices are imposed by the analyst’s intuition and the technology and capacity of the FE software.

Because boundary conditions are applicable to the degrees of freedom at the nodes, rather than to the elements, their prescription is related to the selection of element types. In most commercially available FE software, physical constraints are invoked by zero displacements and/or rotations at the user-defined nodes. For situations of partial restraints, an elastic foundation, which is usually simulated by a series of springs, can also be adapted.

In the subsequent discussion, the longitudinal axis of the nave will be described as the x-axis, the transverse axis of the nave will be designated the y-axis, and the vertical axis will be the z- axis. The recommended boundary conditions vary, depending on the nature of the study being conducted. If the purpose of the study is to determine the thrusts of the vaults on the supporting structure, then the vault is modeled from the springing point upwards. In such cases, the boundary conditions include restraint in the x, y, and z directions at the pier tops. The edges of the vault are subject to symmetry boundary conditions. Generally, at least, the nave wall above the level of the springing of the vault should be modeled. This wall is subject to restraint in the x direction for the full height of any buttressing, including wall buttresses. The wall should also be laterally restrained at the bearing point of any roof trusses (Figure 7) (Erdogmus 2004). These guidelines are developed based on linearly elastic behavior assumption and subject to error when time dependent support settlements or inelastic distortions occur.

In another type of model, the full height of a bay, or 1/4 bay, may be modeled. In this case, the base of the pier and the exterior wall is restrained in the x, y, and z directions. If a 1/4 bay is modeled, symmetry boundary conditions are imposed along the centerlines of the bay along the longitudinal axis, x-translation, and rotation about the x and y-axes, and along the transverse axis, y-translation and rotation about the x and y-axes. When a full bay is modeled, similar conditions are applied to represent the effect of the adjacent bays.

Solution

Depending on the goals of the study, two main solution types are available in most FE software: static and dynamic analysis. Static analysis slows the forces being applied so that no inertia or dissipative resistance is activated, while dynamic analysis examines how the external loads are balanced with inertial, dissipative and elastic resistance of the system within the discretized time domain. Both of the analysis approaches can be solved for linear or nonlinear behavior range.

A model based on linearly elastic principles enables significant progress in the assessment of a masonry structure without the convergence problems associated with non-linear analysis. The behavior of a structure under self-weight and service load can be determined quite accurately as long as the tensile stresses in the structure are low. At the very least, the initial analysis for any masonry structure should include a linearly elastic analysis to determine the presence and magnitude of tensile stresses in the masonry. The assessment of the need for a non-linear analysis is made on the basis of the initial linear analysis. Factors that would indicate the need for a non-linear analysis include tensile stresses exceeding the rupture strength of the masonry in a large zone, the presence of large foundation settlement or the presence of large cracks of unexplained origin in the structure.

Some widely distributed frame analysis programs allow appropriate input for the non-linear analysis of plane arches. The softening behavior under compression and cracking behavior under tension can be input in the form of stress-strain laws, which can be used to define hinges, based on axial force-moment interaction. These hinges can be located at appropriate intervals within the arch, and can assess the behavior of an arch beyond the elastic limit. Non-linear methods of analysis should be undertaken with caution, only when necessary, and only after completion of an accompanying linear analysis.

Calibration of Model and Sensitivity Analysis

The analysis process must continue beyond the solution of the initial model. The analytical estimates must always be verified with known information on the structure, such as the results of in situ testing, or the location of existing damage. Following the initial solution, the model may be refined, the boundary conditions may be altered, material properties may be adjusted and the model is re-run. This process is generally known as tuning or calibration of the model. It is crucial to start with physically reasonable boundary conditions and material properties in the initial FE model.

One widely used method of extracting knowledge from the existing structure is the techniques of experimental modal analysis (EMA). EMA, a form of vibration testing, is concerned with the identification of dynamic parameters of a system, which are used a reference data during the calibration of the FE model. A brief theoretical background to the method along with its application to several monumental masonry buildings is included in the manual (Boothby et al. 2006): the main nave vaults of Cathedral of St. John the Divine (NY) and National Cathedral (DC), and the reading room domes of State Education Building (NY). Much more detail on the use of this method for model validation is available in Atamturktur (2006) and Erdogmus (2004).

It is also prudent to identify which parameters have the greatest influence on the outcome of the model, by varying all of the major parameters, such as material strength or stiffness, and examining their influence on the results. This process is known as sensitivity analysis, and its application is also detailed within the manual.

Conclusions

The assessment of three-dimensional structures, such as domes and vaults, requires a level of sophistication which may be addressed by the implementation of the FE method. The process of building a model involves the steps of the development of a solid model, selection of element types, meshing, and application of boundary conditions, application of loads, solution, and postprocessing. Several cycles of review of the analysis results, updating of the model, and reanalysis are generally required to develop a satisfactory analytical model. Displacement results are significantly more reliable than stress results and separate substructure modeling may be necessary in order to determine more accurate stress results. Simplification of the complex geometry of an actual structure is generally necessary to make the size of the model manageable. In general, walls, buttresses and vault ribs are modeled using solid brick elements. At least three elements through the thickness are necessary in order to capture bending moments and rotations. Areas of complex geometry for which accurate stress results are not necessary, such as rubble filling, may be modeled with tetrahedral solid elements. Shell elements are appropriate for the modeling of masonry vaults and domes, and simplify many of the modeling issues by working explicitly with internal force and moment resultants.

Varied techniques are available for the modeling of actual structures. The geometry of the structure may be input either through the program interface, or by importing a model from a computer-aided drafting program. Decisions on material properties, element types, meshing, boundary conditions and loading conditions result from careful review of the structure in the field and of the literature on the particular construction method employed. The results of the analysis need to be interpreted in conjunction with a detailed inspection of the condition of the structure. If possible, a number of additional validation methods are available to increase the confidence in the quality of the structural model. In any event, a sensitivity analysis of the model should be undertaken.

A responsible assessment of a historic load-bearing masonry structure can be undertaken using modern engineering tools, combined with an understanding of the existing structure, and a strong sympathy for the preservation of structures from the past.

References

1. ANSYS Inc. 2003: ANSYS v. 8.0 Elements Reference.
2. Atamturktur 2006: Atamturktur, S., “Structural Assessment of Guastavino Domes,” M.S. Thesis, The Pennsylvania State University, May 2006.
3. Biemiller 2006: Biemiller L., “Catoctin Creek Aqueduct Restoration”, Society for Industrial Archeology, vol. 35, no. 2 pp. 8-9.
4. Boothby 2006: Boothby, T. E., S. Atamturktur, E. Erdogmus, “Manual for the Assessment of Load Bearing Unreinforced Masonry Structures” (www.arche.psu.edu/masonry-assessment).
5. Erdogmus 2004: Erdogmus, E., “Structural Appraisal of the Florentine Gothic Construction System,” PhD Thesis, The Pennsylvania State University.
6. Heyman 1982: Heyman, J., The Masonry Arch, Chichester, Ellis Horwood, 1982.
7. Heyman 1995: Heyman, J., The Stone Skeleton, Cambridge, Cambridge University Press, 1995.
8. Sansalone 1997: Sansalone, M., “Impact-Echo: The Complete Story,” ACI Structural Journal, 94(6):777-786, November-December 1997.