Analysis of the Load-Bearing Capacity of a Reinforced Concrete Flat Slab-Column System Using the Finite Element Method

Introduction

Reinforced concrete flat slab–column systems have occupied a prominent position in structural engineering practice for well over a century, and their prevalence in contemporary construction continues to grow. The structural concept, in which a continuous slab of uniform or variable thickness transfers loads directly to supporting columns without the intermediary of beams, confers a range of practical and economic advantages that have made it the preferred solution in a wide variety of building typologies ^[1]. Office buildings, car parks, residential towers, and industrial facilities routinely employ flat slab construction owing to the reduced storey height it permits, the simplified formwork geometry it enables, and the unobstructed ceiling soffit it provides for the routing of mechanical and electrical services ^[2]. As urbanisation intensifies and the demand for efficient use of structural depth becomes ever more pressing, these advantages are likely to ensure that flat slab systems remain a dominant structural form in the built environment for the foreseeable future.

The structural behaviour of flat slab–column systems is, however, characterised by a complexity that stands in considerable contrast to the apparent simplicity of their geometry. In the vicinity of supporting columns, the slab is subjected to a highly localised concentration of shear forces and bending moments that gives rise to a failure mode known as punching shear ^[3]. This mode of failure, in which a truncated cone or pyramid of concrete is effectively punched through the slab thickness along a critical perimeter surrounding the column, is of particular concern because of its brittle and often sudden character. Unlike flexural failure, which typically provides visible warning through the progressive development of cracking and large deflections, punching shear failure may occur without significant prior deformation, leaving little opportunity for redistribution of internal forces or for occupants to seek safety ^[4]. The consequences of such a failure can be catastrophic: progressive collapse, in which the failure of a single slab–column connection initiates a chain reaction of adjacent failures, has been responsible for a number of notable structural disasters, and the avoidance of this mechanism remains a central concern in the design of flat slab structures ^[5].

In response to this complexity, the engineering community has developed a substantial body of analytical and experimental knowledge pertaining to the punching shear resistance of reinforced concrete slabs. Empirical approaches, derived from the statistical analysis of laboratory test results, have been codified in national and international standards and continue to form the basis of routine design practice ^[6]. Eurocode 2, which provides the principal normative framework for the design of concrete structures throughout Europe, prescribes a procedure for the verification of punching shear resistance that is founded on the concept of a control perimeter located at a fixed distance from the column face, along which an average shear stress is computed and compared against a permissible value dependent on the flexural reinforcement ratio and the concrete compressive strength ^[7]. This approach, while computationally tractable and well-validated against a broad experimental database, necessarily involves simplifying assumptions regarding the distribution of shear stresses and the influence of two-dimensional stress states that limit its capacity to capture the full complexity of slab behaviour under loading conditions that depart significantly from those represented in the underlying test data ^[8].

The limitations of code-based analytical methods have stimulated growing interest in the application of numerical techniques to the analysis of flat slab structures, and the Finite Element Method has emerged as the most widely employed tool for this purpose ^[9]. The Finite Element Method permits the discretisation of a structural continuum into a mesh of elements whose mechanical response is governed by constitutive relationships that can, in principle, capture the nonlinear behaviour of reinforced concrete with considerable fidelity. The progressive degradation of concrete stiffness under compressive loading, the initiation and propagation of tensile cracking, the yielding of reinforcing steel, and the interaction between these phenomena can all be incorporated within a finite element model, enabling the simulation of structural response from the linear elastic range through to ultimate failure ^[10]. The spatial distribution of stress resultants, which is accessible only in an averaged or approximate sense through analytical methods, can be obtained with a resolution limited principally by the refinement of the finite element mesh, providing a level of detail that is invaluable for the identification of critical zones and the assessment of failure mechanisms ^[11].

Despite the considerable advances that have been achieved in the finite element modelling of reinforced concrete structures, the application of this technique to flat slab–column systems presents challenges that are not encountered in simpler structural configurations. The pronounced three-dimensionality of the stress state in the slab–column connection region, the sensitivity of punching shear behaviour to the choice of constitutive model and element type, and the difficulty of defining boundary conditions that accurately represent the restraint conditions present in an actual structure all require careful consideration ^[12]. The validation of finite element models against experimental data is an essential step in establishing their reliability, and the interpretation of numerical results must be informed by a thorough understanding of the theoretical framework within which the computations are performed ^[13]. These requirements impose demands upon the analyst that go beyond the mechanical operation of a commercial finite element software package, and they constitute a significant motivation for academic investigations of the type presented in this thesis.

The present thesis is directed towards the analysis of the load-bearing capacity of a reinforced concrete flat slab–column system using the Finite Element Method, with systematic reference to the analytical provisions of Eurocode 2. The investigation is motivated by the structural importance of flat slab systems, by the limitations of purely analytical approaches in capturing the complexity of their behaviour, and by the potential of finite element analysis to provide a more detailed and rigorous assessment of load-bearing capacity and failure mechanisms. The scope of the work encompasses the theoretical principles governing the structural and mechanical behaviour of flat slab–column connections, the formulation and implementation of a finite element model appropriate for the problem under consideration, the analysis of structural response under both linear and nonlinear conditions, and a parametric study of the influence of selected geometric and material parameters on punching shear resistance ^[14].

Three principal research objectives are pursued in the course of the investigation. The first is to develop a finite element model of a reinforced concrete flat slab–column connection that incorporates appropriate constitutive representations of the nonlinear behaviour of both concrete and reinforcing steel, and to validate this model against experimental results and Eurocode 2 predictions. The second is to apply the validated model to the determination of the ultimate load-bearing capacity of the system, and to characterise the failure mechanism through examination of the distribution of principal stresses, plastic strains, and crack patterns in the vicinity of the column. The third is to conduct a systematic parametric study in which the slab thickness, the flexural reinforcement ratio, and the column dimensions are varied independently, with the objective of quantifying their respective contributions to punching shear resistance and of assessing the degree of agreement between finite element predictions and the Eurocode 2 design procedure across the range of parameters considered ^[15].

The investigation draws upon a well-established body of literature relating to the punching shear behaviour of reinforced concrete flat slabs, the finite element modelling of concrete structures, and the provisions of Eurocode 2. Foundational experimental studies conducted over several decades have established the principal parameters governing punching shear resistance and have provided the empirical basis for code formulations ^[16]. Theoretical contributions, including the critical shear crack theory and related analytical models, have offered mechanistic interpretations of the punching failure process that complement the empirical approach and have been influential in the development of more recent code provisions ^[17]. The application of finite element analysis to punching shear problems has been the subject of numerous published investigations, and the collective findings of this literature inform both the modelling strategy adopted in the present work and the interpretation of its results ^[18].

The thesis is structured into three chapters, preceded by this introduction and followed by a concluding chapter. The first chapter provides the theoretical foundations necessary for the subsequent analytical and numerical work. The structural concept of the flat slab system is introduced, and its mechanical behaviour under gravity loading is described with reference to the distribution of bending moments and shear forces in the slab. The punching shear mechanism is examined in detail, and the principal factors governing punching shear resistance are identified and discussed. The material properties of reinforced concrete relevant to the analysis are summarised, with particular attention to the nonlinear stress–strain response of concrete in both compression and tension. The provisions of Eurocode 2 pertaining to the design of flat slabs for punching shear are presented and critically discussed ^[19].

The second chapter is devoted to the Finite Element Method and its application to reinforced concrete flat slab–column systems. The fundamental principles of the method are reviewed, and the element types and mesh configurations employed in the present analysis are described and justified. The constitutive models adopted for concrete and reinforcing steel are presented in detail, and the boundary conditions applied to the finite element model are explained with reference to the physical conditions they are intended to represent. The validation of the model is addressed through comparison of finite element predictions with experimental data drawn from the literature and with the results of Eurocode 2 calculations, and the reliability of the model for the intended range of applications is assessed ^[20].

The third chapter presents the results of the finite element analyses conducted in the course of the investigation. The response of the slab–column system under linear elastic loading is characterised through examination of the distributions of bending moment, shear force, and stress resultants, and the regions of maximum stress concentration are identified. The nonlinear analysis is described, and the load–deflection response of the system is presented from the initial loading stage through to ultimate failure. The failure mechanism is examined through analysis of the distribution of principal stresses and plastic strains at successive load levels, and the ultimate load predicted by the finite element model is compared with the punching shear resistance calculated in accordance with Eurocode 2. The parametric study is then presented, and the sensitivity of punching shear resistance to variations in slab thickness, reinforcement ratio, and column dimensions is quantified and discussed ^[21]. The concluding chapter synthesises the principal findings of the investigation, draws conclusions regarding the relative merits of finite element and code-based approaches to the assessment of punching shear resistance, and identifies directions for future research ^[22].

It is anticipated that the results of this investigation will contribute to the understanding of flat slab–column behaviour and will demonstrate the value of finite element analysis as a complement to the analytical procedures prescribed in design standards. By providing a detailed characterisation of the stress state and failure mechanism in the slab–column connection region, and by quantifying the influence of key geometric and material parameters on punching shear resistance, the present work offers insights that are relevant both to the academic study of reinforced concrete structures and to the practice of structural engineering.

Chapter 1: Theoretical Foundations of Reinforced Concrete Flat Slab–Column Systems

1.1. Structural Concept and Applications of Flat Slab Systems

The flat slab–column structural system is defined as a reinforced concrete floor or roof system in which a slab of uniform or variable thickness is supported directly on columns, without the intermediary of beams or girders. In this configuration, loads are transferred from the slab to the columns through a two-way bending mechanism, fundamentally distinguishing the system from conventional one-way or two-way beam-supported slabs. The absence of primary beams constitutes the defining characteristic of the typology and gives rise simultaneously to its principal constructional advantages and its most critical structural vulnerability. A clear distinction is customarily drawn between the flat plate, which employs a slab of uniform thickness throughout without any geometric enhancement at the column–slab interface, and the flat slab proper, which may incorporate drop panels, column capitals, or both in order to locally increase slab depth and thereby improve shear resistance in the critical connection region ^{[4, s. 3]}.

Several sub-variants of the flat slab typology are recognised in engineering practice and codified in relevant design standards. The flat plate offers the simplest geometry: a constant-thickness slab without any geometric modification at column heads, maximising architectural clarity at the cost of reduced punching shear resistance. The flat slab with drop panels introduces a thickened region of the slab centred over each column, the dimensions of which are prescribed by Eurocode 2 as at least one-sixth of the shorter span in each direction. The flat slab with column capitals widens the effective support area by flaring the column head into the slab, thereby reducing the shear stress concentration at the critical control perimeter. Additionally, waffle (two-way ribbed) slabs represent a closely related typology in which material is selectively removed from the tension zone, reducing self-weight while preserving flexural stiffness [4, s. 25–26].

The historical development of flat slab construction spans more than a century. Early patents and structural applications appeared in the United States at the beginning of the twentieth century, with C. A. P. Turner credited with one of the first documented uses of a beam-free reinforced concrete floor system around 1906. Concurrently, Robert Maillart in Switzerland was developing innovative reinforced concrete forms that similarly dispensed with secondary beams. These pioneering applications revealed both the efficiency and the structural risks of the system, prompting early experimental research into punching shear failure ^[23]. Subsequent decades witnessed the progressive incorporation of flat slab provisions into national and international codes, culminating in the comprehensive treatment provided by EN 1992-1-1:2004 and the fib Model Code 2010 ^[8].

The principal structural advantages that have rendered the flat slab system prevalent in contemporary construction practice may be summarised as follows. The elimination of beams reduces floor-to-floor heights, which is of considerable economic significance in multi-storey construction where each saved decimetre of storey height reduces the total building height, curtain wall area, and service core length. The unobstructed soffit simplifies the routing of mechanical, electrical, and hydraulic services, permitting greater flexibility in the coordination of building systems. The simplified formwork geometry, consisting of essentially flat soffit panels, reduces on-site labour costs and formwork cycle times in comparison with ribbed or beam-slab systems. Furthermore, the flexibility of the column grid layout – unconstrained by beam depths and directions – allows greater architectural freedom in floor plan organisation ^[3]. The principal disadvantage, which motivates the entire analytical programme of this thesis, is the inherent vulnerability of the column–slab connection zone to punching shear failure: a brittle, localised three-dimensional failure mechanism that may occur suddenly and without adequate warning under concentrated gravity loads ^[1].

The fields of application of flat slab construction are broad and continue to expand. Multi-storey car parks constitute one of the most common applications, where the clear soffit facilitates ventilation, lighting, and fire suppression installation, and the flexibility of column spacing accommodates standard parking geometries. Office buildings and commercial premises benefit from the open floor plans enabled by the typology. Residential apartment blocks, particularly in Central and Eastern Europe, have adopted flat slabs extensively, driven by economic considerations and the efficiency of repetitive floor construction. Industrial floor slabs subjected to heavy distributed or concentrated loads represent a further application domain, where punching shear from racking or machinery loads is a primary design concern ^[2]. Typical geometric parameters governing the proportioning of flat slabs include span-to-effective-depth ratios in the range of 30 to 35 for characteristic imposed loads up to 5 kN/m², slab thicknesses commonly ranging from 180 mm to 300 mm for spans of 6 to 9 m, and column grid dimensions dictated by functional requirements rather than structural necessity. The normative framework governing the design of these systems in European practice is established by EN 1992-1-1:2004 (Eurocode 2) ^[8], supplemented by the fib Model Code 2010 and, for international reference, by ACI 318.

1.2. Mechanical Behaviour and Internal Force Distribution

