Shear Lag Analysis of Simply Supported Box Girders Considering Axial Equilibrium and Shear Deformation

: The conventional methods used to analyze the shear lag eﬀect in simply supported box girders assume that the neutral axis of the section coincides with the centroid, which does not strictly satisfy the axial equilibrium condition. To address this problem, this study proposes an analysis method in which three independent functions for the shear lag are employed to deﬁne the diﬀerent shear lag strengths of the top slab, the bo�om slab, and the cantilever slab. To fulﬁll the axial equilibrium condition of the box girder and to automatically locate the neutral axis position, the longitudinal displacement of the web is introduced. The shear deformation of the box girder is also considered. The governing diﬀerential equations and corresponding boundary conditions for displacement variables such as deﬂection and rotation of the box girder are derived through the application of the principle of virtual work. The diﬀerential equations are solved by utilizing the boundary conditions to obtain the analytical expressions of the shear lag function, longitudinal displacement of the web, rotation, deﬂection, and neutral axis position. Furthermore, after performing the ﬁnite element analysis, the eﬀectiveness of the proposed method is veriﬁed by comparing the results with those obtained from conventional methods and ﬁnite element analysis. Furthermore, the inﬂuence of the axial equilibrium condition is quantiﬁed on axial force and stress diﬀerence ratios under three methods. Finally, extensive parametric analysis is carried out to investigate the eﬀect of diﬀerent parameter ratios on the ratios of the stress diﬀerence of the ﬂanges. The results show that when the axial equilibrium condition is not considered, the axial stresses in the upper ﬂange of the simply supported box girder are underestimated, especially at the intersection of the top, cantilever slab, and web, and the axial stresses in the lower flange are overestimated. As a result, the method in this study is able to calculate the axial stresses and deflections on simply supported box girders more accurately.


