Multidimensional nonlinear numerical simulation of post‐tensioned concrete girders with different prestressing levels

Most of the Italian infrastructures have been built since the mid‐twentieth century. Prestressed concrete (PC) has been widely used, especially for bridges and viaducts. These structures have proved effective over time, but their performance could be affected by degradation phenomena (i.e., corrosion of the tendons and residual prestressing) or construction defects (i.e., grouting of the ducts), which need to be accurately modeled. This paper focuses on the flexural response of two reduced‐scale, posttensioned, PC bridge girders, prestressed with two different jacking forces and tested under four‐point bending configuration. Four multidimensional numerical models were developed to simulate the experimental behaviors of the selected specimens, using two alternative computational strategies, namely, the finite element method (FEM) and applied element method (AEM). Three FEM modeling were considered, ranging from one‐dimensional lumped‐plasticity models to two‐dimensional and three‐dimensional (3D) spread‐plasticity models. Besides, 3D AEM models were developed into another software package. The accuracy and the computational costs of the two numerical strategies are discussed in this paper. Also, the main features of each numerical modeling technique are addressed, including comparisons between numerical and experimental global response data as well as the statistical characterization of the model error. Based on these different features, the authors also suggest the most suitable numerical strategy for future studies on PC girders.

longitudinal and cross girders of bridge decks.In many countries such as Italy, Germany, France, and United States of America, the road network consists of thousands of PC bridge decks, mostly built starting from the 1950s. 2,3For instance, the Italian bridges stock consist of more than 90% of PC structure dating back to the period between 1960 to 1980. 4 Therefore, assessing the residual load-carrying capacity of aging existing bridges has gradually become an important task. 5,6xperimental tests are the most straightaway method to study both flexural and shear behavior of these structures.5][16][17] These works demonstrated how existing PC bridges with none-to-slight damage or deterioration ensured satisfactory load-bearing capacity, overperforming structural requirements by code provisions.][21][22] In the last decades, numerical methods were adopted to study the response of PC structures focusing on different aspects.Based on the experimental tests on PC samples, the well-known finite elements (FEs) software package Abaqus FEA ® had shown to be able to predict not only the structural behavior at both system and sectional scale (i.e., force-displacement and moment-curvature relationships 23,24 ) but also the local response parameters (i.e., tensile cracking, release of the strands, end anchorages [25][26][27][28] ) of PC members with good accuracy.Experimentally validated FEs codes were also developed to account for long-term effects (creep, shrinkage, and relaxation) on the flexural response of PC beams with either bonded 29 and unbonded tendons. 30,31umerical models of PC structures with unbonded tendons were also proposed.Some authors compared the global response of PC beams with bonded and internally unbonded tendons, finding significant differences in terms of cracking pattern as well as ultimate load capacity and deflection. 32A parametric finite element analysis (FEA) was also carried out to investigate the influence of the prestressing levels on the flexural response of beams with unbonded tendons, revealing that higher effective prestress leads to higher load-carrying capacity at ultimate stage even if with lower deflection (i.e., reduced ductility). 33Some attempts to simplify the complex, member-dependent, stress-strain distribution along unbonded tendons assuming them as constants were also proposed using both FEs 34,35 or applied elements (AEs) 36 methods.However, to the authors' best knowledge, a comprehensive and comparative multidimensional numerical study on the flexural behavior of PC girders is still missing.
This paper investigates the flexural behavior and the load-bearing capacity of two 1/5 scaled posttensioned PC girders with two different levels of initial prestressing.The girders were tested under quasi-static four-point bending up to peak vertical force.Numerical models of the girder specimens were developed according to different computational strategies and using three software packages, in order to accurately reproduce the experimental response and observed damage.Two numerical approaches were implemented, namely the finite element method (FEM) and the applied element method (AEM), using the following types of elements: one-dimensional (1D), two-dimensional (2D), and three-dimensional (3D) FEs and 3D AEs.Comparisons between numerical results and experimental outcomes are herein presented in terms of global (e.g., response curve) and local (e.g., softening response of the materials) response of the girders and their damage, to highlight the pros and cons of the different methods, as well as to characterize the statistical model error for their accuracy assessment and possible probabilistic simulations.The validation of the different numerical models under study against experimental data creates a strong basis for structural assessment of PC girders in existing simply supported bridges.