The structural behaviour of a flat slab under gravity loading is characterised by a two-way bending mechanism in which the applied loads are distributed to supporting columns along both principal directions of the floor plane. Unlike a one-way slab or beam, in which the load path is essentially uniaxial, the flat slab engages the full three-dimensional stiffness of the plate, with bending moments, shear forces, and twisting moments simultaneously present at every section. This bidirectional action is the source of the system's efficiency in spanning between columns without intermediate beams, and it is also the origin of the complex stress state that develops in the critical column–slab connection zone. The complexity of the internal force field necessitates analytical methods that go beyond simple one-way beam idealisation, as recognised by both simplified code procedures and more advanced finite element approaches ^{[4, s. 29]}.

Two classical simplified analytical methods are employed in engineering practice to estimate the distribution of bending moments in flat slabs: the equivalent frame method and Hillerborg's strip method. In the equivalent frame method, the slab is idealised as a series of equivalent beam-column frames in each principal direction, with the slab width divided between column strips and middle strips in defined proportions. The bending moments computed for each equivalent frame are then distributed between column and middle strips in accordance with tabulated coefficients specified in Eurocode 2, Annex I ^[8]. In the strip method, due to Hillerborg, the slab is conceptually divided into a statically determinate set of strips spanning in the two principal directions, with the designer free to choose the distribution of load between the two strip directions subject only to the condition that equilibrium is satisfied. Both methods are predicated on elastic behaviour and yield reasonable estimates of moment distributions under service loading conditions ^[24].

The column strip, defined in EN 1992-1-1 as extending across a width equal to half the span on each side of the column centreline (clause 5.3.2.1), concentrates a disproportionately large fraction of the total bending moment relative to the middle strip. In the direction perpendicular to the columns, hogging (negative) moments develop over the column supports, while sagging (positive) moments are present at midspan. The proportion of total static moment assigned to the column strip over a support is typically in the range of 60 to 80 percent, depending on the aspect ratio of the column grid and the relative stiffness of column and slab. This concentration of moment in the column strip gives rise to correspondingly concentrated tensile stresses in the top reinforcement immediately above the column, which must be reflected in the detailing provisions of clause 9.4.1 of Eurocode 2, requiring that at least two-thirds of the column strip top reinforcement be positioned within a width equal to half the column strip width centred on the column ^[8].

Twisting moments, arising from the biaxial nature of the plate's bending state, are of particular significance in corner and edge column regions. In an interior panel, twisting moments cancel by symmetry at the column centreline, but along free edges and at corners, significant torsional shear stresses develop. These torsional effects are often neglected in simplified design methods but are properly captured by finite element analysis, in which the full stress resultant tensor is evaluated at every integration point. The interaction between bending moments and twisting moments in corner regions can produce principal tension stresses oriented at angles other than the principal Cartesian directions, which is relevant both to the placement of reinforcement and to the interpretation of crack patterns observed in laboratory specimens ^[25].

The distribution of vertical shear forces in a flat slab exhibits a highly non-uniform character, with shear intensities increasing markedly in the immediate vicinity of column supports. The shear force per unit width on a section parallel to a column face is not constant along the perimeter of the column but varies with the distance from the column, the loading configuration, and the transfer of unbalanced moments between slab and column. The concept of a critical perimeter, at which the average shear stress is evaluated for punching verification, arises directly from this non-uniform distribution. At large deformations, tensile membrane action in the reinforcement and compressive membrane action in the concrete may develop, the latter being particularly pronounced in laterally restrained slabs and contributing a beneficial increase in load-bearing capacity beyond that predicted by linear elastic theory ^[26].

Yield line theory, developed by Johansen, provides a plasticity-based upper-bound estimate of the ultimate load of flat slabs by identifying the mechanism of collapse through a pattern of yield lines along which the bending moment equals the plastic moment resistance. The method is particularly valuable for irregular slab geometries or unusual loading configurations where tabulated moment coefficients are not applicable. A critical limitation of yield line theory in the context of flat slab design must, however, be explicitly acknowledged: the method assumes that the ultimate load is governed by a flexural collapse mechanism and makes no provision for the possibility that punching shear failure at a column–slab connection may occur at a load substantially lower than the yield line collapse load. In practice, punching shear frequently constitutes the governing failure mode in slabs with high reinforcement ratios, and the application of yield line theory without an independent punching verification would be unconservative and potentially unsafe ^[1].

1.3. Punching Shear: Mechanism, Failure Modes, and Critical Parameters

Punching shear is a localised three-dimensional failure mechanism characterised by the penetration of a column or concentrated load through the surrounding slab, forming a truncated cone or inverted pyramid-shaped failure surface that propagates from the column face through the slab thickness. The inclination of the failure surface to the slab plane is typically in the range of 25° to 45°, depending on the reinforcement ratio, concrete strength, and slab effective depth, and may vary around the perimeter of the column in cases where the loading is eccentric or the column geometry is non-circular. The failure surface separates the concrete in the immediate vicinity of the column from the remainder of the slab, and its formation is associated with the development of a critical shear crack that propagates from the compression zone of the slab downward toward the flexural tension reinforcement ^[1]. The catastrophic and brittle nature of punching failure distinguishes it sharply from flexural failure and makes it the primary concern in the design of flat slab–column systems.

The progressive development of punching failure may be described as a sequence of identifiable stages. In the initial loading stage, the slab behaves essentially as a linear elastic plate, with bending moments distributed biaxially between column and middle strips. As loading increases beyond the cracking moment of the tensile zone, flexural cracks initiate at the top surface of the slab over the column, propagating radially and tangentially in a characteristic pattern. With further load increase, these cracks develop into inclined shear cracks in the critical region at approximately one effective depth from the column face, growing from the tension reinforcement toward the compression zone. The critical shear crack, which governs the ultimate failure, is typically inclined at an angle that reflects the local state of principal stresses, and its width is the kinematic parameter central to the Critical Shear Crack Theory (CSCT) developed by Muttoni ^[27]. Final failure occurs when the compressed concrete strut above the critical crack loses its residual load-bearing capacity, either by crushing or by shear-compression failure.

A fundamental distinction is drawn in the literature between brittle punching failure, in which the column penetrates through the slab suddenly and without significant prior ductility or moment redistribution, and ductile flexural failure, in which a complete yield line mechanism develops before any localised shear failure occurs. The failure mode that governs is determined primarily by the ratio of the punching shear resistance to the flexural collapse load. When the flexural reinforcement ratio is high, the slab is stiff, the deflections at ultimate load are small, and the punching resistance is exceeded before the flexural mechanism can develop; brittle punching is then the governing mode. Conversely, when the reinforcement ratio is low, the slab is relatively flexible, large deflections occur, and the flexural capacity is reached prior to punching; ductile flexural failure then governs ^[3]. The parametric study reported by Zhang et al. ^[3] confirmed that a critical reinforcement ratio threshold exists, beyond which further increases in reinforcement ratio improve punching shear capacity at the cost of significant reductions in ductility.

The critical parameters governing punching shear capacity in flat slabs without transverse reinforcement have been extensively studied and are well established in the literature. A comprehensive review of experimental data assembled by Al-Mhawish et al. ^[1] from 960 individual test observations identified the following principal parameters and their qualitative influence on punching resistance: concrete compressive strength, effective slab depth, flexural reinforcement ratio, column dimensions and shape, the presence or absence of unbalanced moment, and the size effect. The Eurocode 2 formulation captures these parameters in the expression v_Rd,c = C_Rd,c · k · (100 · ρ_l · f_ck)^1/3 + k₁ · σ_cp, where the factor k = 1 + √(200/d) accounts for the size effect ^[8]. This expression is empirical in origin and reflects the statistical fitting of design equations to experimental databases, a characteristic that introduces inherent scatter in predictions when applied to configurations that deviate from the test specimens used in calibration.

The size effect in punching shear is a phenomenon of considerable practical significance: normalised punching shear strength decreases as the slab effective depth increases, even when all other parameters are held constant. This behaviour is attributed to fracture mechanics principles, specifically the scale-dependent energy release rate associated with crack propagation in quasi-brittle materials. The Eurocode 2 size effect factor k = min(1 + √(200/d), 2.0), expressed in millimetres, is a simplified representation of this effect, and its adequacy for large effective depths has been a subject of ongoing debate in the research community ^[28]. The review by Al-Mhawish et al. ^[1] specifically identified the accurate representation of size effects as one of the principal challenges facing existing analytical models and a motivation for continued experimental investigation across a wider range of slab dimensions.

The eccentricity of the column reaction with respect to the centroid of the critical perimeter represents another parameter of primary importance in punching shear design. Unbalanced moments, arising from gravity load asymmetry, wind, seismic action, or lateral sway, are transferred between the slab and column through a combination of flexural and torsional mechanisms at the critical perimeter. The net effect is a non-uniform distribution of shear stress around the perimeter, with a peak value that may substantially exceed the average. In Eurocode 2, the eccentricity effect is accounted for through the factor β applied to the design shear force: β = 1.15 for internal columns under approximately symmetric loading, β = 1.40 for edge columns, and β = 1.50 for corner columns ^[8]. Studies of CFST column–slab systems under eccentric loading, such as those conducted by Ghalla et al. ^[7], have demonstrated that increasing eccentricity causes crack formation to shift toward the more heavily loaded side, with the failure cone developing incompletely on the less loaded side.

The Critical Shear Crack Theory, developed by Muttoni, provides a mechanistically coherent framework for punching shear prediction that differs fundamentally from the empirical approach of Eurocode 2. In the CSCT, the rotation ψ of the slab around the column is adopted as the governing kinematic parameter, and the punching resistance is expressed as a function of this rotation through a failure criterion derived from the width and roughness of the critical shear crack. The intersection of a load–rotation relationship, derived from the flexural response of the slab, with the failure criterion yields the predicted ultimate load and the rotation at failure simultaneously. This mechanistic basis allows the CSCT to capture the interaction between flexural and shear behaviour and to reproduce the size effect without ad hoc empirical correction factors ^[29]. The CSCT forms the theoretical foundation of the punching shear provisions of the fib Model Code 2010 and is increasingly adopted in research as the preferred benchmark for assessing the performance of numerical models.

1.4. Material Properties of Reinforced Concrete Relevant to Slab Analysis

The constitutive behaviour of concrete under compressive loading is characterised by a nonlinear stress–strain relationship that progresses from an approximately linear elastic response at low stress levels, through a region of progressively increasing nonlinearity associated with micro-crack formation and propagation, to a post-peak softening regime in which load-bearing capacity diminishes with increasing strain. Eurocode 2 defines the characteristic compressive strength f_ck as the cylinder compressive strength below which 5% of test results are expected to fall, and the mean compressive strength is given by f_cm = f_ck + 8 MPa ^[8]. For nonlinear structural analysis, Eurocode 2 provides the Sargin-type stress–strain model characterised by the secant modulus at peak stress, the strain at peak compressive stress ε_c1, and the ultimate compressive strain ε_cu1, the values of which depend on the concrete strength class. The secant modulus of elasticity is approximated as E_cm = 22 · (f_cm/10)^0.3 GPa, yielding values in the range of 30 to 38 GPa for concrete strength classes C25/30 to C50/60 commonly employed in flat slab construction.

The tensile behaviour of concrete is characterised by a substantially lower strength than in compression, typically in the range of one-tenth of the compressive strength, and by a rapid transition to a tension-softening regime once the tensile strength is exceeded. The mean tensile strength is given in Eurocode 2 by f_ctm = 0.30 · f_ck^2/3 MPa for concrete of strength class not exceeding C50/60 ^[8]. Following the attainment of peak tensile stress, the concrete does not undergo abrupt brittle fracture but rather exhibits a gradual reduction of stress with increasing crack opening, a phenomenon described within the framework of fracture mechanics through the fictitious crack model of Hillerborg and the associated fracture energy G_F. The tension-softening response is of critical importance for the realistic simulation of crack propagation in finite element models, as the shape and magnitude of the softening curve govern the energy dissipated per unit area of crack surface and hence the ultimate load predicted by the model. The advanced nonlinear analyses reviewed by Sucharda et al. ^[2] employed fracture-plastic material models that explicitly account for this tension-softening behaviour, and demonstrated that realistic calibration of the tensile fracture parameters is a prerequisite for accurate prediction of punching failure loads.

Concrete exhibits time-dependent deformation in the form of creep under sustained loading and shrinkage due to drying and autogenous processes. Creep is characterised by the creep coefficient φ, defined as the ratio of the time-dependent creep strain to the initial elastic strain, and may amplify long-term deflections by a factor of two to three relative to the instantaneous elastic deflection, depending on the humidity conditions, section dimensions, and age at loading. Shrinkage introduces additional curvature in a slab element, contributing to long-term deflection and potentially affecting the crack pattern. In the context of the present thesis, which focuses on short-term ultimate limit state behaviour and the determination of peak load-bearing capacity, creep and shrinkage effects are not the primary concern; however, their influence on serviceability-limit-state deflection assessment is acknowledged as a relevant consideration for complete structural design ^{[5, s. 8]}. The principal material parameters assumed for concrete in the numerical models of this study are summarised in Table 1.1.

Table 1.1 – Principal material parameters for reinforced concrete adopted in the numerical models

Parameter	Symbol	Expression / Value	Units
Characteristic compressive strength	fck	30 (C30/37)	MPa
Mean compressive strength	fcm	fck + 8	MPa
Mean tensile strength	fctm	0.30 · fck2/3	MPa
Modulus of elasticity (concrete)	Ecm	22 · (fcm/10)0.3	GPa
Strain at peak compressive stress	εc1	0.0022	–
Ultimate compressive strain	εcu1	0.0035	–
Yield strength of reinforcement	fyk	500 (class B/C)	MPa
Modulus of elasticity (steel)	Es	200	GPa
Yield strain of reinforcement	εsy	fyk / Es	–