Introduction
The uneven distribution of bending stress along any cross section of a thin-walled box girder and the non-uniform longitudinal stress distribution in the flange at the webflange connection due to shear interaction result in what is known as the shear lag effect [1,2].The shear lag effect is widely observed in various engineering structures, for example, beam structures, bridge engineering, and concrete structures [3][4][5][6][7][8].Neglecting to consider the shear lag effect can lead to two detrimental consequences: an underestimation of vertical deformation in the box girder and an underestimation of peak stress in the wing plate, significantly impacting the applicability and safety of the structure.
In previous studies, the shear lag effect of box girder had been extensively studied.Reissner [9] introduced the orthotropic plate method for the analysis of the shear lag effect in cantilever box girders with corrugated wing plates.Zhang [10] proposed a simpler improvement of the segment method to analyze the shear lag in cantilever box beams, and its correctness was verified through finite element analysis.Abdel-Sayed [11] further applied the folding plate orthotropic plate method for effective analysis of steel bridge panels.The result showed that the sectional area of the longitudinal bracing in relation to the sectional area of the skin plate is also significant.However, these research methods are limited to simple box-girder structures and fail to address the requirements of complex engineering scenarios.Reissner [12], in a pioneering effort, employed the energy variational method to investigate the shear lag effect in rectangular box girders.As a result, the shear lag effect was found to be the reason why the deflection of beams increased.Subsequently, according to Reissner's method, Zhou [13] established a finite element method that took into account the interaction between the bending deformation and the shear lag deformation of a box girder and defined a shear-lag stiffness matrix.Moreover, Wu and Liu et al. [14] established the associated finite element stiffness matrix and equivalent nodal forces vector.Moffa and Dowling [15] adopted the finite element method to examine the impact of the shear lag effect on steel box girders.Lee et al. [16] investigated the root causes of the negative shear lag effect in cantilever box girders using plate and shell elements.Prokic [17] proposed a finite element method based on the warpage function for shear lag effect analysis of box girders.Foutch [18] and Chang and Zheng [19] analyzed the negative shear lag effect using the energy method.Cheung [20] introduced the finite strip method and applied it to the study of the shear lag effect in separated and multicompartment box girders.Cheung gave a detailed description of the application of the finite strip method in practical engineering.Chang et al. [21] utilized the finite difference method to analyze the shear lag effect in variable-section cantilever box girders.Furthermore, Yang et al. [22] introduced the random wired segment method for the analysis of the shear lag effect in box girders.However, these studies assume that the neutral axis coincides with the centroid axis of the cross section.Consequently, the axial equilibrium condition of the box girder section cannot be strictly satisfied, resulting in significant errors in the analysis results.
For the purpose of fulfilling the axial equilibrium condition, Ni et al. [23] constructed a deflection function that achieves axial self-equilibrium of normal stresses in the crosssection and analyzed the shear lag effect in wide-flange beams.Xie et al. [24] assumed a new longitudinal displacement function and derived the differential equations for the bending deformation of thin-walled box beams, as well as the formulation for calculating the element stiffness coefficients.Based on the energy variational method, Chen and Shen et al. [25] proposed a closed-form solution for the shear lag function, and the results indicate that this method can accurately predict the normal stress in the flange.Subsequently, Lin et al. [26] used an infinite series of higher-order polynomials to describe the longitudinal displacement of the flanges of the box beam and introduced additional terms in the longitudinal displacement to consider the axial equilibrium conditions.However, the significance of the additional items in these studies is difficult to understand, and they ignore the shear deformation of box girders.Zhou et al. [27] derived the fundamental expression of the shear lag deflection function by considering the relationship between the in-plane shear deformation and the longitudinal displacement of each flange plate in the box girder and proved that there is a discrepancy between the neutral axis and the centroid axis throughout the cross-section.Therefore, Li and Wan [28] introduced a novel and precise theoretical analysis approach to forecast the shear lag phenomenon in thin-walled singlebox girders, which demonstrates remarkable effectiveness and accuracy.They made sufficient adjustments to the coefficients at the flanges and compared the results with the finite element method.We propose a new formulation of the beam finite element method (B3S) for predicting the shear lag and shear deformation effects of thin-walled singlechamber and multi-chamber box girders [29].However, the additional terms are very complex for practical applications in the existing research methods.Meanwhile, the axial force is not considered, and the calculation process is complicated.For box girders with different dimensions, further comprehensive studies should be conducted.Therefore, Tan and Zhang et al. [30] studied the impact of the shear lag effect on ultra-high-performance concrete structures, and a parametric investigation of various factors has been conducted to further examine the impact of different parameters on the shear lag effect of the structures.In addition, Zhao and Liu et al. [31] proposed a beam finite element model considering the shear lag effect of steel-concrete composite box beams and found that the warping intensity function of the concrete slab and the steel beam at the end induced by shear lag exhibited a respective increase of 111.64% and 7.01%, respectively.
This study proposes an analytical method to overcome the limitations of conventional methods, which consider both axial equilibrium and shear deformation.During the analysis and derivation process, three independent functions for shear lag are employed to characterize the varying intensities of shear lag in each slab.Additionally, the axial equilibrium condition of the box girder is satisfied by incorporating the longitudinal displacement of the web, which automatically determines the position of the neutral axis.Using the principle of virtual work, obtain the corresponding control differential equations and boundary conditions.The governing differential equations were solved by utilizing the boundary conditions to obtain the analytical expressions of the shear lag function, longitudinal displacement of the web, rotation, deflection, and neutral axis position.Furthermore, the finite element analysis was performed, comparing the results obtained from the finite element method, conventional methods, and the proposed method.Moreover, the influence of the axial equilibrium condition on the calculation results of a simply supported girder is quantified by calculating the stress, axial force, and stress ratio under three different methods.Finally, investigate the influence of several different parameters on the stress difference ratio by changing the model parameters and give some reasonable suggestions.

Box Girder Analysis Model
Figure 1 shows the cross-section of the RC box girder, where

Strain and Displacement Analysis of Box Girders
To develop the shear lag analysis model for box girders, the following four assumptions have been made: (1) Uniformly distributed loads are symmetrically applied on the top surface of the web to prevent cross-sectional torsion, distortion, and lateral bending.