| Description of specimens
The experimental program deals with quasi-static four-point bending tests on two girder specimens, which were designed in a reduced scale compared to an existing PC bridge, 37 using a geometric scale factor S L = 5.Assuming an acceleration scale S A = 1, a consistent model is obtained with a mass scale equal to S 3 L ¼ 125 and the girder self-weight reduced accordingly compared to the prototype bridge. 38he specimens share the same geometry: a T-shaped cross section with a total length of 6600 mm and a total height of 440 mm (Figure 1).The size of the flange and web was equal to 480 Â 60 and 150 Â 380 mm 2 , respectively.
Mild steel rebars consisted of three levels of two 8-mm diameter longitudinal bars, vertical stirrups with 8-mm nominal diameter and 100-mm spacing, and a welded wire mesh in the flange with nominal diameter and mesh size equal to 6 and 200 mm, respectively.
Two seven-wire mono-strand tendons, each of them having a gross area of 150 mm 2 , were used to prestress the girders after concrete curing.Each tendon had a symmetrical layout with respect to the entire specimen, with 1000-mm-long straight configuration in the central part of the girder and 5600-mm-long parabolic configuration in the lateral parts of the girder.
The two specimens were designed with two different levels of initial prestress.The jacking force of the specimens was set equal to P HP = 150 kN (high prestress [HP], corresponding to σ p,HP ¼ 1000 MPa) or P HP ¼ 75 kN (low prestress [LP], corresponding to σ p,HP ¼ 500 MPa) for both tendons of the same girder in order to provide different levels of prestress.Accordingly, the two specimens are herein identified as S1-HP (Specimen 1, HP) and S2-LP (Specimen 2, LP), respectively.
Prestressing losses at the time of testing due to initial deformation of the strand anchorages, friction, creep, and shrinkage phenomena were measured using load cells at the end of the girders.The effective prestressing forces (P e ) in bottom and upper tendons were found equal to 111 kN (percentage reduction of À26%) and 104 kN (À31%) in S1-HP specimen (P e,HP ), and 36.3 kN (À52%) and 33.4 kN (À55%) in S2-LP specimen (P e,LP ), respectively.
For each girder, six cubic concrete samples were casted to be tested up to failure, providing a mean cubic compressive strength R cm = 30 MPa at 28 days of curing.Considering the concrete curing between casting and testing of samples (i.e., 228 and 263 days for S1-HP and S2-LP, respectively), the cylindrical compressive strength at the day of testing (f cm (t)) was estimated equal to 30 MPa according to EN-1992-2-1. 39This value is in line with a recent study from these authors based on data mining of mechanical properties from existing PC bridge decks. 40Regarding prestressing steel, the stress at 1% of total strain f p1 , ultimate tensile stress f pu , and fracture elongation ε pu were found to be 1782 MPa, 1969 MPa, and 0.07, respectively, whereas Young's modulus E p ¼ 204, 215 MPa, as per minimum requirements confirmed by acceptance tests provided by the supplier. 41he mild steel was designed as per B450C class according to the Italian Building Code, 42 with nominal values of yield strength f sy ¼ 450 MPa, ultimate tensile strength f su ¼ 540 MPa, and Young's modulus E s ¼ 210, 000 MPa.

| Test setup
The experimental setup under the universal loading machine and one of the specimens before testing at the laboratory of the Department of Structures for Engineering and Architecture at University of Naples Federico II are shown in Figure 2a,b, respectively.Each girder specimen was quasi-statically tested under four-point bending.The loading procedure consisted of two phases, as follows: a cyclic load protocol with force control (Phase 1) and monotonic load protocol with displacement control and constant rate of 0.05 mm/s up to failure (Phase 2).The cyclic protocol was designed to investigate the performance of PC girders under service loads, whereas the monotonic protocol was used to assess the ultimate response. 38In Phase 1, the specimen attained three increasing levels of peak vertical force, that is, F e,1 = 33.1 kN (serviceability level), F e,2 = 48.3kN (design level) and F e,3 = 72.4kN (1.5 times design level), according to ACI 437.1R-07. 43I G U R E 1 Schematic of the loading scheme with half longitudinal section, midspan cross-section, and end section of specimen (dimensions in mm).
Each girder was simply supported by two squareshaped rubber bearings (150 Â 150 mm 2 in plan) to allow free rotations, while laterally restraining the girder to avoid torsional rotations.The vertical force F e was applied by means of a rigid steel beam with two loading points having a clear distance of 850 mm, thus resulting into a shear-span length (measured starting from the mid part of the support) L v = 2800 mm (Figure 1).
The girder deflection at midspan cross section δ e was measured using a vertical linear variable differential transducer while a load cell was installed between the vertical actuator and the rigid steel beam to measure the actual load imposed to the specimen.Further details on the experimental tests are reported in Losanno et al. 38