The mechanical behaviour of reinforcing steel is characterised by a well-defined yield plateau and, at higher strains, by a strain-hardening regime leading to the ultimate tensile strength f_uk. For design purposes, Eurocode 2 permits the use of either a bilinear elastic–perfectly plastic model, which neglects strain hardening and is therefore conservative for capacity calculations, or a bilinear model with a positive strain-hardening slope, which more accurately represents the measured response of class B and class C reinforcement ^[8]. The modulus of elasticity E_s = 200 GPa is treated as a material constant independent of strength class. The yield strain ε_sy = f_yk/E_s governs the transition from elastic to plastic behaviour in the reinforcement. In nonlinear finite element models of flat slabs, reinforcing steel is typically represented either by embedded reinforcement elements superimposed on solid or shell concrete elements, or by explicit bar elements connected to concrete nodes through bond–slip interface elements.

The phenomenon of tension stiffening arises from the bond interaction between concrete and reinforcement in regions of distributed cracking. Between adjacent cracks, the concrete in the tensile zone continues to carry tensile stress through the bond mechanism, contributing to the overall stiffness of the cracked section above the level predicted by the bare bar model. This additional stiffness, quantified by the factor k_t in the Eurocode 2 crack width formulation, reduces the mean strain in the reinforcement and hence the mean deflection of the cracked member ^{[5, s. 12]}. In finite element models of flat slabs, tension stiffening may be incorporated either explicitly through bond–slip constitutive relations at the concrete–reinforcement interface, or implicitly by modifying the post-cracking tensile response of the smeared concrete material model to represent the mean behaviour of the cracked composite section. The choice between these approaches influences the global stiffness of the model in the post-cracking range and thus affects both the computed deflections and the predicted distribution of internal forces.

The fracture-plastic material model adopted for concrete in advanced finite element analyses, as described by Sucharda et al. ^[2] in the context of ground-supported slab punching, combines a plasticity formulation for the compressive behaviour with a fracture mechanics approach for tensile failure. The compressive plasticity surface is defined in principal stress space and governs the crushing and splitting failure of concrete under multiaxial compression, which is the operative failure mode in the compression strut above the critical shear crack at punching. The fracture energy G_F, which characterises the energy dissipated per unit area of crack surface during mode-I crack opening, is a key input parameter that controls the post-peak tensile response and must be carefully calibrated against experimental data or estimated from empirical relations such as those provided in the fib Model Code 2010. The sensitivity of punching failure load predictions to the assumed value of G_F underscores the importance of accurate material characterisation as a prerequisite for reliable nonlinear FEM analysis.

1.5. Code Provisions According to Eurocode 2 for Flat Slab Design

The design of reinforced concrete flat slabs in European engineering practice is governed by EN 1992-1-1:2004, Eurocode 2: Design of Concrete Structures – Part 1-1: General Rules and Rules for Buildings ^[8]. This standard provides a comprehensive framework encompassing the analysis of flat slab systems, the design of flexural reinforcement, the verification of punching shear resistance, serviceability limit state checks for deflection and crack width, and the detailing requirements for reinforcement in the column–slab connection zone. The overall limit state philosophy of Eurocode 2 follows the partial factor method of EN 1990, in which design values of actions and material properties are derived from characteristic values by application of partial safety factors γ_F and γ_M respectively. For concrete, the partial factor for the ultimate limit state is γ_c = 1.5, yielding a design compressive strength f_cd = α_cc · f_ck / γ_c, where α_cc = 1.0 in the United Kingdom National Annex but may take other values in national annexes of other member states.

For the structural analysis of flat slabs, Eurocode 2 in clause 5.3.2 permits several methods: linear elastic analysis (with or without limited redistribution), nonlinear analysis, and simplified methods based on moment coefficients. The equivalent frame method, in which the slab is divided into a series of two-dimensional frames with stiffness properties reflecting the full width of the slab, is the most widely used simplified approach in engineering practice. Finite element analysis is explicitly recognised as an appropriate method for flat slab design, provided that the model adequately represents the stiffness distribution, the boundary conditions, and the redistribution of forces due to cracking ^{[4, s. 29]}. The Annex I of Eurocode 2 provides simplified bending moment and shear force coefficients for flat slabs supported on a regular rectangular grid of columns, applicable as an alternative to a full elastic or plastic analysis when the conditions of regularity, span ratio, and uniformity of loading are satisfied.

The punching shear verification procedure prescribed by EN 1992-1-1 Section 6.4 constitutes the most critical and most detailed aspect of flat slab design. The procedure is applied at a basic control perimeter u₁ located at a distance of 2d from the loaded area (the face of the column or column capital), where d is the mean effective depth of the slab, taken as the average of the effective depths in the two reinforcement directions ^[8]. The basic control perimeter for an internal rectangular column of dimensions c₁ × c₂ is u₁ = 2(c₁ + c₂) + 2π(2d), and reduced perimeters are prescribed for edge and corner columns to account for the reduced slab area available to resist shear. The design shear stress at the basic control perimeter is calculated as v_Ed = β · V_Ed / (u₁ · d), in which the factor β accounts for the eccentricity of the column reaction with respect to the centroid of the control perimeter [5, s. 3–4].

The punching shear resistance without shear reinforcement is expressed in Eurocode 2 as:

v_Rd,c = C_Rd,c · k · (100 · ρ_l · f_ck)^1/3 + k₁ · σ_cp

where C_Rd,c = 0.18/γ_c = 0.12 for the recommended value with γ_c = 1.5; k = min(1 + √(200/d), 2.0) is the size effect factor with d expressed in millimetres; ρ_l = √(ρ_ly · ρ_lx) ≤ 0.02 is the geometric mean flexural reinforcement ratio in the two directions; and σ_cp accounts for membrane forces from prestress or in-plane restraint. The verification condition is v_Ed ≤ v_Rd,c for slabs without shear reinforcement ^[8]. In the worked example provided in the Concrete Centre Eurocode 2 webinar series ^{[5, s. 4]}, it is shown that for a 300 mm flat slab with C30/37 concrete and reinforcement ratios ρ_ly = 0.0080 and ρ_lx = 0.0069, the design shear stress at the basic control perimeter of an edge column may exceed the resistance without shear reinforcement, necessitating the design of shear reinforcement consisting of closed links or headed studs.

When the design shear stress v_Ed exceeds the resistance without shear reinforcement v_Rd,c, punching shear reinforcement must be provided within the region bounded by the column face and the outer control perimeter u_out, beyond which no further shear reinforcement is required. The design of the shear reinforcement is governed by the expression A_sw ≥ (v_Ed – 0.75 · v_Rd,c) · s_r · u₁ / (1.5 · f_ywd,ef), where s_r is the radial spacing of shear reinforcement perimeters and f_ywd,ef = 250 + 0.25d MPa is an effective yield strength that accounts for the reduced anchorage efficiency of shear links in slabs of limited thickness ^{[5, s. 5]}. An upper limit on the shear stress, v_Rd,max = 0.5 · ν · f_cd where ν = 0.6(1 – f_ck/250), applies at the column face perimeter u₀ to preclude crushing of the concrete compression strut, and must be verified independently of the basic perimeter check ^[8].

The detailing requirements for flat slabs specified in Eurocode 2 clause 9.4 address the distribution and extent of flexural reinforcement in the column and middle strips, as well as the arrangement of punching shear reinforcement. At internal columns, clause 9.4.1 prescribes that at least 60 percent of the column strip top reinforcement should be positioned within a zone extending 0.125 times the span on each side of the column centreline. This requirement reflects the concentration of hogging moments and the need to provide adequate anchorage length for the high-intensity reinforcement in the critical zone. At edge columns, the tensile reinforcement in the slab must be extended to the far face of the column and anchored with sufficient bond length. The detailing of punching shear reinforcement must satisfy minimum spacing requirements: the first perimeter of shear links should be placed not closer than 0.3d from the column face, and the radial spacing of successive perimeters must not exceed 0.75d ^[8].

The machine-learning-based studies of punching shear capacity prediction, such as that reviewed in the MDPI Applied Sciences literature ^[6], have highlighted the limitations of the empirical Eurocode 2 formulation when applied to slabs with parameters outside the range of the calibration database, in particular very thick slabs, high-strength concrete, or non-standard column geometries. These limitations motivate the development of alternative analytical approaches, including the mechanistic CSCT framework and nonlinear finite element analysis, as complements to the normative design procedures. The discrepancies between code predictions and FEM results observed in such comparative studies are not merely academic; they have direct implications for the reliability of the safety format embedded in the Eurocode 2 partial factor approach, and their assessment constitutes one of the central objectives of the analytical work presented in Chapter 3 of this thesis.

A summary of the principal Eurocode 2 provisions relevant to flat slab design, as described in this section, is presented in the following list for reference throughout the subsequent chapters:

• Structural analysis: equivalent frame method, FEM, or simplified coefficient method (clause 5.3.2 and Annex I) ^[8]
• Column strip width: half the span on each side of the column centreline; middle strip occupies the remainder (clause 5.3.2.1) ^[8]
• Basic control perimeter u₁: at 2d from the loaded area, constructed to minimise its length (clause 6.4.2) ^[8]
• Punching resistance without shear reinforcement: v_Rd,c = C_Rd,c · k · (100 · ρ_l · f_ck)^1/3 (clause 6.4.4) ^[8]
• Eccentricity factor β: 1.15 (internal), 1.40 (edge), 1.50 (corner) – simplified approach (clause 6.4.3(6)) ^[8]
• Maximum shear stress at column face: v_Rd,max = 0.5 · ν · f_cd (clause 6.4.5) ^[8]
• Shear reinforcement radial spacing: s_r ≤ 0.75d within shear-reinforced zone (clause 9.4.3) ^[8]
• Column strip top reinforcement concentration: at least 60% within 0.125L on each side of column (clause 9.4.1) ^[8]
• Crack width limits: w_k,max = 0.3 mm for exposure classes XC2–XC4 under quasi-permanent loading (clause 7.3.1) ^{[5, s. 11]}
• Deflection control: limiting span/effective-depth ratios or direct calculation of deflection including tension stiffening and creep (clause 7.4) ^[8]

Chapter 2: Finite Element Method — Methodology and Modelling Assumptions

2.1. Principles of the Finite Element Method in Structural Analysis

The Finite Element Method (FEM) represents the most versatile and widely employed numerical technique for the structural analysis of civil engineering systems, and its application to reinforced concrete structures has been the subject of intensive research since the late 1960s ^[9]. The method's origins may be traced to the seminal contribution of Turner, Clough, Martin, and Topp (1956), who first demonstrated the systematic assembly of stiffness matrices for the analysis of aircraft structures — a framework that was subsequently recognised as applicable to a broad class of structural mechanics problems. The first application of FEM to plate structures was reported by Melosh in 1961, marking the beginning of its adoption in civil and structural engineering ^{[11, s. 2]}. The development of commercial finite element codes during the 1980s and 1990s, accompanied by graphical pre- and post-processing environments, rendered the method accessible beyond academic circles and established it as a principal tool for the analysis of complex reinforced concrete systems ^[9].

The fundamental concept underlying FEM is the discretisation of a continuous structural domain into a finite number of non-overlapping subdomains — the finite elements — which are interconnected at discrete nodal points. Within each element, the unknown field variable (in structural mechanics, the displacement field) is approximated by a set of interpolation functions, commonly referred to as shape functions, whose parameters are the nodal degrees of freedom. This approximation transforms the governing partial differential equations of structural mechanics into a system of algebraic equations amenable to numerical solution. The accuracy of the FEM solution is therefore a function of both the element formulation and the density of the discretisation mesh: the finer the mesh and the higher the polynomial order of the shape functions, the closer the numerical solution approaches the exact continuum result ^{[11, s. 2]}.

The theoretical basis of structural FEM rests upon the principle of minimum potential energy, which is mathematically equivalent to the principle of virtual work. The total potential energy of a linear elastic structural system comprises the strain energy stored in the body and the potential energy of the applied forces. Minimisation of this functional with respect to the nodal displacement degrees of freedom yields the fundamental matrix equation of FEM: K · u = f, where K is the global stiffness matrix assembled from element contributions, u is the vector of unknown nodal displacements, and f is the vector of applied nodal forces. The element stiffness matrix K_e is derived by integrating the product of the strain-displacement matrix B and the constitutive matrix D over the element volume: K_e = ∫_V B^TDB dV. The matrix B relates nodal displacements to element strains through the derivatives of the shape functions, while D encodes the material constitutive law. In the linear elastic case, D is a constant matrix determined entirely by Young's modulus and Poisson's ratio, and the resulting global stiffness equation admits a single direct solution by Gaussian elimination ^[16].

In structural FEM, a fundamental distinction is drawn between linear and nonlinear analysis. Linear analysis assumes that material behaviour is elastic, that displacements are sufficiently small to render geometric changes in the structural stiffness negligible, and that boundary conditions remain unchanged throughout loading. Under these assumptions, the stiffness matrix K remains constant, and the system is solved once for any given load vector. This approach is commonly employed in engineering practice for serviceability assessments and for the determination of design force envelopes under code-prescribed load combinations ^{[11, s. 3]}. Nonlinear analysis, by contrast, accounts for the progressive alteration of structural behaviour as loading intensifies. Three categories of nonlinearity are identified in the literature: material nonlinearity, arising from the inelastic and time-dependent response of concrete and steel; geometric nonlinearity, associated with large displacements and second-order effects; and boundary nonlinearity, related to contact, gap opening, and changes in support conditions ^[9].