Reinforcing
(2) The shear deformation in the box girder is analyzed using the Timoshenko beam theory.
(3) The vertical and transverse strain of flanges and out-of-plane shear strain are negligible.(4) The reinforcement and concrete remain in a linear elastic state, with no occurrence of bond slip between them.
Assumptions play a crucial role in any analysis as they simplify the problem.However, it is important to acknowledge that assumptions may introduce some level of approximation and may affect the accuracy of the analysis.Especially for the assumption of no. 4, in the actual situation, the bond slip between the steel bar and the concrete is difficult to avoid.The bond slip phenomenon causes the relative displacement between reinforcement and concrete, which affects the stiffness, deformation, and bearing capacity of the structure.Therefore, if the adhesive slip is ignored in the analysis, it may lead to an inaccurate prediction of the structural response.In some cases, this can lead to underestimating the deformation and load-bearing capacity of the structure.
Based on the above assumptions, the displacement function of each slab in the box girder can be established according to Figure 2, as follows: The meaning of 1, 2,3, 4 refers to each slab of box girder; w is the longitudinal dis- placement of the web;  is the angle of the section; i f and i  are the shear lag func- tions and shear lag displacement functions, respectively; and v is the deflection.
Based on Equations ( 1)-( 3), the strain in each slab of the box girder can be expressed as follows: where i  and i  are the strain in each slab of the box girder.
The selection of the shear lag displacement function is crucial in shear lag analysis.Therefore, the following principles should be followed: (1) For thin-walled structures, the shear lag displacement function is a continuous function solely dependent on variable x and symmetric around the Y-axis [32].(2) ( ) 0 i x    ensures that the shear stresses are zero at the center of the top slab, the center of the bo om slab, and the edges of the cantilever slab.
(3) ( ) 0 i x   ensures the continuity of deformation at the junction of the flange and the web.
In the analysis of shear lag in box girders, the use of a quadratic parabola as the warpage displacement function has proven adequate [33].Therefore, i  can be expressed in the following form:

Governing Differential Equations and Boundary Conditions
According to Equations ( 1)-( 10), following Xiang and He [34], and applying the principle of virtual work, the global equilibrium condition can be expressed in the following form: The virtual work carried out by the internal forces of the box girder is shown on the left-hand side of the equation, where N, V, and M are respectively the axial force, shear force, and bending moment borne by the entire cross-section of the box girder, and Wic and Qic are respectively the double moment and shear force borne by each slab of the box girder.The virtual work carried out by the external forces of the box girder is shown on the right-hand side of the equation, where , , , z y i p p q g are the external loads acting along the longitudinal direction of the box girder, and N * , V * , M * , and Wi * are the external loads acting on the beam end.The internal forces borne by the cross-section of the box girder, included in the lefthand side of Equation (11), can be defined as follows: In the above formula, Equation ( 14) is multiplied by y on the basis of Equation ( 12), and Equation ( 15) is multiplied by i  based on Equation (12), where EA and EB are the tensile section stiffnesses that are independent of the shear lag displacement function, EI is the flexural section stiffness that is independent of the shear lag displacement function, subscript  relates to the parameters of the shear lag displacement function, , ,  are the cross-sectional stiffness of each slab of the box girder, and sh GA is the shear stiffness.
sh GA and id G  can be expressed as follows: where c G is the shear modulus; sh A is the shear area of the section; and ic A is the con- crete integration area.
By performing divisional integration on Equation (11), the following control differential equations are obtained: The boundary conditions at the beam end can be expressed as follows: Equation ( 19) represents the longitudinal balance of the box girder.Equation (20) represents the vertical balance of the box girder.Equation (21) represents the balance of the box girder in the rotational direction around the X-axis.Equation (22) represents the balance between shear and axial stresses in each slab of box girder.

Analytical Solution of the Governing Differential Equations
Adding the subscripts h and p to represent, respectively, the homogeneous and particular solutions of the generalized displacement equation yields the following:

Homogeneous Solution
According to the homogeneous form of Equations ( 19)-( 22), the second-order differential equations of h w and h  with respect to ih f are obtained as follows: , 1, 2,3 where By substituting Equations ( 28) and ( 29) into Equation ( 22), the homogeneous differential equations can be expressed in matrix form as follows: where The elements of the matrix above are defined as follows: , 1, 2, 3 , , 1, 2, 3 Based on Equation (32), the homogeneous solution of ih f can be wri en as follows: where