| Experimental results
The monotonic envelopes of the force-displacement response curves of the two specimens are shown in Figure 3, where the bending moment M e associated with each force level F e through shear-span length L v is also highlighted under varying deflection δ e .The total deflection (δ) can be obtained by combining δ e with prestressing-induced upward displacement (δ p < 0), the latter measured equal to À2.5 and À 0.5 mm under HP and LP levels, respectively.
During the experimental tests, the cracks pattern of each girder was also monitored by visual inspection.Figure 4 shows the front view of the midspan of the girders at two significant performance levels: (i) cracking onset F e,cr (Figure 4a,b) and (ii) vertical force at maximum deflection F e δ e, max ð Þ(Figure 4c,d).A first crack was observed around the midspan cross section of both girders, measuring the following forces and deflections: F HP e,cr ¼ 49:0 kN and δ HP e,cr ¼ 5:54 mm in case of S1-HP specimen; and F LP e,cr ¼ 17:4 kN and δ LP e,cr ¼ 1:74 mm in case of S2-LP specimen.It can be noted how the ratio between vertical forces at cracking initiation (F HP e,cr =F LP e,cr ≈ 2:82) is similar to the ratio between prestress levels (P HP =P LP ≈ 2:00), as expected in linear response range.Figure 3 shows the influence of the initial prestress level on the flexural response of each girder.Regardless of the prestress level, a similar initial vertical stiffness was found for the two specimens as well as both girders tend to attain the same peak vertical force although under different imposed displacements.Higher displacements were imposed to S2-LP due to additional belts pulling down the girder allowing the transverse head of the actuator rigidly translating downwards.The following values of initial stiffness (defined as the ratio between a vertical force of 10 kN and the corresponding imposed vertical displacement 38 3. Considering a vertical force level slightly higher than F e,cr for both specimens, S2-LP returns higher deformations than S1-HP.Accordingly, the initial prestress level significantly affects both the cracking threshold and the post-cracked response while playing a minor role on the uncracked phase and peak response.

| NUMERICAL SIMULATION
A multidimensional comparative numerical study is illustrated in the following sections.Numerical models of both PC girders were developed using two main approaches: (i) the FEM and (ii) the AEM.Using the first method, three types of FEs were used, namely, 1D, 2D, and 3D elements, whereas 3D AEs were adopted through the second numerical method.Accordingly, numerical analyses are herein referred as 1D, 2D, and 3D finite element analyses (FEAs), and 3D applied element analyses (AEAs).The 1D and 2D FEAs were run using the structural analysis software SAP2000©, 44 whereas 3D FEAs and AEAs using Abaqus FEA© 45 and Extreme Loading for Structures© (ELS), 46 respectively.

| Concrete
To perform 3D numerical analyses, concrete was modeled in both Abaqus and ELS according to fib Model Code 2010. 47The compressive and pre-cracking tensile response of the material are described in terms of stressstrain relationship as follows: where σ c is the compressive stress, σ ct is the tensile stress, 47 ), ε c,1 is the strain at maximum compressive stress, E c1 is the secant modulus from the origin to the peak compressive stress, f cm is the cylindrical compressive strength, ε ct is the tensile strain at cracking onset, and f ctm is the corresponding mean tensile strength.For the post-cracking response of concrete, a bilinear equation describes the tensile stress-crack opening (w) behavior as follows: being The numerical values of each of the parameters listed in Equations ( 1) and ( 2) are shown in Figure 5a,b, according to Section 2.1.
In 1D and 2D FEAs with SAP2000, the mechanical behavior of concrete was simulated through the Mander constitutive model. 48The compressive behavior is described by the following stress-strain equation: Regarding the tensile behavior, two linear stress-strain branches are used for the pre-cracking response (i.e., up to f ctm ¼ 0:62 ffiffiffiffiffiffiffi f cm p ) and post-cracking response (i.e., up to ε ct, lim ¼ αf ctm =E c , being α a constant parameter dependent on the ratio f cm =f ctm ).The complete model is illustrated in Figure 5c with the numerical notation for each of the parameters.
Concrete was modeled in Abaqus according to the concrete damage plasticity model (CPDM), which is a continuum, plasticity-based damage model. 45In CPDM, two main failure mechanisms are assumed together with default setting of parameters 45 : (i) tensile cracking and (ii) concrete crushing.The evolution of the yield (or failure) surface is controlled by two hardening variables, namely, the tensile equivalent plastic strain and the compressive equivalent plastic strain.Under uniaxial tension, a linear elastic stress-strain behavior is assumed up to the value of tensile strength (f ctm in Figure 5a), when micro-cracking in the concrete material arises.Beyond the failure stress, the formation of micro-cracks is represented macroscopically with a softening stressstrain or stress-crack opening relationship according to a smeared crack approach (Figure 5b).Under uniaxial compression, the behavior is linear up to yielding stress threshold (0:4 Â f cm in Figure 5a 47 ).The plastic response typically shows stress hardening up to compressive strength (f cm in Figure 5a, corresponding to ε c1 ), followed by a softening branch up to failure strain (ε c, lim in Figure 5a).Scalar degradation variables d c (in compression) or d t (in tension) are used to compute the reduction of the elastic modulus, as 5a for E ci ).These parameters are activated when the postpeak compressive/tensile response is triggered.

| Reinforcing steel
Prestressing and mild steel were considered as isotropic materials in each of the numerical simulations.Simplified numerical stress-strain models were used.As for the prestressing steel, three types of stress-strain models were used: (i) elastic-perfectly plastic model in 1D and 2D FEAs, based on the software parameters settings for the element adopted (see Section 3.2.1),(ii) elastic-plastic with hardening models in 3D AEAs, and (iii) Ramberg-Osgood model 49 in 3D FEAs.For the first type of constitutive model, the ultimate stress (i.e., f 1,2D py ¼ f pt ¼ 1969 MPa) was adopted as equivalent yielding, according to fib Model Code. 47These three assumptions for the constitutive laws of the prestressing steel were used to investigate the different results that can be obtained using simplified rather than accurate stressstrain models.
For each numerical analysis, mild steel was modeled using an elastic-plastic with hardening material model.The diagrams for both mild and prestressing steel are F I G U R E 5 Constitutive models adopted for numerical analyses: (a) compressive behavior of uncracked concrete, 47 (b) tensile softening behavior of concrete according to fib Model Code, 47 (c) Mander compressive and tensile, 48 (d) bilinear behavior of prestressing steel, and (e) mild steel.AEA, applied element analyses; FEA, finite element analysis.
shown in Figure 5d,e, which also include the values of the main mechanical parameters.