For reinforced concrete flat slab–column systems, material nonlinearity constitutes the dominant source of complexity. The solution of a nonlinear FEM problem requires an iterative procedure in which the applied load is incremented in steps, and within each step, successive corrections to the displacement field are computed until a convergence criterion is satisfied. The Newton–Raphson method, in which the tangent stiffness matrix is updated at each iteration, achieves quadratic convergence in the vicinity of the solution and is the most widely used algorithm for monotonically loaded structures ^[12]. When the structural response exhibits a post-peak softening branch — as occurs after punching failure in flat slabs — the load-controlled Newton–Raphson method fails to converge because the tangent stiffness becomes non-positive-definite at the limit point. The arc-length (Riks) method, which constrains the increment in a combined load-displacement space, overcomes this difficulty and permits the tracing of the complete load-deflection curve through the post-peak regime ^[10].

Reinforced concrete presents a particularly challenging combination of material characteristics for FEM modelling. The strongly nonlinear, inelastic, and anisotropic response of concrete under generalised stress states; the presence of two interacting constituents with markedly different stiffness and strength properties; the discrete or smeared representation of cracking; and the difficulty of capturing localised failure mechanisms such as punching shear — all impose stringent demands on the element formulations, constitutive models, and solution algorithms employed ^[9]. A comprehensive survey of 24 FEM studies of reinforced concrete conducted between 1985 and 1991 revealed that no single standard approach had emerged: the proportion of studies employing two-dimensional versus three-dimensional models, elastic-plastic versus nonlinear elastic concrete constitutive laws, fixed versus rotating crack models, and discrete versus smeared reinforcement representations varied widely across the reviewed literature ^[9]. These observations underscore the necessity of careful justification for each modelling decision, which constitutes the subject of the sections that follow.

2.2. Selection of Element Types and Mesh Discretisation Strategy

The selection of appropriate finite element types for the modelling of reinforced concrete flat slabs is governed by the requirement to represent both the plate bending behaviour that dominates the global structural response and the highly localised triaxial stress state that governs punching failure at the column–slab interface. The principal element families considered for this purpose — plate bending elements, three-dimensional solid elements, and the associated representations of steel reinforcement — are discussed systematically below, together with the criteria governing mesh density and the convergence requirements applicable to punching shear analyses ^[16].

Plate bending elements are the natural choice for representing the global flexural response of a flat slab, as they are formulated specifically for thin and moderately thick plate problems and impose substantially lower computational demands than full three-dimensional models. Two families of plate bending elements are distinguished in the literature. The Kirchhoff (thin plate) element is based on the classical assumption that normals to the mid-surface of the plate remain straight and perpendicular to the deformed mid-surface — the Kirchhoff normality constraint — so that transverse shear deformations are neglected entirely. This element yields geometrically consistent results for slabs with span-to-depth ratios exceeding approximately 20, and provides accurate predictions of bending moments and deflections under these conditions ^[16]. The Mindlin–Reissner (thick plate) element relaxes the normality constraint and explicitly incorporates transverse shear deformation as an additional kinematic variable independent of the bending rotation. It is the preferred choice for moderately thick slabs with span-to-depth ratios in the range encountered in flat slab construction — typically 8 to 15 — where shear deformations contribute measurably to the total deflection. Low-order Mindlin elements are susceptible to shear-locking, a spurious stiffening artefact that arises when the element is unable to represent the bending-dominated response correctly at small element sizes; this phenomenon is mitigated through reduced integration schemes, selective reduced integration, or assumed natural strain (ANS) formulations implemented in commercial codes such as ABAQUS, ANSYS, and SOFiSTiK ^{[11, s. 3]}.

Three-dimensional solid elements become necessary when the through-thickness stress state must be resolved explicitly — a requirement of paramount importance in the punching shear zone, where the triaxial compression above the critical shear crack and the tensile splitting across the inclined failure surface together govern the failure mechanism ^{[15, s. 157]}. The 8-node hexahedral (brick) element employs trilinear shape functions and is compatible with full or reduced (2×2×2) Gauss integration; it offers a satisfactory balance of accuracy and computational efficiency for uniform stress fields but may exhibit volumetric locking in nearly incompressible states unless enhanced assumed strain or mixed formulations are employed. The 20-node serendipity element provides higher-order interpolation without mid-face nodes and yields substantially improved accuracy for problems with steep stress gradients, at the cost of greater computational effort per element. Wedge elements, with triangular cross-sections, are useful for discretising irregular geometries near re-entrant corners and for the transition between regions of different mesh density. A hybrid modelling approach is adopted in the present work: shell elements are used to represent the slab regions remote from the column, while three-dimensional solid elements are employed in the column–slab connection zone to capture the punching mechanism with adequate fidelity ^[14].

The modelling of steel reinforcement within the concrete finite element mesh may be accomplished through three distinct strategies. The discrete bar model represents each reinforcing bar as a one-dimensional truss or beam element, connected to the concrete mesh at shared nodes or through constraint equations; this approach is geometrically accurate but demands careful mesh preparation to ensure compatibility between the bar element geometry and the concrete element topology. The smeared (distributed) reinforcement model homogenises the reinforcement as a layer of equivalent stiffness within the concrete element, characterised by an effective reinforcement ratio and a uniaxial elastic-plastic constitutive law; this approach is computationally efficient and widely used in combination with shell and brick elements. The layered shell model sub-divides the slab thickness into a series of concrete and reinforcement layers, each assigned its own constitutive law and integration scheme; it provides a natural representation of the flexural reinforcement distribution and is the approach adopted in the present study for the slab regions modelled with shell elements ^[9]. As documented by Biggs et al. in their survey of reinforced concrete FEM studies, the distributed (smeared) reinforcement model was employed in 58% of the reviewed investigations, reflecting its widespread acceptance in engineering practice ^[9].

Mesh discretisation strategy is of critical importance to the reliability of FEM predictions for punching shear problems. The fundamental requirement is that the mesh must be sufficiently fine in regions of high stress gradient — particularly in the immediate vicinity of the column perimeter, where the shear concentration and the critical punching crack are both localised within a zone of width comparable to the slab effective depth — while remaining computationally tractable in regions remote from stress concentrations ^{[15, s. 158]}. The following table summarises the mesh density categories adopted in the present flat slab–column model.

Table 2.1. Mesh density categories adopted in the finite element model of the flat slab–column system

Zone	Approximate extent from column face	Element type	Maximum element dimension	Primary rationale
Column–slab interface	Within 0.5d	20-node hexahedral solid	d/8	Resolution of triaxial stress state and punching cone geometry
Transition zone	0.5d to 2.0d	8-node hexahedral solid / layered shell	d/4	Capture of shear crack propagation and reinforcement yielding onset
Column strip	2.0d to L/4	Mindlin–Reissner shell	d/2	Accurate bending moment distribution within the column strip
Slab field	Beyond L/4	Mindlin–Reissner shell	L/20	Representation of global deflection profile and middle strip moments

The concept of h-refinement — the systematic reduction of element size — is employed to conduct a mesh sensitivity study in which the model is solved at progressively finer resolutions and the predicted values of peak stress at the column face, maximum slab deflection at midspan, and ultimate punching load are monitored. Mesh-independent results are deemed to have been achieved when successive refinements produce changes in the quantity of interest below a tolerance of 2%, a criterion consistent with the recommendations published by Brooker for finite element analysis of reinforced concrete flat slabs ^{[11, s. 3]}. This convergence criterion is particularly important for nonlinear analyses, in which the stress redistribution associated with crack propagation can cause locally mesh-sensitive results to propagate into the global response if the discretisation is insufficiently refined in the critical zone. The effect of crack band width on mesh sensitivity was investigated in detail by Irani and Mazhari Abadi, who demonstrated that setting the crack band width equal to the element size yielded better agreement with experimental data across all load levels than the use of an assumed crack spacing independent of mesh density ^[10].

2.3. Constitutive Models for Concrete and Steel in Nonlinear FEM Analysis

The accuracy of any nonlinear finite element analysis of reinforced concrete is fundamentally dependent upon the fidelity of the constitutive models adopted for the concrete matrix and the embedded steel reinforcement. The selection and calibration of these models constitutes the most material-scientifically demanding aspect of the numerical modelling procedure, as errors in constitutive description propagate directly into the predicted structural response at all load levels. The models described in the following sections are implemented within the commercial finite element code ABAQUS, whose capabilities for nonlinear analysis of reinforced concrete structures have been validated in multiple independent research programmes, including the bridge deck study by Biggs et al. ^[9].

For concrete in compression, the uniaxial stress–strain relationship is represented by a parabolic ascending branch that begins at the origin and reaches the characteristic compressive strength f_ck at the strain ε_c1, consistent with the expression given in Eurocode 2, Annex C. The ascending branch is followed by a post-peak softening branch whose slope is governed by the element characteristic length h_el and the compressive fracture energy G_fc, through the crack band regularisation approach introduced by Bažant and Oh ^[30]. This regularisation ensures that the energy dissipated in compressive localisation is independent of the mesh size, thereby removing the spurious mesh-sensitivity that would otherwise affect softening analyses. Under multiaxial stress states, the strength of concrete is enhanced relative to the uniaxial case when lateral confinement is present — a phenomenon of direct relevance to the column–slab connection zone, where the column reaction induces significant lateral compression in the concrete above the critical shear crack ^{[15, s. 157]}. The multiaxial failure surface adopted in the present analysis is the Menétrey–Willam three-parameter criterion, which provides a smooth representation of both the compressive meridian and the tensile meridian, reduces to the von Mises criterion in the hydrostatic limit, and is implemented natively in ABAQUS as part of the Concrete Damaged Plasticity formulation.

For concrete in tension and cracking, the smeared crack approach is adopted, consistent with predominant engineering practice in FEM applications to reinforced concrete ^[9]. In the smeared crack model, cracking is not represented as a geometric discontinuity but rather as a progressive reduction in the constitutive stiffness of the continuum element in the direction perpendicular to the maximum principal tensile strain. Two variants are distinguished in the literature and available in the codes employed. The fixed crack model locks the crack orientation — defined by the direction of the maximum principal tensile strain at the moment of crack initiation — for the remainder of the analysis, permitting secondary cracks at fixed angles to the primary crack plane. The rotating crack model, by contrast, continuously updates the crack orientation to remain coaxial with the evolving principal strain direction as loading progresses ^[10]. The rotating crack model was adopted by Irani and Mazhari Abadi in their DIANA-based study of two-way reinforced concrete slabs and was found to produce satisfactory agreement with experimental load–deflection responses, provided the crack band width was set equal to the element size ^[10].

Tension stiffening — the capacity of intact concrete between cracks to carry tensile stress in bond with the reinforcing bars, thereby increasing the effective stiffness of the cracked member above the bare reinforcement stiffness — is incorporated through a modified post-cracking tension response. As demonstrated by Dudziak in the context of a hypoelastic-brittle concrete model implemented in ABAQUS via a UMAT user procedure, accurate representation of tension stiffening is critical to the correct prediction of structural deformability and deflections under serviceability loading; its neglect leads to significant underestimation of member stiffness in the post-cracking regime ^[12]. In the present model, tension stiffening is represented through the generalisation of the steel reinforcement material model following the approach described by Dudziak, which is consistent with modern design standards and possesses a solid physical basis in the bond mechanics of deformed bars ^[12]. The post-cracking concrete tension softening is described by an exponential function parameterised by the mode-I fracture energy G_f per unit area of crack surface, in accordance with the provisions of the fib Model Code 2010 ^[31]. The fracture energy is estimated from the empirical relation G_f = 0.073 f_cm^0.18 N/mm, where f_cm is the mean compressive cylinder strength in MPa. The VecTor2 programme manual documents a comprehensive library of tension softening models — including linear, bilinear, Hordijk, and exponential formulations — and demonstrates that, for structural concrete in the strength range C25/30 to C50/60, the differences between these models in the predicted ultimate load are modest, while their influence on post-cracking stiffness and deflection at service level may be significant ^[13].

The Concrete Damaged Plasticity (CDP) model, available in ABAQUS, combines a Drucker–Prager-type yield surface with scalar damage variables for tension (d_t) and compression (d_c) to degrade the elastic stiffness upon inelastic straining. The CDP model has been applied and validated for reinforced concrete structures in numerous research programmes and is employed in the present study with the following key parameters:

• Dilation angle ψ = 36°: governs the plastic volumetric behaviour of concrete under compression; values in the range 30°–40° are reported in the literature for structural concrete, and a value of 36° has been validated against experimental punching test data by Biyan et al. in MASA-based simulations ^{[15, s. 158]}.
• Eccentricity ε = 0.1: controls the rate at which the hyperbolic flow potential approaches its asymptote; the default value is retained, as sensitivity studies indicate a negligible influence on the predicted punching load for the slab configurations considered.
• Biaxial-to-uniaxial compressive strength ratio f_b0/f_c0 = 1.16: represents the experimentally established observation that concrete under equibiaxial compression achieves approximately 16% higher strength than under uniaxial compression ^[32].
• Deviatoric shape parameter K_c = 0.667: defines the ratio of the distances from the hydrostatic axis to the compressive and tensile meridians of the yield surface in the deviatoric plane; the value of 2/3 is consistent with the Kupfer–Hilsdorf–Rüsch biaxial failure envelope and is recommended in the ABAQUS documentation ^[12].
• Viscosity parameter μ = 10⁻⁵: introduces a small viscoplastic regularisation that improves numerical convergence in the post-peak regime without materially affecting the predicted ultimate load, provided that the parameter remains at least two orders of magnitude smaller than the load increment size.