A Particular Solution
Assuming that the box girder is subjected to a uniformly distributed vertical load y p , the specific solution of the shear lag function is: where d is an integral constant that can be determined by the boundary conditions as follows: Based on Equation ( 2), the coordinates of the neutral axis y can be wri en as:

Validation of the Method
The finite element method (FEM) can provide abundant data for comparison with the proposed method (PM).By establishing a finite element model and comparing the results obtained from finite elements with the results obtained from the proposed method and the conventional method, the curve of the results obtained from the proposed method is closer to the curve of the finite element results, so the results of the proposed method are closer to the finite element results, which verifies the validity of the proposed method.In this paper, according to the derivation process of the proposed method, when the variable w is eliminated, the resulting analytical solution can be approximated as the ana- lytical solution of the conventional method.The model of a simple box beam was referenced in the literature [35].Consequently, the data analysis is conducted using the model shown in Figure 3. Figure 3a  Finite element analysis is performed using ABAQUS 2021.The concrete element is simulated by using C3D8 solid elements (Figure 4a), while the steel reinforcement is simulated by using T3D2 beam elements (Figure 4b), with relevant properties assigned.To avoid stress concentration at the supports, a suitable thickness of pads is added at both ends of the supports.Reference points are established to apply the boundary conditions to simply supported conditions.The left support is modeled as a hinged support, while the right support is modeled as a link support.The embedded region functionality is used to embed the reinforcement within the concrete, and the tie functionality is used to connect the pads with the concrete.The reference points are coupled with the pads.Figure 5 shows the finite element operation results.Figure 6 shows the mid-span axial stress curves of each slab.Obviously, the results obtained from the proposed method are closer to those obtained from the finite element method than those obtained from the conventional method.Figure 7 shows the deflection curve.The deflection curve obtained from the proposed method in this study closely aligns with the deflection curve resulting from the conventional method, with a negligible difference compared to the finite element deflection curve.By comparing the stresses obtained from the above analysis, the proposed method provides stress values that are closer to the finite element results compared to the conventional method.This validates the effectiveness of the proposed shear lag analysis method that considers the variation in neutral axis height.To facilitate a comprehensive analysis of such an influence, the differences in axial stress and deflection are defined as follows: where p s is the stress of the proposed method, and c s is the stress of the conventional method.
where p v is the deflection of the proposed method, and c v is the deflection of the con- ventional method.
Table 1 shows the ratios of the stress difference and the ratios of the deflection difference in the simply supported box girder.It is evident that the change in neutral axis height has li le influence on the deflection but significantly impacts the stress.Therefore, the subsequent sections focus on analyzing the influence of the change in neutral axis height on the axial stresses of the box girder's cross-section and of each slab.

Example Analysis
Figure 8 shows the details of the simply supported concrete box girder considered in this example.The length is L = 4.0 m, a uniformly distributed load of 14.2k , which is applied to the top surface of the web along the length direction.The elastic modulus of concrete is 38.1 GPa, and its Poisson's ratio is 0.2.The weight of the box girder is not considered.
(a) Figure 9 shows the axial force along the length direction of a simply supported concrete box girder.The result shows that when the conventional method is employed for calculation, the box girder experiences relatively large axial forces.In contrast, the axial force obtained through the method proposed in this study is almost zero.A, B,  C, D, and E represent the joint between the top slab and the web, the midpoint of the top slab, the edge of the cantilever slab, the joint between the bo om slab and the web, and the midpoint of the bo om slab, respectively.In Figure 10a, the maximum stresses of the A, B, and C obtained from the conventional method are −0.653Mpa, −0.568 Mpa, and −0.601 Mpa, respectively, while the maximum stresses of the A, B, and C obtained from the proposed method are −0.685Mpa, −0.596 Mpa, and −0.631 Mpa, respectively.As shown in Figure 10b, the maximum stresses of the D and E obtained from the conventional method are 1.17 Mpa and 1.08 Mpa, respectively, while the maximum stresses of the D and E obtained from the proposed method are 1.14 Mpa and 1.05 Mpa, respectively.As shown in Figure 11, the maximum stress difference ratio of the upper flange along the cross-section is 4.73%, whereas the maximum stress difference ratio of the lower flange is −2.54%.In Figure 11b, in the upper flange, the position with the highest stress is A, reaching a maximum of 7.7%; this ratio at other positions is approximately 4.7%.Consequently, neglecting the axial equilibrium condition results in an underestimation of stress on the upper flanges while overestimating stress on the lower flanges.