| One-and two-dimensional FE modeling
The development of the 1D FEMs of the PC girders relied on the use of two types of FEs, that is, the frame element for concrete components and the tendon element for prestressing steel (Figure 6a).Two-node linear elements were used to mesh both the frame and the tendon.In 2D FEMs, the strands were modeled using the same tendon object used in 1D analyses, while the concrete was modeled using 2D rectangular-shaped (B base Â B height ) four-node FEs.Both in 1D and 2D models, the tendons are embedded into the concrete frames/shells at each node (Figure 6a), reproducing the perfect bonding of the experimental tests.
Section designer from SAP2000 allows drawing the T-shaped cross section and the longitudinal rebars within.The parabolic-linear path of the prestressing reinforcement can be perfectly reproduced by the design menu of SAP2000 (Figure 6a).Cross-sectional area and material properties shown in Figure 5d are then assigned to both tendons.In both 1D and 2D models, the values of effective prestressing forces (Section 2.1) were imposed as internal force values.
The nonlinear response of the structure can be simulated via 1D FE modeling using specific plastic hinges (i.e., lumped plasticity) to be defined either in terms of bending moment-rotation of the element or bending moment-curvature of the generic cross section.Where no plastic hinges are defined, the numerical model develops a linear elastic response with no degradation of the flexural response regardless of the materials modeling.To define the number and position of the plastic hinges, the girder was divided into 13 frame elements, with variable lengths included in the range 548 mm (end frame elements) to 425 mm (midspan frame elements).The plastic hinges were then located at the midpoint of each frame element (Figure 7).Plastic hinges were defined starting from the analytical study of the bending moment-curvature response of the corresponding cross sections, that is, taking into account the eccentricity of the tendons at different positions along the girder.The following material models were used in the analytical study of the cross sections 39 : • Parabola-rectangular uniaxial stress-strain behavior for concrete in compression, with cylindrical compressive strength equal to 30 MPa and elastic modulus of 30,000 MPa.Conventional yielding strain and ultimate strain were assumed equal to 0.002 and 0.0035, respectively.• Linear uniaxial stress-strain behavior for uncracked concrete in tension, with tensile strength equal to 2.9 MPa.• No tensile resistance of cracked concrete.
• Elastic-plastic uniaxial stress-strain behavior with hardening for mild steel, with yielding stress and ultimate stress equal to 450 and 540 MPa, respectively.The elastic modulus was set to 210,000 MPa.Accordingly, the yielding and ultimate strains were respectively set to 0.0021 and 0.07.• Elastic-plastic stress-strain constitutive law for prestressing steel, using a yield stress of 1969 MPa. 47ssuming an elastic modulus of 204,215 MPa, yielding and ultimate strains were set equal to 0.0087 and 0.070, respectively.
The moment-curvature diagrams used for 1D FEAs are shown in Figure 8, where a very limited influence of the initial prestress level on the ultimate response of the girders is observed as per experimental observations.Four significant thresholds were defined according to software default settings, as follows: (B) first cracking, (C) first yielding of rebar, (D) first yielding of prestressing steel, and (E) ultimate compressive response of concrete.
Boundary conditions were assigned to specific joints of the girder to model both restraints and loads.Given the dimension of the analysis in terms of geometry and loading conditions, the 1D model of the girder was restrained by assigning a pin to the left-hand side end joint and a roller to the other end joint.Two joints located at L v from the end sections were also implemented to assign the vertical loads, while monitoring the output force and displacement at the midspan cross section (Figure 6a).
In the 2D models developed in SAP2000, a four-node multilayered nonlinear shell element was used to model the concrete girder (Figure 6b).Two different layers were assigned as concrete and longitudinal rebars.The material properties given in Figure 5c allows the nonlinear response of the single concrete shell element to be taken into account, whereas the layers of rebars are modeled using the elastic-plastic with hardening models of Figure 5e.Given the T-shaped cross section, the flange and web were separately modeled, meshed, and merged with each other.Prestressing strands were modeled using the same tendon object and material properties used in 1D analyses (Figure 6b).Boundary conditions were defined by assuming the web shell to be simply supported at the two bottom joints of the end sections, hence restraining them against both vertical and horizontal translations as well as out-of-plane (i.e., torsional) Moment-curvature diagrams of plastic hinges used in 1D finite element analyses.S1-HP, Specimen 1, high prestress; S2-LP, Specimen 2, low prestress.rotations.The vertical loads were applied to two joints of the flange, allowing nonlinear response assessment of the girder as load versus midspan cross section deflection (Figure 6b).