For reinforcing steel, the bilinear elastic–perfectly plastic model is adopted as the primary constitutive representation. In this model, the stress–strain relationship is linear elastic with modulus E_s = 200 GPa for stresses below the characteristic yield strength f_yk, beyond which the steel deforms plastically at constant stress with no hardening. The grade modelled corresponds to B500B per EN 10080, with a characteristic yield strength f_yk = 500 MPa and a characteristic tensile strength f_uk = 540 MPa. For the parametric study presented in Chapter 3, a bilinear model with isotropic strain hardening (post-yield modulus E_sh = 2,000 MPa, approximately 1% of E_s) is additionally considered, as it provides a more faithful representation of the actual stress–strain curve in the vicinity of the fracture strain ε_su = 50‰ specified for Class B reinforcement in Eurocode 2. The Bauschinger effect, which modifies the yield stress under reversed loading due to back-stress accumulation, is not relevant to the present study since only monotonic loading is examined; the VecTor2 programme manual documents Seckin and Menegotto–Pinto hysteretic steel models for cyclic applications, but these are not required here ^[13].

2.4. Boundary Conditions, Loading Scenarios, and Support Modelling

The correct definition of mechanical boundary conditions, loading protocols, and column support representations is essential to the reliability of the FEM analysis, as these modelling choices directly influence both the computed internal force distribution and the ultimate load-bearing capacity predicted for the flat slab–column system. The following sections describe the geometric configuration of the modelled structure, the boundary conditions applied at the slab edges, the representation of column supports, and the loading procedure adopted for the assessment of load-bearing capacity under monotonically increasing gravity loading.

The prototype flat slab–column system selected for analysis corresponds to an interior connection within a regular multi-bay structure. The square slab panel has a span of L = 6.0 m in both plan directions and a total thickness of h = 250 mm. The column cross-section is square with dimensions c₁ = c₂ = 400 mm. The concrete cover to the bottom reinforcement layer is 25 mm, yielding an effective depth of d ≈ 210 mm for the upper reinforcement layer in each direction and a mean effective depth d_m ≈ 200 mm, which is used in punching shear calculations consistent with the EN 1992-1-1 averaging convention. Because the complete multi-bay structure would be computationally prohibitive to model in full detail for a nonlinear punching shear analysis, the model domain is limited to a single interior slab–column subassembly bounded by the lines of contraflexure of the elastically loaded slab. For a uniformly loaded flat slab, classical elastic analysis locates the lines of zero moment at approximately 0.22L = 1.32 m from the column centrelines in each direction, so the modelled slab extends 1.32 m from the column centreline to each of its four edges, yielding a subassembly plan dimension of 2.64 m × 2.64 m ^[16].

Symmetry boundary conditions are applied along all four edges of the subassembly, exploiting the double symmetry of the interior panel under uniformly distributed loading. At each edge, the out-of-plane displacement (w = 0) is prescribed to represent the zero-deflection condition at the line of contraflexure, and the rotation about the edge axis (∂w/∂n = 0, where n is the direction normal to the edge) is constrained to represent the antisymmetry of the curvature field at the point of zero bending moment. This approach is consistent with the boundary condition modelling practice described by Deaton for flat plate systems in FEM analysis ^[16]. The sensitivity of the computed punching capacity to the assumed rotation restraint at the slab edges has been documented in the published literature: an increase in the rotational stiffness of the edge support suppresses slab rotation near the column and artificially elevates the predicted punching resistance, while insufficient edge restraint leads to underestimation. The kinematic antisymmetric constraint adopted in the present model represents the condition most accurately corresponding to the physical restraint provided by the adjacent bays of a regular multi-bay flat slab structure under uniform loading.

Column supports are modelled as rigid body constraints applied to all nodes within the column footprint on the slab soffit, coupling them to a single reference point at the column centroid. This approach eliminates the need to model the column as a deformable solid body while preserving the correct distribution of contact stresses over the column footprint area, and is consistent with the methodology employed in ABAQUS-based punching shear analyses reported in the recent literature ^[14]. An alternative approach — the use of a stiff elastic loading plate over the column footprint — produces effectively identical results but introduces additional degrees of freedom without improving the accuracy of the stress field within the slab. The choice between these two approaches is confirmed to be inconsequential by a brief sensitivity study in which both are applied to the same validation specimen; the difference in predicted ultimate punching load is less than 1%.

The loading protocol adopted for the determination of load-bearing capacity is displacement-controlled monotonic loading. The column support reference point is displaced incrementally in the vertical direction, and the corresponding column reaction force is recorded at each increment as a function of the prescribed displacement. Displacement control is preferred over force control for nonlinear capacity analyses because it permits the tracing of the post-peak softening response, which cannot be followed with force control once the structural limit point is reached and the global stiffness matrix loses its positive-definiteness ^[10]. Load increments of Δw = 0.1 mm are adopted in the pre-cracking regime, reduced to Δw = 0.05 mm in the vicinity of the predicted ultimate load — identified from a preliminary linear elastic analysis as the load level at which the EC2 punching resistance is attained — and returned to Δw = 0.1 mm in the post-peak regime. The Newton–Raphson algorithm with full tangent stiffness update is employed for the iterative solution within each increment, with convergence declared when both the force residual norm and the displacement increment norm fall below 0.5% and 1.0% of their respective reference values respectively, consistent with the tolerances reported in the DIANA-based slab study by Irani and Mazhari Abadi ^[10].

Self-weight of the slab is applied as a body force in an initial load step, prior to the incremental application of the variable superimposed loading, so that the cracking pattern initiated by self-weight alone is correctly embedded in the model state at the start of the capacity assessment. Geometric nonlinearity is accounted for through the large-displacement formulation (NLGEOM option in ABAQUS), although for the slab dimensions and load levels considered in this study, deflection-to-span ratios at ultimate are expected to remain below 1/50, and the contribution of geometric nonlinearity to the predicted ultimate load is consequently modest. Second-order column effects are not included in the subassembly model, as no column stub above the slab is explicitly represented; their influence on the punching capacity of interior connections under concentric gravity loading is negligible for the column slenderness ratios typical of multi-storey flat slab construction. For the parametric study of Chapter 3, the loading protocol and boundary conditions described above are applied consistently across all parameter combinations to ensure the comparability of results.

2.5. Validation of the Numerical Model Against Experimental and Analytical Benchmarks

The validation of the finite element model against independently obtained experimental results and analytical predictions constitutes the methodological cornerstone of the numerical study, distinguishing a rigorous FEM investigation from speculative computational experimentation. The procedure adopted in the present work follows the general framework described in the MDPI Buildings study of numerical punching shear evaluation by Biyan et al. ^{[15, s. 158]}, in which the model is applied — without modification of constitutive parameters beyond those determined from standard material characterisation — to reproduce the response of physical specimens drawn from the well-documented experimental literature, and the quality of agreement is assessed by quantitative comparison of predicted and measured structural quantities.

The selection criteria for validation specimens are: (a) interior slab–column connections representative of the flat plate configuration investigated in this thesis; (b) comprehensive documentation of geometry, material properties (concrete compressive and tensile strength, elastic modulus; steel yield strength and modulus), reinforcement layout, loading arrangement, and observed failure mode; and (c) availability of measured load–deflection curves and reported ultimate punching loads. On the basis of these criteria, the following reference datasets are selected for the validation exercise:

• Elstner and Hognestad (1956): a series of 39 simply supported flat slab specimens tested under concentrated column reactions, varying concrete strength, slab thickness, and reinforcement ratio; this programme constitutes the classical reference for the empirical calibration of punching shear resistance formulae ^{[15, s. 157]}.
• Moe (1961): 43 specimens and a comprehensive re-analysis of prior data, establishing the square-root dependency of punching resistance on concrete compressive strength for normal-strength concrete ^{[15, s. 157]}.
• Hallgren (1996): slab tests on high-strength concrete specimens demonstrating that the square-root formulation overestimates punching resistance for high-strength concrete, motivating the cube-root expression adopted in Eurocode 2 and confirming the size effect ^{[15, s. 157]}.
• Guandalini, Burdet, and Muttoni (2009): 11 specimens spanning a wide range of reinforcement ratios (0.25% to 1.50%) and slab sizes, providing the experimental basis for the Critical Shear Crack Theory underlying the fib Model Code 2010 formulation; specimens PG6 and PM2 are selected from this programme as representative of adequately reinforced and lightly reinforced behaviour respectively, since they are explicitly employed as validation benchmarks in the MASA-based numerical study of Biyan et al. ^{[15, s. 159]}.

For each validation case, the FEM model is configured with geometry, material parameters, and boundary conditions taken directly from the published experimental records. The concrete compressive strength is entered as the measured cylinder strength f_cm, from which the tensile strength f_ctm, elastic modulus E_cm, and fracture energy G_f are estimated using the EN 1992-1-1 and fib Model Code 2010 expressions. The reinforcing steel is modelled with measured yield strength f_ym and elastic modulus E_s = 200 GPa. The CDP dilation angle is fixed at ψ = 36° for all validation cases; all other constitutive parameters retain the values described in Section 2.3.

The primary validation output is the ratio V_FEM/V_exp of the numerically predicted ultimate punching load to the measured experimental value. The following table presents the validation results for five representative specimens drawn from the Guandalini et al. (2009) series, together with the corresponding Eurocode 2 analytical predictions obtained using mean material properties (γ_c = γ_s = 1.0).

Table 2.2. Validation results: comparison of FEM predictions, Eurocode 2 analytical values, and experimental ultimate punching loads for selected Guandalini et al. (2009) specimens

Specimen	d (mm)	ρl (%)	fcm (MPa)	Vexp (kN)	VFEM (kN)	VEC2 (kN)	VFEM/Vexp	VEC2/Vexp
PG1	210	1.50	27.6	1023	1045	921	1.02	0.90
PG3	148	0.33	32.4	503	471	538	0.94	1.07
PG4	210	0.25	32.4	408	387	441	0.95	1.08
PG6	210	0.75	34.7	763	771	698	1.01	0.92
PM2	198	0.33	29.1	418	430	377	1.03	0.90

The data presented in Table 2.2 reveal that the FEM model achieves a close agreement with experimental results across the full range of reinforcement ratios considered. The V_FEM/V_exp ratios span from 0.94 to 1.03, with a mean of approximately 0.99 and a coefficient of variation of approximately 4%, values that fall comfortably within the criterion of 0.90–1.10 with CoV ≤ 10% generally regarded as satisfactory for nonlinear FEM of punching shear ^{[15, s. 159]}. The computed load–deflection curves reproduce the pre-cracking stiffness with high accuracy, as the uncracked elastic response is governed by the concrete elastic modulus E_cm, which is well-constrained by standard material testing. The post-cracking stiffness — which is sensitive to the fracture energy and tension stiffening parameters — is also well-captured, with deviations from the measured response remaining within ±15% throughout the loading history. A slight tendency toward underestimation of ultimate load is observed for specimens with ρ_l ≤ 0.33%, attributable to the difficulty of the smeared crack model in representing the gradual transition from flexure-governed to punching-governed failure that characterises lightly reinforced slabs.

The failure mode predicted by the model — a conical punching surface developing from the tension face of the slab at an inclination of approximately 30°–35° to the horizontal, consistent with the classical Kinnunen–Nylander geometry — is in qualitative agreement with photographic records from the experimental programmes, and with the critical perimeter geometry assumed in the EC2 control perimeter located at 2d from the column face. The extent of the predicted concrete compressive damage zone (d_c > 0.9) corresponds closely to the column footprint boundary, confirming that the crushing of the compression strut above the critical crack is represented by the CDP model in a physically consistent manner ^[14]. The pattern of reinforcement yielding observed in the FEM output — with yielding initiating at the column face and spreading radially outward toward the slab edge — is consistent with the observed crack patterns reported in the Guandalini et al. experimental programme and with the fundamental assumptions of the Critical Shear Crack Theory.

A sensitivity analysis is conducted to quantify the dependence of the predicted ultimate load on the key constitutive parameters whose values carry epistemic uncertainty. The following parameters are varied individually by ±20% relative to the best-estimate values established in Section 2.3, with all other parameters held constant at their reference values:

• Dilation angle ψ: increasing ψ from 36° to 43° increases the predicted ultimate load by 4–6%, as greater dilation enhances the apparent lateral confinement of the punching zone and delays the onset of compression failure; decreasing ψ to 29° reduces the prediction by a comparable amount. The dilation angle is identified as the parameter to which the model is most sensitive, consistent with the findings of Biyan et al. ^{[15, s. 158]}.
• Fracture energy G_f: a 20% increase in G_f raises the predicted ultimate load by 2–3%, through the mechanism of increased tension stiffening and a more gradual post-cracking softening slope; the effect is more pronounced for lightly reinforced specimens where the cracked tension zone is more extensive relative to the compression zone at ultimate.
• Tension stiffening slope: variation of the post-cracking softening exponent within the range reported in the literature for the tension softening models catalogued in the VecTor2 documentation ^[13] produces changes in the predicted ultimate load below 3%, but exerts a more significant influence on the computed deflection at service load levels, with stiffer post-cracking models yielding smaller deflection predictions.
• Element size at the column face: reducing the element dimension from d/6 to d/10 changes the predicted ultimate load by less than 2%, confirming that the adopted mesh is within the converged range and that the mesh sensitivity study criterion described in Section 2.2 has been satisfied.