Parametric Study
For a be er analysis of the influence of the change in neutral axis height on the stress of simply supported concrete box girders, this section focuses on analyzing the following factors: the ratio between the half-width of the top slab and the span of the box girder Based on the analysis of the above section, the mid-span axial stress difference ratios between the upper and lower flanges of the box girder section, as well as their variation with different parameter ratios, are discussed below.

Effect of 1 / b L
Taking 1 b = 55 cm, L is 110, 55, 36.7, 27.5, and 22 cm, respectively.The 1 / b L ratios are 0.05, 0.1, 0.15, 0.2, and 0.25, respectively.As shown in Figure 12, the stress difference ratios between the upper and lower flanges increase with an increase in 1 / b L .The stress difference ratio of the upper flange reaches a maximum value of 8.1%, and the numerical value of the lower flange does not exceed 5%.As shown in Figure 13, the stress difference ratios of the five positions gradually increase with the increase in 1 / b L , reaching 8.1%, 6.6%, 6.8%, −4.3%, and −3.2%, respectively.Consequently, the change in 1 / b L has a greater effect on the stress difference ratios of the upper flange and a lesser effect on the lower flange.Therefore, it is necessary to focus on the stresses of the upper flange of the box girder, especially at the junction between the top slab and the web (position A).   14a, the stress difference ratios of the upper flange increase with the increase in / , reaching a maximum value of 5.6%.As shown in Figure 14, when / ranges from 0.25 to 0.75, the stress difference ratios of the lower flange decrease with the increase in 3 1 / b b , and when 3 1 / b b ranges from 0.75 to 1.25, the difference ratios of the lower flange increase with the increase in 3 1 / b b , but the numerical value does not exceed 3%.As shown in Figure 15, the stress difference ratios at five positions increase with the increase in 3 1 / b b , with the stress dif- ference ratios at A, B, and C being more than 5%.Therefore, the change of / affects the stress difference ratios of the upper flange more than the lower flange./ b b ranges between 0.8 and 1.0, the stress difference ratio of the upper and lower flanges increase with an increase in 2 1 / b b , with the stress differ- ence ratio of the upper flange reaching a maximum value of 12.5% and the lower flange not exceeding 5%.As shown in Figure 17, the stress difference ratios at five positions first decrease and then increase with the increase in 2 1 / b b .The stress difference ratios at A, B, and C are very high, whereas those at D and E are relatively small, not exceeding 5%.Therefore, the change in 2 1 / b b has a greater effect on the stress difference ratios of the upper flange than the lower flange./ b b on the stress difference ratio is larger.Therefore, when analyzing the stress of a box girder, it is necessary to consider the axial equilibrium condition, particularly the box girders that have short spans, long cantilevers, and short bo om slabs.