| Three-dimensional modeling according to FEM and AEM
In 3D analyses, brick element with cubic shape is used to model the concrete girder.The 3D FE modeling of the girders was performed through the following elements: eight-node first-order elements with full integration rule (named C3D8 in Abaqus library) for concrete, and twonode truss elements (named T3D2 in Abaqus library) for both mild and prestressing steels (Figure 6c).The embedded technique was used to model the interaction between longitudinal bars, stirrups, and tendons, assuming a perfect bonding condition with the surrounding concrete.
The initial prestress was applied to the tendons using a temperature-based approach. 50,51For 1D rod element, the effective tendon initial prestress σ p,e (Section 2.1) can be obtained using the following equation: being α p the thermal expansion coefficient of prestressing steel and ΔT the applied thermal variation corresponding to prestress level at the day of testing.
The applied element modeling of the girders was carried out by discretizing the structural element into a number of discrete elements, which were connected to each other through a predefined number of normal and shear nonlinear springs assigned to the elements' surfaces.The strain level of each spring defines the damage reached in a single part of the element.The behavior of the defined springs depends on the relative distance on the same surface of the element; moreover, stresses and deformations of each spring are related to the volume of the elements connected by the spring.The stiffness matrix of each spring is defined starting from the element centroid displacements.Then, the global stiffness matrix is obtained by combining different stiffness contributions.
The AEM approach allows accounting for the failure of any springs based on stress-or strain-criteria.When the stresses and/or strains of the single spring reach their cracking limits, cracks are formed (i.e., the so-called "open cracks" in ELS).Further or wider cracks develop under increasing load, while they may close under cyclic loading (i.e., loading/unloading pattern).When maximum strain attains the ultimate deformation capacity (i.e., the so-called "separation strain" in ELS) the springs become inactive, and elements are disconnected from each other.Thus, a permanent crack is formed, and adjacent elements are separated.
Reinforcing bars and strands were modeled as specific springs.The "custom RFT" tool of ELS was used to generate user-defined reinforcement layout, including customization of prestressing tendons and mild steel (Figure 6d).In ELS, the prestress was applied as force along the strands.
In 3D modeling, the boundary conditions were applied to the girder through suitable model items, namely, bearings and loading crossbeam/pins (Figure 6c,d).Given their very limited influence on test results, these parts were modeled via linear-elastic materials with elastic stiffness set equal to 2000 and 210,000 MPa, respectively.In Abaqus, a contact interaction between two items (i.e., between girder and supports, or girder and crossbeam) was used, assuming a friction coefficient of 0.6 for simulation of shear behavior over the interaction surface. 52In ELS software, the interaction between two items was automatically managed, controlling the contact with respect to the mechanical properties of the two contact bodies.

| Nonlinear analysis procedures
Nonlinear static analysis is implemented in each of the numerical method through two loading steps: (i) combined self-weight and prestressing and (ii) four-point bending test according to the experimental set-up (Section 2.2).The force-displacement response curves are then obtained for the midspan cross section.
The mixed iterative event-to-event nonlinear solution control method has been used in SAP2000. 44Using a specified iteration convergence tolerance equal to 0.0001 (set as default), iteration is used to reach equilibrium at each step and substep of the analysis.In each step, first a constant-stiffness iteration is tried; if no convergence is obtained with this attempt, the Newton-Raphson (tangent-stiffness) iteration is then tried.When both of these methods fail, the step size is reduced, repeating the process.With the event-to-event strategy, an increment of load is applied until significant change in the stiffness of the structure occurs (i.e., in one of the elements a nonlinear "event" occurs).The stiffness matrix is then reformed, then to apply a further increment of load until the next "event".
In 3D FEAs, the Abaqus/standard nonlinear implicit solver is used. 53Using implicit method, equilibrium is enforced between the external applied loads and the internal generated reaction forces at every solution step using Newton-Raphson method.This method is both incremental and iterative and ensures an unconditional stability.
Two types of solver options are available in 3D AEAs software 54 : (i) the exact solver and (ii) the iterative solver.The first type of solver can be used both in static and dynamic analysis when the estimated time stepping is relatively large (i.e., >0.001 s), while the iterative solver needs to be used in complex nonlinear cases (e.g., penetration or blast analyses) where the time stepping is very little (i.e., <0.0001 s).For the purposes of this work, the exact solver was used.