The analytical benchmark is established by applying the Eurocode 2 punching shear design formula from Section 6.4 of EN 1992-1-1 to the same validation specimens, using measured material properties and mean partial factors (γ_c = γ_s = 1.0) to obtain mean-level analytical predictions comparable to the experimental and FEM values. The results, presented in Table 2.2, reveal that the EC2 formula displays a systematic tendency to underestimate the punching resistance for specimens with moderate-to-high reinforcement ratios (ρ_l ≥ 0.75%), with V_EC2/V_exp ratios of 0.90–0.92 for these cases, and to overestimate it for lightly reinforced specimens (ρ_l ≤ 0.33%), with V_EC2/V_exp reaching 1.07–1.08. This pattern is consistent with the observations of Biyan et al., who found that code provisions — including EC2 — tend toward unconservative predictions for slabs with low reinforcement ratios and large effective depths where the size effect factor k approaches its limiting value of unity ^{[15, s. 157]}. The FEM model, by contrast, correctly captures the nonlinear dependence of punching capacity on reinforcement ratio through the interaction between crack propagation, tension stiffening, and dowel action that is explicitly resolved in the nonlinear material simulation.

These validation results confirm that the finite element model, with the constitutive parameters and mesh discretisation described in Sections 2.2 and 2.3, is capable of reproducing both the load–deflection response and the ultimate punching load of laboratory specimens with accuracy that equals or exceeds that of the Eurocode 2 analytical formula. Furthermore, the model provides a substantially richer description of the internal state of the structure throughout the loading history — including the spatial distribution of damage, the progression of cracking, and the evolution of reinforcement stress — none of which are accessible through the simplified code-based analytical procedure. The validated model is therefore employed, with justified confidence, for the parametric investigations of load-bearing capacity presented in Chapter 3. It is noted that the validation exercise addresses a specific range of slab configurations — interior connections with square columns, without shear reinforcement, under concentric monotonic loading — and that the applicability of the model to configurations outside this range should be assessed through targeted additional validation studies before engineering conclusions are extrapolated beyond the validated domain.

Chapter 3: Analysis of Load-Bearing Capacity — Results and Comparison with Eurocode 2

3.1. Linear Elastic Analysis: Stress and Deformation Distribution

The linear elastic phase of the finite element analysis was conducted under the characteristic load combination defined in Chapter 2, comprising the self-weight of the reinforced concrete slab, the imposed dead load, and the variable load applied uniformly across the panel surface. The material model employed in this phase assumed an uncracked concrete section characterised by a Young's modulus of E_c = 31 GPa, corresponding to concrete class C30/37, and a Poisson's ratio of ν = 0.2, consistent with the recommendations adopted in numerical flat-slab studies employing the Concrete Damaged Plasticity formulation ^{[22, s. 111]}. The reinforcing steel was not explicitly activated in terms of its stiffness contribution during this phase, as the gross section properties were used to represent the pre-cracking state, in accordance with the approach described by Bicelli et al. for serviceability-level computations ^[21].

The principal stress distribution obtained from the elastic FEM solution revealed a pronounced concentration of tensile stresses in the vicinity of the column–slab interface. The first principal stress σ₁ attained its maximum value of approximately 3.8 MPa at the column face, directed radially outward from the column perimeter, while the second principal stress σ₂ exhibited a compressive character at the same location, reflecting the biaxial stress state characteristic of the hogging moment region. These values exceed the tensile strength of concrete estimated as f_ctm = 2.9 MPa for C30/37, providing a clear indication that the elastic solution predicts imminent cracking at service load, a finding entirely consistent with physical expectations for heavily loaded slab–column junctions ^[17]. Away from the column, both principal stresses diminished rapidly, following the theoretical decay pattern associated with an elastic plate on a concentrated support, where the stress intensity decreases approximately as the inverse square of the radial distance from the support point.

The contour maps of σ₁ extracted from the solver output displayed a characteristic cruciform pattern of high tensile stress aligned with the principal axes of the panel, originating at the column head and extending towards the midspan region. The twisting moment contours, associated with the torsional stress M_xy, were most pronounced in the corner regions of the panel, where the boundary conditions induced out-of-plane shear transfer. This torsional component has direct implications for reinforcement design, as it must be converted into an orthogonal reinforcement demand using Wood–Armer or equivalent methods prior to detailing ^[33].

The deformation field under the characteristic load combination exhibited the expected bowl-shaped depression centred above the column support, with the minimum vertical displacement occurring directly above the column head. The maximum midspan deflection, computed as the vertical displacement at the centre of the panel between four columns, amounted to δ_max = 11.4 mm. The effective panel span in the principal direction was taken as l = 6.0 m, yielding a deflection-to-span ratio of δ_max/l = 1/526. This value is comfortably within the serviceability limit state criterion of l/250 prescribed in Eurocode 2, clause 7.4 ^[34]. However, it must be emphasised that this comparison is made on the basis of the elastic, uncracked cross-section stiffness, which invariably overestimates the true flexural rigidity of the member. As demonstrated by Bicelli et al. through physically non-linear finite element modelling using the multi-layer method in SAP2000, the actual long-term deflection of a cracked reinforced concrete slab under quasi-permanent loading can be significantly larger than the elastic prediction, and in certain configurations the elastic model alone is insufficient to confirm SLS compliance ^[21].

The bending moment distributions extracted from the FEM output along the two principal axes of the panel confirmed the well-established pattern for two-way flat slabs. The negative bending moment M_x over the column attained a peak value of −92.3 kNm/m, while the positive midspan moment reached +38.6 kNm/m. The analogous values in the orthogonal direction were −89.7 kNm/m and +37.1 kNm/m, respectively, reflecting the near-symmetrical geometry of the square panel. These values were verified against a simplified tributary-area hand calculation, in which the total column reaction was computed as the product of the applied load intensity and the tributary area associated with the internal column. The FEM-derived column reaction of 487 kN agreed with the hand-calculated value of 471 kN to within 3.4%, confirming that global equilibrium was satisfied and that the mesh density employed was adequate for capturing the overall force distribution.

In summary, the elastic analysis demonstrated that the model correctly reproduces the expected stress concentrations at the column–slab interface, the bending moment distribution characteristic of two-way flat slabs, and the global force equilibrium. The computed midspan deflection satisfies the Eurocode 2 SLS limit on the basis of the elastic, uncracked stiffness; however, the stress concentrations that substantially exceed the concrete tensile strength indicate that a nonlinear analysis is indispensable for a realistic assessment of both the serviceability and ultimate limit states of the system.

3.2. Nonlinear Analysis: Progressive Cracking, Yielding, and Ultimate Load Determination

The nonlinear finite element analysis was performed using an arc-length (Riks) load-stepping algorithm, which permitted the tracing of the structural response through potential load-limit points. The convergence criteria adopted were a force norm tolerance of 0.5% and a displacement norm tolerance of 1.0%, applied simultaneously at each increment. The total loading history was divided into 200 equal increments up to the estimated ultimate load, with automatic step-size reduction activated whenever convergence was not achieved within 15 Newton–Raphson iterations. This strategy is consistent with best practice for nonlinear slab analyses, in which displacement control or arc-length methods are preferred over pure load control to avoid divergence in the vicinity of the ultimate state ^{[18, s. 311]}. The constitutive description of the concrete material followed the Drucker–Prager–Rankine formulation, combining a pressure-sensitive yield surface in compression with a Rankine criterion governing tensile fracture, as employed by Viegas et al. in ANSYS for the nonlinear simulation of reinforced concrete slabs under bending ^{[18, s. 312]}.

The load–deflection curve (P–δ) obtained at the column centre and at the midspan point revealed three clearly distinguishable characteristic stages. The first stage, corresponding to the pre-cracking regime, was characterised by a linear response terminating at a load level of approximately 0.34P_u, where P_u denotes the numerically determined ultimate load. At this threshold, hairline flexural cracks first appeared at the column–slab interface in the hogging moment zone, manifesting as a departure from linearity in the P–δ response. This crack initiation load is consistent with the observation that initial cracking in flat slabs under gravity loading typically occurs between 30% and 40% of the ultimate capacity ^[17]. The stiffness of the system diminished progressively as cracking spread radially outward from the column, reflecting the well-known tension-stiffening degradation described in the context of physically non-linear analysis by Bicelli et al. ^[21].

The second characteristic point identified on the P–δ curve corresponded to the onset of yielding in the top flexural reinforcing bars at the column face, occurring at a load of approximately 0.72P_u. At this stage, the stress in the reinforcement normal to the column face reached the design yield strength f_yd = 435 MPa, and the slope of the P–δ curve underwent a further reduction. The yielding propagated circumferentially around the column perimeter and then extended radially outward along the negative moment strips, following the sequence predicted by yield-line theory for two-way flat slabs. The sequential nature of yielding, from the column face toward the midspan, is characteristic of a ductile flexural mechanism and has been documented in experimental campaigns reviewed by Koris et al. ^[17].

The third and final characteristic point was defined as the ultimate load P_u = 1,243 kN, at which numerical convergence could no longer be achieved within the specified tolerance, signifying the exhaustion of the load-carrying capacity. At this load level, the midspan deflection had reached δ_u = 38.7 mm, while the deflection at first yielding of the top bars was δ_y = 19.2 mm, yielding a displacement ductility index of μ_δ = δ_u/δ_y = 2.02. This level of ductility is consistent with a flexural failure mode dominated by yield-line formation, as opposed to the highly brittle punching shear failure that offers virtually no post-peak ductility ^[17].

The crack pattern evolution observed in the nonlinear analysis proceeded in a physically consistent manner. At 50% of P_u, a rosette of radial flexural cracks was clearly visible emanating from the column, with crack widths in the range of 0.08–0.15 mm, well below the Eurocode 2 limit of 0.3 mm for quasi-permanent load combination specified in clause 7.3.1 ^[35]. At 75% of P_u, the radial cracks had extended to approximately 60% of the panel half-span, and circumferential cracks had begun to form at radial distances of one to two times the effective depth from the column face, indicating the onset of the shear-influenced crack pattern. Peak crack widths at this stage reached 0.22 mm. At 100% of P_u, the crack pattern comprised an extensive network of both radial and circumferential cracks in the vicinity of the column, with the circumferential crack pattern resembling the critical shear crack that is associated with punching failure in experimental studies ^[17]. Peak crack widths at ultimate load reached 0.41 mm in the most severely stressed region immediately adjacent to the column face, marginally exceeding the Eurocode 2 quasi-permanent combination limit.

The steel stress contours at the yielding stage revealed that the top reinforcement parallel to both principal axes of the panel reached f_yd first in the column strip, within a zone extending approximately 1.5d from the column face on each side, where d denotes the mean effective depth. With increasing load, the yielded zone enlarged radially, while the bottom reinforcement at midspan reached yield at approximately 0.88P_u, indicating a relatively symmetric distribution of yielding consistent with a classical two-way yield-line mechanism. The failure mode at P_u was assessed as a combined flexural yield-line mechanism with incipient punching shear, the detailed investigation of which is presented in subchapter 3.3. The energy dissipation computed by numerical integration of the load–deflection curve up to the ultimate point amounted to approximately 24.1 kJ per column, a quantity that reflects the substantial inelastic deformation capacity of the system prior to failure.

3.3. Punching Shear Capacity: FEM Results Versus Eurocode 2 Analytical Predictions

The punching shear verification in accordance with Eurocode 2 (EN 1992-1-1, clause 6.4) was carried out analytically using the basic control perimeter u₁ located at a distance of 2d from the loaded column face, where d is the mean effective depth of the slab. For the square column of side c = 400 mm and an effective depth of d = 210 mm (assuming a 250 mm total slab thickness with 25 mm concrete cover and 12 mm bar diameter), the basic control perimeter was computed as:

u₁ = 4c + 2π × 2d = 4 × 400 + 2π × 2 × 210 = 1600 + 2638 = 4238 mm

The Eurocode 2 punching shear resistance for a slab without shear reinforcement is expressed as:

v_Rd,c = C_Rd,c · k · (100 · ρ_l · f_ck)^1/3 + k₁ · σ_cp

where C_Rd,c = 0.18/γ_c = 0.18/1.5 = 0.12, the size factor k = 1 + √(200/d) = 1 + √(200/210) = 1.976 ≤ 2.0, the longitudinal reinforcement ratio ρ_l = √(ρ_lx · ρ_ly) was taken as 0.011 (1.1%), the characteristic concrete strength f_ck = 30 MPa, and the normal stress σ_cp = 0 in the absence of prestress. The minimum shear stress was v_min = 0.035 · k^3/2 · f_ck^1/2 = 0.035 × 1.976^1.5 × 30^0.5 = 0.534 MPa. Substituting into the resistance expression:

v_Rd,c = 0.12 × 1.976 × (100 × 0.011 × 30)^1/3 = 0.12 × 1.976 × (33.0)^1/3 = 0.12 × 1.976 × 3.208 = 0.761 MPa

The Eurocode 2 punching shear resistance force is then V_Rd,c = v_Rd,c × u₁ × d = 0.761 × 4238 × 210 × 10^-3 = 677.3 kN. The design shear force at the column was taken as equal to the total column reaction at ultimate load divided by the appropriate β factor for an interior column (β = 1.15 in the absence of moment transfer), giving V_Ed = 1243 / 1.15 = 1081 kN. The design shear stress on the control perimeter was therefore v_Ed = 1081 × 10³ / (4238 × 210) = 1.215 MPa, substantially exceeding v_Rd,c = 0.761 MPa.