Conclusions
The method proposed in this study is able to calculate the axial force and axial stress of simply supported box girders more accurately, considering the effects of axial equilibrium and shear deformation, to solve the shortcomings of the conventional method.According to the aforementioned research, the following conclusions can be drawn: (1) The proposed method in this paper is able to calculate more accurately the axial force, axial stresses, and deflections on box girders.The deflection difference ratio in simply supported RC box girders is a maximum of 0.60%, while the stress difference ratio is a maximum of 7.72% and a minimum of 2.51%, so considering the axial equilibrium conditions has more effect on the stress than on the deflection.(2) The maximum stresses obtained from the method in this study are −0.69Mpa in the upper flange and 1.14 Mpa in the lower flange, while the maximum stresses obtained from the conventional method are −0.65 Mpa in the upper flange and 1.17 Mpa in the lower flange.The maximum ratios of the difference in stresses obtained from the two methods for the upper flange and the lower flange are 7.72% and 3.91%, respectively.Therefore, if the axial equilibrium condition is not considered, the conventional method underestimates the stresses in the upper flange and overestimates the stresses in the lower flange.(3) The influence of the height-to-width ratio and the box girder thickness-to-web thickness ratio has a small effect on the stress difference ratio, which is found to be negligible.The stress difference ratio is significantly influenced by the ratio of the box girder length to the half-width of the top slab, the ratio of the half-width of the bo om slab to the top slab, and the ratio of the width of the cantilever slab to the top slab.
Particularly in the upper flange, with a maximum of 8.1%, 5.7%, and 12.5%, respectively.Therefore, for the short span length, the width of the bo om slab is much smaller than the top slab, and for the long cantilever slab of the simply supported box girder, we should particularly consider the axial equilibrium condition and calculate their stresses more accurately to ensure the safety of the engineering quality.
In the following study, the method can be applied to the box girder structure with a more complex structure and force conditions, which has great practical significance.Foundation of China (No. 52239009), and the Natural Science Foundation of Hubei Province (No. 2022CFD026).
are the widths of each slab, w t is the thickness of the web, f t is the thickness of the flange, andL and H are the length and height of the box girder, respectively.In terms of the Cartesian coordinate system {O, X, Y, Z}, the box girder is symmetric around the Y-and Z-axes, and the origin O can be set at any height of the web.

Figure 1 .
Figure 1.Cross-section of a simply supported RC box girder.

Figure 2 .
Figure 2. Components of the cross-section displacement.

 ; and 2 n
are the integral constants;  is a characteristic matrix com- posed of characteristic terms of matrix1 T t is the eigenvalue of matrix 1 T t  .
shows the cross-section of the box girder considered for validation.The simply supported RC box girder has a length of L = 4.2 m and supports a vertically uniformly distributed load of 16.8 y p  kN/m, as shown in Figure 3b.The load is symmetrically applied to the top surface of the web of the box girder.The layout of the reinforcing steel bars is shown in Figure 3c.The elastic moduli of concrete and steel are 34.5 GPa and 200 GPa, respectively.The weight of the box girder is not considered in the calculations.

Figure 3 .
Figure 3. Simply supported RC box girder used for validation: (a) cross-section dimensions; (b) longitudinal view of applied loads; and (c) reinforcement layout (units: mm).

Figure 4 .Figure 5 .
Figure 4. Discretization of the finite element model for the simply supported box girder: (a) concrete and (b) reinforcement.

Figure 6 .
Figure 6.Axial stresses of the simply supported RC box girder flanges: (a) upper flange and (b) lower flange.

Figure 7 .
Figure 7. Deflection of the simply supported RC box girder.

Figure 9 .
Figure 9. Axial force along the box girder.

1 / 1 /
b L ; the ratio between the cantilever slab and the half-width of the top slab 3 1 / b b ; the ratio between the half-width of the bo om slab and the half-width of the top slab 2 1 / b b ; the ratio between the height and the half-width of the top slab h b ; the ratio between the thickness of the flanges and the width of the web / f w t t .To avoid interdependence between different parameters, each analysis involves only one parameter variation while keeping the others constant.Based on the above example analysis, the elastic modulus of concrete is 38.1 GPa, and Poisson's ratio is 0.2.Taking 1 b = 55 cm, w t = 8 cm, f t = 7 cm, 2 b = 40 cm, and 3 b = 45 cm, the length of the box girder is L = 4 m, and the height of the box girder is H = 45 cm, respectively.The relative position parameters of each slab of the box girder are defined as follows:

Figure 13 . 1 b
Figure 13.Effect of / on of five positions.

Figure 15 .
Figure 15.Effect of / on of five positions.

Figure 17 . 1 b 1 / h b , and when 1 /then decrease with the increase in 1 / 1 /Figure 18 .
Figure 17.Effect of / on of five positions.

Figure 19 .
Figure 19.Effect of ℎ/ on of five positions.

Table 1 .
η and η1 along the length direction of a simply supported box girder.