| Mesh sensitivity analysis
For each numerical model, the optimal size of the single mesh element of the concrete girder was studied.Two outputs were used for this scope, namely the maximum (F e, max ) and ultimate (F e,ult ) vertical force, corresponding to peak and post-peak flexural response of the girders, respectively.More specifically, the ultimate force was associated with either numerical instability or ultimate stress/strain of materials.In order to have a measure of the post-peak branch, the mesh size with the maximum ratio F max =F ult was selected as the optimum.
Given the frame partition according to the locations of the plastic hinges (Section 3.2.1),no sensitivity analysis was performed on the 1D models as the size of the single linear finite element within the frame object does not affect the results under proper calibration of plastic-hinge parameters.Thus, the SAP2000 default auto-mesh tool was used for 1D FEA.
For 2D FEAs, the aspect ratio (i.e., ratio between base and width) of the single rectangular element was assumed variable between 0.55 and 8.7.In this range, a maximum percentage variation of the control ratio F max =F ult equal to 0.21% was obtained, resulting in a slight influence of the mesh size on the results of interest in this paper.Thus, for the web 3 elements along the height and 96 along the width (i.e., aspect ratio equal to 0.543) were set, while for the flange 4 elements along the width and 96 along the length of the girder (i.e., rectangular element with base to width ratio equal to 1.7) were chosen.
In 3D AEAs, base sides in the range 80-120 mm were tested.It was seen how the post-peak and collapse response can be obtained regardless of the size of the single discrete element, with a minimum F max =F ult equal to 3.74.Thus, 100 Â 100 Â 100 (mm Â mm) brick element was chosen both for S1-HP and S2-LP.
3D FEAs results show more significant meshdependency, as can be seen from Figure 9.The maximum F max =F ult ratios were obtained using 50 Â 50 Â 50 (mm Â mm) brick elements for both S1-HP and S2-LP, while increasing the mesh size would return unstable analyses after peaking the vertical force.
It is worth pointing out how the sensitivity mesh analyses returned significantly higher values of F max =F ult through 3D AEAs (minimum ultimate-to-maximum force ratio equal to 3.74) compared to 3D FEAs ( 4 | DISCUSSION OF RESULTS

| Force-displacement curves, damage distributions, and model error estimation
The experimental-numerical comparison regarding S1-HP and S2-LP specimens are shown in Figure 10 in terms of force-displacement curves, assuming the total deflection as δ ¼ δ p þ δ e (see Section 2.3).Similar to Figure 10, the maximum bending moment associated with force/deflection level is also given along the secondary vertical axis.The numerical results include the postpeak branch prior to numerical instability (i.e., softening of the concrete).Thus, Figure 10 also reports a prediction of the response of the girders beyond the maximum experimental displacement, where the test was stopped due to stroke limitations of the actuators.An excellent experimental-numerical agreement was found through all selected numerical strategies, with 1D FEA only limited to peak response with no softening according to plastic hinge parameters.The response provided by each of the numerical models fully reproduces the experimental curves, ranging from the initially linear response to the ultimate condition.Therefore, the nonlinear behavior of the specimens was well simulated accounting for both concrete cracking and reinforcement yielding.
The initial upward deflection due to prestressing action (i.e., negative deflection at zero external load) is well matched by the numerical models.Table 1 outlines the numerical deflections compared to the experimental values, showing maximum absolute differences equal to 0.43 and 0.12 mm among the different strategies for HP and LP levels, respectively.
The first cracking was obtained from numerical models corresponding to the value of vertical force where the tensile strength is attained at midspan.Values of F e,cr are reported in Table 1, that shows a maximum percentage difference in the range [0%, 10%] among the different strategies.A better accuracy is obtained considering the peak load of the girders, observing a maximum percentage difference of 2%.
A comparison between numerical and analytical predictions corresponding to moment-curvature response of midspan cross section in Figure 8 is also reported in Table 2.An excellent matching is obtained for each of the selected performance indicator, with maximum absolute differences equal to 0.34 mm (HP) and 0.19 mm (LP) on the prestressed-induced upward deflection, and maximum percentage differences in the range [À5%; +8%] for the cracking and peak vertical forces.Such an outcome allows validating numerical models against well-known analytical models and demonstrating a satisfactory accuracy of hand calculations.
Figures 11-14 report the response of the girders at peak vertical force through the different models.In 1D FEAs, the incremental vertical displacement gradually activates the plastic hinges at different locations of the girder.In the case of S2-LP specimen, the two plastic hinges located around the midspan cross section attained the maximum bending moment (Figure 11b), whereas ultimate nonlinear response of S1-HP specimen propagated outside the loaded part (Figure 11a), as the analytical ultimate curvatures related to the ultimate bending moment of HP girder are smaller compared to the LP case.
The 2D FEAs provide the axial stresses across each shell (Figure 12), giving information on the tensile and compressive response of every section.Throughout the   In Abaqus, contour plots in terms of concrete tensile damage (Figure 13) provide an excellent simulation of the cracking pattern of the elements.Both for HP (Figure 13a) and LP (Figure 13b), the 3D FEAs were able to return location and extension along the height of the girder of the single cracks that develop at different vertical load levels.
Finally, contours plots in terms of axial tensile strains in 3D AEAs showed larger concentration of the maximum values at the midspan of S2-LP specimen (Figure 14b) than S1-HP specimen (Figure 14a).Even if both specimens attained similar peak loads, the different modeling strategies confirmed that the initial prestress level affects the local response of the girder in terms of damage distributions.
The global accuracy of each numerical model is herein computed using the following error coefficient: where F num,i and F exp,i are the numerical and experimental forces at the i-th level of total deflection δ i .The model error (E n ) for each numerical strategy and specimen was thus statistically characterized as normally distributed random variable with the following mean value and standard deviation: where n is the number of deflection levels used to discretize the overall displacement range.The accuracy of a numerical analysis increases when E num reduces, that is, the numerical curve matches the experimental one with better approximation.Figure 15 shows the 5 mm discretization intervals of the experimental curves while Table 3 outlines the mean value and standard deviation for each numerical model and girder specimen (or, equivalently, prestress level).As can be seen, the maximum percentage error is equal to 2.22% and is related to 1D FEAs of S2-LP with a standard deviation equal to 0.006; the best accuracy is obtained for 3D FEAs of S2-LP, with an error of around 0.3% and a standard deviation equal to 0.006.

| Stress distributions within tendons
The tendons internal stresses (primary axis) and corresponding internal forces (secondary axis) at midspan cross section were plotted against the imposed vertical displacement in Figure 16, addressing the two different steps of the load pattern during the analysis.In the first step, the prestress was applied up to the measured value at the day of testing (i.e., P e,HP or P e,LP , see the slope to the left-hand side of Figure 16; see Section 2.1), and these values were compared in Figure 16 against the target ones (i.e., jacking prestress, dashed lines) both for HP and LP.Loading-induced increasing stresses were then obtained in both tendons following their stress-strain constitutive model (slope to the right-hand side in Figure 16).The initial prestress level affects the yielding point of the tendons as a lower imposed displacement is required to yield the strands in case of S1-HP specimen compared to S2-LP specimen.In no 3D, cases the ultimate tensile strength was reached at ultimate conditions (f py ≤ σ p < f pu ) in the range of the deflections investigated.These results show how using accurate stress-strain models (i.e., Ramberg-Osgood or bilinear with hardening) rather than simplified constitutive laws (elastic-plastic), provides measurements of the stresses trend in the strands in line with the experimental results, as no tendons' rupture was attained during the experimental tests as confirmed by peak strains lower than ultimate values.
Table 4 reports on the numerical values of the displacement at the yielding of tendons (δ yp ) and maximum stress obtained in the analyses (σ p δ max ð Þ).Tendon yielding at midspan was achieved under a vertical displacement ranging in the intervals [58, 77 mm] and [79, 131 mm] in the case of S1-HP and S2-LP specimens, respectively.

| Comparative assessment and discussion
The results of the numerical simulations demonstrated how the flexural response of PC girders can be predicted with different numerical strategies resulting into different modeling effort, computational burden, and response accuracy.A number of features were considered for an exhaustive comparison between different analyses, with pros and cons that are herein discussed and summarized in Table 5.
i. Computational burden.As known, the computation time of the numerical models depends on different factors, such as the mesh size of the individual elements, the constitutive models of the materials, and boundary conditions.The analyses were performed using an Intel(R) Core(TM) i7-4790 CPU, with 3.60 GHz clock frequency and 8.00 GB RAM.
With the given material properties, boundary conditions and mesh size, for FEAs less than 1, 5, and 60 min are needed to completely solve 1D, 2D, and 3D FEs models, while surprisingly ELS completely run the job in 5 min only.ii.Nonlinear global response.In 1D analyses, the nonlinear response of the girder can be obtained solely starting from a proper definition of plastic hinges properties and their discrete locations.This implies the analytical study of different sections of the PC girder and prior knowledge of the global response of the same girder for the correct positioning of each plastic hinge.Considering 2D and 3D analyses, the nonlinear response can be directly implemented using appropriate nonlinear elements combined with nonlinear material behavior.iii.Local response.Both 1D and 2D FEAs allow a satisfactory simulation of the global response of PC girders with lack of information in terms of local response, such as onset and progressive cracking, steel yielding, and tendon-concrete interaction.Such an information can be obtained from 3D analyses with plenty of possible output parameters, such as local stress-strain response of concrete and steel reinforcing elements.iv.Damage.The simulation of damage to concrete elements can be carried out using 3D FEAs.The CDPM used in Abaqus provides detailed response in terms of incremental concrete damage, both in compression and tension.Damage is implicitly considered in ELS through the failure of connection springs between adjacent elements, but explicit information on the concrete progressive damage can be solely estimated in terms of axial strains.v. Peak response.As observed in Section 4, the ultimate response (including peak flexural force/ bending moment and vertical deformations) can be predicted with good accuracy using each of the numerical strategies implemented in this study.Peak response was found to mainly depend on materials' strength, particularly tendons yielding and ultimate stress values (Figure 16).vi.Softening.The post-peak flexural response of the PC girders could be predicted by 2D and 3D analyses.In 1D models, the post-peak response would  demand a softening phase when modeling plastic hinges, but this would limit the ease of modeling.In 2D analyses, additional softening in the forcedisplacement curve can be obtained under larger deflections with further computational burden.vii.Progressive collapse.The collapse behavior of the girder can be studied in 3D analyses.This can be addressed as the most promising feature of AEA implementing progressive collapse condition (i.e., element separation) both for static and dynamic analyses.Using the implicit solver in 3D FEAs as done in this work would limit the post-peak response of the models, and the element deletion can only be implemented using an explicit dynamic solver.However, in this type of analysis the computational burden would increase significantly.viii.Interaction between model items.The 3D FEA in Abaqus allows for an accurate modeling of the T A B L E 5 Pros and cons of the numerical strategies.
• Ease of modeling (including sections with complex geometry).• Large number of analyses (i.e., results) in short time.
• Accurate prediction of the global response up to peak flexural force.
• Approximate and average local response parameters • Only in-plane response can be studied.
• No interaction between independent objects.
• Need for analytical study of different sections for plastic hinge setting.
• Ease of modeling.
• Accurate results of the global response up to peak flexural force.• Distributed plasticity.
• Softening included under large displacements.
• No interaction between independent objects • Tricky to model sections with complex geometry.

3D FEA
• Accurate modeling of each object.
• Damage modeling of the concrete elements.
• Accurate prediction of the peak and post-peak flexural response.
• Mesh dependent when using implicit solver.
• Collapse response dependent on the numerical solver.
• Accurate prediction of both peak and post-peak responses.• Progressive collapse implemented in both static and dynamic analyses.• Very little mesh dependent.
• No manageable interaction between different independent objects.• Rebar modeling according to default settings of the software.
interaction between two or more model items; thus, the boundary conditions (e.g., supports and loading crossbeam, adhesion of concrete with reinforcing steel and tendons) can be modeled with great accuracy, setting the appropriate interface behavior.In 3D AEA, ELS automatically adopts a perfect bond between items.
In light of the above considerations, some useful suggestions can be drawn for practical applications.The best features of 1D and 2D modeling techniques are the minimum computational burden and the fast running time.As these models are able to reproduce the global response of a PC member with good-to-excellent agreement, they can be used to obtain a prediction of the global flexural capacity of beam-type bridge girders.The 2D model can be adopted when plastic hinge properties are not available.
When considering 3D FEA, the increase in modeling accuracy produces a rising computational burden.This technique is recommended when local response is of greatest interest: to the best of the author's knowledge, this is the only modeling option to account for local defects in PC structures (e.g., unbonded and partially bonded tendons, deteriorated steel-concrete bond properties) or complex interaction mechanisms between solids.
Finally, 3D AEA provides excellent global results with a rather low computational effort compared to other equidimensional numerical strategies.The best feature of this technique can be found in the post-peak and progressive collapse analysis, thus it can also be used for large scale PC facilities.

| CONCLUSIONS
This paper has investigated the flexural response of two 1/5-scaled posttensioned PC girders, comparing experimental results to the simulation outputs of different computational strategies.The specimens were designed with two different levels of initial prestress (high and low), then four-point bending tests were conducted up to peak flexural response of both girders.On the basis of the monotonic envelope of the flexural response, a comprehensive comparative multidimensional numerical study was developed.Each PC girder was modeled using different analysis dimensions and numerical strategies: 1D and 2D finite element analyses using SAP2000, 3D finite element analyses using Abaqus, and 3D applied element analyses using ELS.Taking into account the different types of numerical models of PC girders, this study provides an unprecedented comparative assessment between distinct numerical simulations.
At both global and local scales, simulation results have been presented in terms of force-displacement curves, vertical displacements due to initial prestressing, first cracking, and stress pattern in the tendons against the vertical displacement imposed to the specimens.Even if the cracking and peak load were approximated with good accuracy through all the methods, local phenomena such as cracking and damage pattern as well as reinforcement-concrete bond could only be detected at the scale of 3D models.While 3D FEA requires a very high numerical burden, 3D AEA also resulted in stable softening behavior beyond peak response with a very limited computational effort.The numerical errors were evaluated for each numerical strategies as a measure of the accuracy of the numerical predictions.It was seen how each of the numerical methods matched the experimental force-displacement curves with great accuracy, as a maximum percentage difference equal to 2.22% (corresponding to 1D FEAs of the lowly prestressed girder) was obtained.Also, the values of the standard deviations computed on this parameter are lower than 0.06.
In the last part of the paper, the main features of each numerical strategy have been addressed to highlight the main pros and cons of the different analysis methods.The 1D FEAs have been issued as the easiest and fastest types of modeling, however at the expense of the accuracy at large vertical displacements.To model the post-peak response, at least 2D elements are needed.This numerical strategy provides quite accurate results, but the modeling of the geometry could be elaborate.Both 3D software are equipped with the best features, as they allow for an accurate prediction of the whole flexural response of the girders, starting from initial behavior to the post-peak and collapse responses.In 3D FEAs, the numerical models can be customized, managing the contact interaction between different items and implementing user-defined material models; however, the largest computational burden is associated to this software.The 3D AEAs returned the most accurate predictions of the girders' responses with a rather low computational burden.However, some setting could not be managed in a direct way and the local response is not as accurate as 3D FEAs.
Further investigations will consider the multidimensional modeling of special issues related to prestressed structures, such as grouting defects, damage, and corrosion of prestressing steel.The model error that was characterized through this study might be considered in future probabilistic investigations, allowing the accuracy of each numerical model to be taken into account when assessing the structural performance and safety of posttensioned PC girders.In addition, the computational strategies herein investigated for longitudinal girders might be used to evaluate the nonlinear behavior of bridge deck up to collapse.

F
I G U R E 2 (a) Experimental four-point bending test setup (front view) and (b) specimen before testing with view of the LVDT, linear variable differential transducer.

F I G U R E 3
Experimental envelope curves of specimens: force and bending moment versus deflection with identification of cracking force levels.S1-HP, Specimen 1, high prestress; S2-LP, Specimen 2, low prestress.F I G U R E 4 Midspan crack patterns of Specimen 1, high prestress (S1-HP) specimen (a, c) and Specimen 2, low prestress (S2-LP) specimen (b, d) at cracking initiation and peak load.

F
I G U R E 6 Numerical models: (a, b) 1D and 2D FEM with SAP2000, (c) 3D FEM with Abaqus, and (d) 3D AEM with Extreme Loading for Structures.AEM, applied element method; FEM, finite element method.

7
Position and nomenclature of the plastic hinges in 1D numerical models.

F I G U R E 1 1
Plastic hinges status at peak vertical forces in 1D finite element analyses.S1-HP, Specimen 1, high prestress; S2-LP, Specimen 2, low prestress.central loaded part of the girder, tensile stresses lower than the tensile strength are associated with the softening branch of Figure 5c.

F I G U R E 1 4
Concrete axial tensile strains in 3D applied elements models at peak vertical force.S1-HP, Specimen 1, high prestress; S2-LP, Specimen 2, low prestress.FI G U R E 1 5 Discretization of the experimental curves for the model errors evaluation.S1-HP, Specimen 1, high prestress; S2-LP, Specimen 2, low prestress.
Experimental-numerical comparison in terms of prestressing-induced deflection, cracking force, and force at maximum experimental deflection.