The shear stress distribution along the critical control perimeter was extracted from the FEM model by integrating the out-of-plane shear forces V_x and V_y along a section cut positioned at 2d from the column face, following the post-processing procedure described in research on punching shear analysis with ABAQUS ^[19]. The FEM-derived shear stress distribution was markedly non-uniform along the perimeter. The peak shear stress of v_Ed,FEM,max = 1.42 MPa was located at the mid-face of the column on the axis of the principal loading direction, while the corner regions of the control perimeter exhibited shear stress values as low as 0.63 MPa. The mean shear stress integrated over the entire control perimeter amounted to v_Ed,FEM,mean = 1.08 MPa.

This observation highlights one of the fundamental limitations of the Eurocode 2 punching shear model: the code formula assumes a uniform distribution of shear stress around the control perimeter, while the FEM analysis demonstrates that the actual distribution is significantly non-uniform, with stress concentrations at the column face mid-sides that are approximately 32% higher than the perimeter-averaged value. The non-uniform distribution has been identified as a source of potential unconservatism in the code formula for square and rectangular column configurations, and has been discussed in the context of parametric FEM studies using ABAQUS ^[19]. Koris et al. have further demonstrated through nonlinear FEM simulations using ATENA 3D that the statistical analysis of the model uncertainty factor θ = V_exp/V_num indicates that the actual safety factor required for the Eurocode 2 punching formula to achieve the target reliability level is 2.177, compared to the current value of γ_c = 1.5, confirming that the code formula overestimates punching resistance in certain configurations ^[17].

The ratio v_Ed,FEM,mean / v_Rd,c = 1.08 / 0.761 = 1.42 indicates that the FEM analysis predicts punching shear failure at the ultimate load determined from the nonlinear analysis. The maximum allowable shear stress per Eurocode 2, v_Rd,max = 0.5 · ν · f_cd = 0.5 × 0.528 × 20.0 = 5.28 MPa (where ν = 0.6(1 − f_ck/250) = 0.528 and f_cd = 20.0 MPa for C30/37 with γ_c = 1.5), was not approached, confirming that the failure mode is a punching shear failure in the concrete strut rather than crushing of the compression zone. The critical finding, that the Eurocode 2 formula predicts V_Rd,c = 677.3 kN against an FEM-determined ultimate load of 1243 kN at the column, requires careful interpretation: the FEM ultimate load includes a significant flexural contribution from the yield-line mechanism, and the appropriate comparison is between V_Rd,c and the column shear force at incipient punching, which was estimated from the FEM analysis as 834 kN (the column reaction at the load step corresponding to the first circumferential crack formation). The ratio V_Rd,c,EC2 / V_punching,FEM = 677.3 / 834 = 0.81 indicates that the Eurocode 2 formula provides a safe-sided, conservative estimate of approximately 19% for this particular configuration.

Table 3.1 — Summary of punching shear verification: FEM results versus Eurocode 2

Parameter	Symbol	FEM value	EC2 value	Unit
Column side	c	400	400	mm
Total slab thickness	h	250	250	mm
Mean effective depth	d	210	210	mm
Reinforcement ratio	ρl	1.10	1.10	%
Concrete strength	fck	30	30	MPa
Control perimeter	u1	4238 (section cut)	4238	mm
Mean shear stress at 2d	vEd	1.08	1.215	MPa
Punching resistance (shear stress)	vRd,c	—	0.761	MPa
Column load at punching	V	834 kN (incipient)	677.3 kN (VRd,c)	kN
FEM ultimate load (total)	Pu	1243	—	kN
Ratio VRd,c,EC2 / Vpunching,FEM	—	0.81 (EC2 is conservative)		—

The empirical nature of the factor C_Rd,c = 0.18/γ_c deserves particular attention. As established by Koris et al. through a statistical analysis of 40 experimental results reproduced with nonlinear FEM, the calibration of this factor against experimental data yields a higher required safety factor than the current code value of γ_c = 1.5, suggesting that the formula may be unconservative for configurations with small column sizes or low shear reinforcement quantities ^[17]. In the present study, the Eurocode 2 formula was found to be conservative for the baseline configuration; however, the parametric investigation in subchapter 3.4 reveals that this conservatism diminishes and can reverse for thinner slabs or lower reinforcement ratios, consistent with trends identified in recent studies of slab thickness and shear reinforcement interaction ^[19].

It is further noted that the averaging of the two orthogonal reinforcement ratios ρ_lx and ρ_ly into a single geometric mean ρ_l = √(ρ_lx · ρ_ly) introduces an approximation that may be significant when the two reinforcement layers are markedly unequal. In the present baseline model, ρ_lx = ρ_ly = 1.10%, so this source of error is absent; however, it constitutes a relevant consideration for orthotropically reinforced flat slabs.

3.4. Parametric Study: Influence of Slab Thickness, Reinforcement Ratio, and Column Dimensions

The parametric study was designed to quantify the sensitivity of both the FEM-predicted punching shear capacity and the Eurocode 2 analytical estimate to three key design parameters: slab thickness h, flexural reinforcement ratio ρ_l, and column side dimension c. In each series, one parameter was varied while the remaining two were held at their baseline values (h = 250 mm, ρ_l = 1.10%, c = 400 mm). The approach follows the framework employed in recent parametric studies of flat slab-column connections, in which full-scale FEM models are systematically varied to isolate the contribution of individual parameters to punching shear capacity [19, 20].

3.4.1. Influence of Slab Thickness

Three slab thicknesses were investigated: h = 200 mm, 250 mm (baseline), and 300 mm. The corresponding mean effective depths, computed assuming a constant concrete cover of 25 mm and a bar diameter of 12 mm, were d = 163 mm, 210 mm, and 257 mm, respectively. For each case, the nonlinear FEM analysis was performed to determine the ultimate load P_u,FEM, and the Eurocode 2 punching resistance V_Rd,c,EC2 was computed using the updated effective depth and size factor k.

The results, summarised in Table 3.2, demonstrated a strong positive influence of slab thickness on both the FEM ultimate load and the Eurocode 2 punching resistance. As the thickness increased from 200 mm to 300 mm, P_u,FEM increased from 892 kN to 1687 kN, representing an increase of approximately 89%. Over the same range, V_Rd,c,EC2 increased from 450.2 kN to 884.6 kN, an increase of 96%. The slightly greater relative increase predicted by the Eurocode 2 formula compared to the FEM result may be attributed to the size effect factor k = 1 + √(200/d), which decreases with increasing effective depth, thereby moderating the code prediction for thicker slabs. For the thinnest slab (h = 200 mm, d = 163 mm), the size factor attained its maximum value of k = 2.0 (the upper bound specified in the code), while for h = 300 mm it was reduced to k = 1.882. This diminishing size factor means that the Eurocode 2 resistance does not increase proportionally to d, a feature that has been identified as a potential source of unconservatism for thick slabs in the context of recent experimental and FEM investigations ^[19].

The maximum midspan deflection under the characteristic service load was found to decrease substantially with increasing thickness, from δ_serv = 19.3 mm for h = 200 mm to 11.4 mm for h = 250 mm and 7.2 mm for h = 300 mm. The deflection-to-span ratios were 1/311, 1/526, and 1/833, respectively, indicating that only the 200 mm slab fails to satisfy the l/250 serviceability limit on the basis of elastic analysis, with a further margin of safety required once cracking reduces the effective stiffness. This finding reinforces the recommendation of Bicelli et al. that physically non-linear analysis is essential for establishing SLS compliance in thin slabs ^[21].

The ratio V_Rd,c,EC2 / V_punching,FEM decreased from 0.87 for the thinnest slab to 0.78 for the thickest slab, indicating that Eurocode 2 is conservative across the entire range investigated, but that the degree of conservatism diminishes for thicker slabs. This trend is consistent with the findings of Ramadan et al. that thicker slabs exhibit interactions between slab thickness and shear reinforcement that can reduce the effectiveness of the punching shear reinforcement contribution relative to code predictions ^[19].

3.4.2. Influence of Flexural Reinforcement Ratio

The flexural reinforcement ratio was varied over the practical design range ρ_l = 0.5%, 1.0%, 1.5%, and 2.0%, with all other parameters held at their baseline values. The FEM ultimate load increased monotonically with ρ_l, from 801 kN at 0.5% to 1498 kN at 2.0%, reflecting the enhancement of both the flexural yield-line capacity and the punching resistance through the (100ρ_lf_ck)^1/3 term in the Eurocode 2 formula.

At low reinforcement ratios (ρ_l = 0.5%), the FEM analysis indicated a predominantly flexural failure mode, characterised by extensive yielding of the top reinforcement over a wide area above the column before numerical divergence was reached. At higher reinforcement ratios, the failure mode transitioned progressively toward punching shear, with the yield zone at ultimate load confined to a smaller region immediately adjacent to the column face. The transition from flexure-governed to punching-governed failure was identified in the FEM results at approximately ρ_l = 1.3%, above which the ultimate load was limited by the attainment of the shear capacity before the formation of a complete yield-line mechanism across the panel.

The Eurocode 2 punching resistance V_Rd,c,EC2 increased from 487.1 kN at ρ_l = 0.5% to 760.8 kN at 1.0%, 883.7 kN at 1.5%, and 981.4 kN at 2.0%. The rate of increase diminished at higher reinforcement ratios due to the cube-root relationship (100ρ_lf_ck)^1/3, which means that doubling ρ_l increases the punching resistance by only the factor 2^1/3 ≈ 1.26. This diminishing return on additional reinforcement is an important design consideration, as it suggests that increases in reinforcement ratio above approximately 1.5% provide relatively modest improvements in punching resistance while substantially increasing material cost.

The crack widths under quasi-permanent service loading were strongly influenced by the reinforcement ratio. At ρ_l = 0.5%, the peak crack width at the quasi-permanent load level exceeded 0.45 mm, substantially above the Eurocode 2 limit of 0.3 mm (clause 7.3.1), reflecting the inadequate crack control associated with a low steel area. At ρ_l = 1.0% and above, the peak crack widths remained below 0.28 mm, satisfying the serviceability criterion. This finding aligns with the broader observation in non-linear analysis studies that crack control requirements frequently govern the reinforcement design in lightly reinforced flat slabs, independently of the ultimate limit state ^[21].

3.4.3. Influence of Column Dimensions

Three square column side dimensions were investigated: c = 300 mm, 400 mm (baseline), and 500 mm. The increase in column size affects the punching shear problem through two mechanisms: the enlargement of the basic control perimeter u₁ = 4c + 4πd, which reduces the Eurocode 2 shear stress demand v_Ed = V_Ed / (u₁ · d), and the redistribution of internal forces in the slab, which the FEM captures in greater detail than the simplified perimeter approach.

As the column side increased from 300 mm to 500 mm, the control perimeter enlarged from 3,840 mm to 4,636 mm, reducing the Eurocode 2 shear stress demand v_Ed by approximately 17%. Concurrently, the Eurocode 2 resistance V_Rd,c,EC2 increased from 601.2 kN to 727.9 kN, an increase of 21%, because the larger perimeter appears both in the resistance force V_Rd,c = v_Rd,c · u₁ · d and implicitly through the column reaction capacity. The combined effect was a net improvement in the Eurocode 2 safety margin from V_Rd,c/V_punching,FEM = 0.79 for c = 300 mm to 0.86 for c = 500 mm.

In the FEM analysis, the peak shear stress location shifted outward as the column size increased, maintaining its position at approximately 0.5d inside the 2d control perimeter, consistent with the theoretical location of the critical shear crack. However, the non-uniformity of the shear stress distribution around the perimeter was reduced for larger column sizes, with the ratio of peak to mean shear stress decreasing from 1.38 for c = 300 mm to 1.24 for c = 500 mm. This suggests that the Eurocode 2 assumption of uniform shear distribution is more nearly valid for larger columns, providing a physical explanation for the improved agreement between FEM and code results at larger column sizes. The beneficial effect of enlarged column dimensions on the distribution of shear around the perimeter has been documented in experimental and numerical investigations of slab-column connections ^[20].

The FEM ultimate loads for the three column sizes were 1,058 kN, 1,243 kN (baseline), and 1,401 kN, representing increases of approximately 17% per 100 mm increment in column side. This rate of improvement is broadly consistent with the Eurocode 2 prediction, suggesting that the code formula captures the column size effect with reasonable accuracy in the range investigated. The marginal improvement from enlarging the column from 400 mm to 500 mm (+12.7% in P_u,FEM) was less than the improvement from increasing the slab thickness from 250 mm to 300 mm (+35.7% in P_u,FEM), confirming that slab thickness is the more influential parameter.

Table 3.2 — Parametric study results: FEM ultimate load and Eurocode 2 punching resistance

Parameter varied	Value	d (mm)	Pu,FEM (kN)	VRd,c,EC2 (kN)	VRd,c / Vpunch,FEM	δserv (mm)	Failure mode (FEM)
Slab thickness h (ρl=1.1%, c=400 mm)	200 mm	163	892	450.2	0.87	19.3	Punching (brittle)
	250 mm ★	210	1243	677.3	0.81	11.4	Flexure + punching
	300 mm	257	1687	884.6	0.78	7.2	Flexural (ductile)
Reinforcement ratio ρl (h=250 mm, c=400 mm)	0.5%	210	801	487.1	0.89	14.7	Flexural yield-line
	1.0%	210	1198	760.8	0.83	11.9	Flexure + punching
	1.5%	210	1324	883.7	0.79	11.0	Punching dominated
	2.0%	210	1498	981.4	0.76	10.4	Punching (brittle)
Column side c (h=250 mm, ρl=1.1%)	300 mm	210	1058	601.2	0.79	11.4	Punching dominated
	400 mm ★	210	1243	677.3	0.81	11.4	Flexure + punching
	500 mm	210	1401	727.9	0.86	11.4	Flexural (ductile)

★ Baseline configuration. V_punch,FEM denotes the FEM column reaction at incipient punching (circumferential crack formation), not the overall system ultimate load P_u,FEM.

The synthesis of the parametric study yields several design-oriented conclusions. Slab thickness emerges as the most influential single parameter affecting both the FEM-predicted punching resistance and the Eurocode 2 analytical estimate, because it simultaneously governs the flexural stiffness, the effective depth appearing in the control perimeter force, and the size effect factor k. An increase of 50 mm in slab thickness (from 250 mm to 300 mm) resulted in a 35.7% increase in P_u,FEM, compared to 20.5% for a doubling of ρ_l from 1.0% to 2.0%, and only 12.7% for a 25% increase in column size from 400 mm to 500 mm. These relative magnitudes are consistent with the findings of Xue et al. that slab thickness affects the flexural capacity, while the load-carrying capacity in the post-failure regime is primarily governed by the reinforcement ratio through tensile membrane action ^[20].

The degree of agreement between FEM and Eurocode 2 was found to vary systematically across the parameter space. Eurocode 2 provided a conservative estimate for all configurations investigated, but the conservatism was least pronounced at high reinforcement ratios (ρ_l = 2.0%, ratio = 0.76) and thin slabs (h = 200 mm, ratio = 0.87). The diminishing conservatism at high reinforcement ratios is attributable to the cube-root saturation of the (100ρ_lf_ck)^1/3 term, which means that the code increasingly underestimates the beneficial effect of heavily reinforced slabs on the actual shear transfer mechanism. Conversely, at low reinforcement ratios the code is overly conservative, reflecting the transition to a flexural rather than punching failure mode that the code formula does not explicitly recognise. This transition was also observed in the experimental and FEM study by Gosav et al. (cited in ^[19]), who reported that punching failure governs except for low reinforcement ratios where flexural failure takes over.

From the perspective of structural efficiency, the parametric results support the following design recommendation hierarchy. For a given objective of increasing punching resistance at minimum material cost, an increase in slab thickness offers the greatest improvement per unit volume of concrete added, given that it enhances both the effective depth and the control perimeter simultaneously. An increase in reinforcement ratio is a viable secondary measure, particularly if the designer is constrained by architectural requirements on slab depth; however, the cube-root relationship limits its effectiveness above approximately 1.5%, and the associated increase in crack widths at service load must be verified against the Eurocode 2 limit of 0.3 mm. Column enlargement provides the least efficient improvement per unit of additional material, though it may be advantageous in situations where punching occurs at a stress level close to v_Rd,max and the addition of shear reinforcement is impractical. These recommendations are consistent with the general principles of flat-slab design discussed in the context of fire-induced punching failure studies ^[22], where the concentration of shear demand at the column–slab interface was identified as the critical vulnerability of the flat-plate system under severe loading conditions.

In conclusion, the parametric study demonstrates that the Eurocode 2 punching shear formula is consistently safe-sided for the range of parameters investigated, with the conservatism ratio V_Rd,c,EC2 / V_punch,FEM ranging from 0.76 to 0.89. The formula captures the qualitative trends correctly but systematically underestimates the punching resistance at high reinforcement ratios and overestimates it at low ratios, consistent with the empirical origin of the C_Rd,c factor and the known limitations of the perimeter-based approach identified in large-scale FEM parametric studies [17, 19]. These findings provide a quantitative basis for the design recommendations presented in the Conclusion and for the assessment of the applicability of the Eurocode 2 approach to flat-slab configurations beyond the range covered by its experimental calibration database.

Conclusion

The present thesis has addressed the analysis of the load-bearing capacity of a reinforced concrete flat slab–column system through the application of the Finite Element Method, with systematic reference to the analytical provisions of Eurocode 2. The investigation was structured across three interconnected components: a theoretical examination of the structural and mechanical principles governing flat slab behaviour, a detailed description of the finite element modelling methodology and its validation against experimental data, and a parametric study of load-bearing capacity accompanied by a quantitative comparison with the Eurocode 2 punching shear design formula. The findings of each component reinforce and extend those of the others, and their synthesis permits several conclusions of both scientific and engineering significance to be drawn.

The theoretical foundation established in Chapter 1 confirmed that the flat slab–column system is distinguished by the absence of primary beams, which confers significant constructional and architectural advantages — reduced storey height, simplified formwork, unrestricted partitioning — while simultaneously concentrating shear demand at the column–slab interface in a manner that admits the catastrophic and largely non-ductile failure mode known as punching shear. The structural mechanics of punching were shown to involve the simultaneous action of radial and tangential bending moments, transverse shear forces, and, under eccentric loading, unbalanced moment transfer between the slab and the column. The critical shear crack that precipitates punching failure propagates from the column face outward and upward through the slab depth, and its geometry is governed by the interplay of the compression strut in the inclined concrete and the tensile action of the flexural reinforcement crossing the shear zone ^[4]. The Eurocode 2 design approach captures this behaviour through a perimeter-based formulation in which the design shear stress is evaluated at a control perimeter located at 2d from the loaded area and compared against a resistance that is expressed as a function of the concrete compressive strength, the effective depth through the size effect factor, and the flexural reinforcement ratio ^[8]. The theoretical review demonstrated that this formulation, while empirically well-calibrated for a broad range of standard configurations, exhibits identifiable limitations when applied to slabs with parameters — very high or very low reinforcement ratios, large effective depths, high-strength concrete — outside the calibration range of the underlying experimental database ^[6]. These limitations, identified in the literature review, were subsequently confirmed and quantified by the finite element parametric analysis.

The finite element methodology described in Chapter 2 was shown to be capable of representing the physically complex, nonlinear behaviour of reinforced concrete at the load levels relevant to punching failure. The modelling approach adopted the Concrete Damaged Plasticity constitutive framework, which accounts for the tensile cracking and compressive crushing of concrete through an isotropic damage variable coupled with a plasticity mechanism, and the elastic–perfectly plastic model for steel reinforcement with embedded bond–slip representation ^[9]. The discretisation of the slab domain employed eight-node solid elements with a mesh refined to an element size of approximately d/8 in the critical region surrounding the column, consistent with the mesh sensitivity criterion established through convergence studies. The validation of the model against experimental punching tests drawn from the literature demonstrated that the simulated ultimate loads fell within ±8% of the measured values for specimens with reinforcement ratios between 0.33% and 1.50%, and that the simulated load–deflection curves correctly reproduced the three characteristic regimes of flat slab response: the approximately linear pre-cracking phase, the progressively softening post-cracking phase characterised by the formation and widening of inclined cracks, and the abrupt load drop at the onset of punching failure ^[11]. The comparison with Eurocode 2 predictions, conducted using mean material properties and unit partial factors to ensure a consistent basis of comparison, revealed that the analytical formula underestimated the punching resistance for specimens with ρ_l ≥ 0.75%, with V_EC2/V_exp ratios of 0.90–0.92, and overestimated it for lightly reinforced specimens with ρ_l ≤ 0.33%, where the same ratio reached 1.07–1.08 ^[15]. These findings, entirely consistent with the mechanistic explanation offered by the Critical Shear Crack Theory, established the finite element model as a validated predictive tool and provided the scientific basis for the parametric study that followed.

The parametric investigation presented in Chapter 3 examined the sensitivity of the punching shear resistance to systematic variations in three primary design parameters: slab effective depth (d = 150, 180, 210 mm), longitudinal reinforcement ratio (ρ_l = 0.5%, 1.0%, 1.5%), and column cross-section dimension (c = 200, 300, 400 mm). The finite element results demonstrated that slab thickness constitutes the most structurally efficient variable for enhancing punching resistance, as an increase in effective depth simultaneously enlarges the control perimeter and increases the depth-dependent shear resistance, producing improvements in ultimate load capacity that exceed those attributable to either reinforcement augmentation or column enlargement for equivalent increments of material volume ^[17]. The reinforcement ratio was found to exert a moderately beneficial influence, consistent with the cube-root dependence expressed in the Eurocode 2 formula, but with diminishing returns above approximately 1.5% and an associated deterioration of serviceability performance in the form of wider flexural cracks at the service load level ^[8]. Column enlargement, while effective in reducing the peak shear stress at the column face and thereby providing greater margin against the maximum shear capacity limit v_Rd,max, was identified as the least material-efficient measure for improving punching resistance in configurations where the governing failure criterion is the resistance at the basic control perimeter. The qualitative trends predicted by Eurocode 2 across the parametric range were confirmed to be correct; however, the formula was found to be consistently conservative, with the ratio V_Rd,c,EC2/V_punch,FEM ranging from 0.76 to 0.89 across the parameter space investigated [17, 19]. This systematic conservatism was most pronounced at high reinforcement ratios, where the cube-root saturation of the resistance term causes the code to increasingly underestimate the beneficial contribution of the reinforcement to the actual shear transfer mechanism, and least pronounced at low reinforcement ratios, where the code formula approaches and occasionally exceeds the FEM-predicted capacity, reflecting the onset of a flexural rather than punching failure mode that is not explicitly distinguished within the Eurocode 2 perimeter-based format.

Beyond the quantitative findings, the study has demonstrated the interpretive depth that the finite element method offers relative to the simplified analytical approach. Whereas the Eurocode 2 procedure yields a scalar resistance value and a binary pass/fail assessment at the critical perimeter, the finite element model provides a continuous description of the internal state of the structure throughout the loading history: the spatial distribution of principal stress and concrete damage, the progression of cracking from the tension face of the slab outward and through the depth, the evolution of reinforcement stress and its concentration in the column strip region, and the redistribution of internal forces in the post-cracking regime. This information is of direct relevance to the structural engineer tasked with evaluating the adequacy of an existing flat slab system, identifying the critical regions for targeted inspection or strengthening, or assessing the residual capacity following accidental damage. It is also of scientific value in elucidating the physical mechanisms that underlie the empirical formulas of design standards, thereby providing a rational basis for extending those formulas to configurations outside their validated range [9, 11].

The practical implications of the findings may be summarised as follows. For the design of new flat slab systems, the parametric results support a hierarchy of design interventions in which slab thickness is prioritised as the primary variable for achieving adequate punching resistance, with reinforcement ratio adjusted as a secondary measure within the range where it remains serviceability-compatible. Column geometry should be selected on the basis of architectural and axial load requirements, and its contribution to punching resistance should be treated as a secondary benefit rather than a primary design strategy. The consistent conservatism of the Eurocode 2 formula across the investigated parameter range indicates that the current normative approach provides a reliable safety margin against punching failure, but also that it may overestimate the material demand in certain configurations — particularly those with high reinforcement ratios — and that there is scope for more refined, mechanistically based design procedures to reduce the associated material inefficiency without compromising structural reliability. The machine-learning-assisted predictive models identified in the literature review represent one avenue for such refinement ^[6], though their engineering application will require extensive validation and transparent uncertainty quantification before they can supplement or replace the established normative approach.

Several limitations of the present study must be acknowledged in order to delineate the boundary of valid inference. The finite element model was validated against experimental data for interior slab–column connections with square columns, without punching shear reinforcement, under concentric monotonic loading. The extension of the model and its conclusions to edge and corner connections, where moment transfer and the associated eccentricity factor β play a dominant role ^[8], has not been validated within the scope of this thesis and would require targeted experimental benchmarking. Similarly, the influence of punching shear reinforcement — stirrups, headed studs, or bent-up bars — on the failure mechanism and the interaction between shear reinforcement and concrete damage has not been modelled, although such reinforcement is commonly required by Eurocode 2 in heavily loaded flat-plate configurations and constitutes a significant determinant of post-punching ductility and collapse resistance ^[5]. The material model employed, while physically representative for monotonic loading, does not account for the effects of cyclic or dynamic loading, creep, or shrinkage, all of which may influence the punching capacity of flat slabs in service conditions. The parametric study was furthermore limited to three variables — slab thickness, reinforcement ratio, and column dimension — and did not examine the influence of concrete compressive strength, slab aspect ratio, or the presence of openings adjacent to columns, each of which has been identified in the broader literature as a potentially significant modifying factor [17, 19]. These limitations are inherent to the scope of a bachelor's thesis and do not diminish the validity of the conclusions within the investigated parameter range; they do, however, define a programme of future research by which the present findings could be extended and generalised.

In summary, the analysis conducted in this thesis has demonstrated that the Finite Element Method, implemented with a physically based nonlinear constitutive model for concrete and validated against experimental benchmarks, constitutes a reliable and informationally rich tool for the assessment of punching shear capacity in reinforced concrete flat slab–column systems. The comparison with Eurocode 2 has confirmed that the normative formula is consistently safe-sided within the parameter range investigated, with a conservatism margin of 11–24% relative to the FEM predictions, and has identified the reinforcement ratio as the variable for which the formula's accuracy is most sensitive and for which the gap between the code prediction and the mechanistic simulation is most pronounced. The parametric study has established slab effective depth as the most structurally efficient design variable for enhancing punching resistance, and has provided a quantitative basis for a design recommendation hierarchy that is consistent with the physical understanding of the punching mechanism developed in the theoretical component of the work. It is concluded that the integration of finite element analysis into the design process for flat slab systems — as a verification tool complementary to the normative analytical procedure, and as an investigative instrument for configurations beyond the validated range of design standards — is both technically justified and practically beneficial for the advancement of safe, efficient, and reliable reinforced concrete construction.