CABLE-STAYED BRIDGE LOADS CAUSED BY TRAFFIC CONGESTION ON THE DECK MEASURED WITH BRIDGE MONITORING SYSTEM

. Structural safety of a bridge depends, among other things, on the number of vehicles passing on its deck, their weights and distribution of loads to their axes. A large number of vehicles can accumulate on the bridge in a controlled state, i.e., during an acceptance test including load testing, and during traffic congestion on the bridge, which is a fortuitous event addressed in this paper. The paper deals with the analysis of load intensity on one bridge carriageway when it is fully and randomly filled during traffic congestion. The influence functions of the forces in the cables are used to determine the effects of loads exerted by the vehicles moving at very low speed. Effects of such loads are studied considering an exemplary cable-stayed bridge. Since the measurement basis was limited, the iterative algorithm was used. It consists in shortening the girder sections under analysis to the area appropriate for determining the load in each successive step of iteration. Ineffectiveness of the traditional analysis show that the load intensity caused by traffic congestion is relatively high. It has been demonstrated that the matrix method may be successfully used to determine the real load of bridges on the basis of selected parameters measured in the bridge structure, also for complex scheme bridges, including the cable-stayed bridges.


Introduction
The loads assumed in the structural analysis during the design stage of a bridge structure are usually of simple homogeneous form (uniformly distributed load combined with some concentrated forces). For large bridges, the analysis of loads with special (military) vehicles is also used to assess their exact load capacities. Over-sized vehicles, such as special transport, are also considered as a separate group of live loads on the bridges.
The intensities of the aforementioned groups of variable loads (LS1 … LS5) are compared for a cable-stayed bridge over the Vistula River in Płock, Poland (Figures 1 and 2). It is assumed that the load occurs on one of the two carriageways, each 8 m wide. The review of various types of loads presented above allows concluding that the nature and intensity of actions on the supporting structures of the bridge can be quite diverse. Real loads on road bridges, by their very nature, include random features, as they are composed of actions caused by many vehicles of various sizes and weights. This paper deals with the analysis of a special case of bridge operation -traffic congestion on the roadway caused by road collision near the bridge. Under these circumstances, traffic congestion occurs when the vehicles of unknown masses and geometry accumulate on one carriageway of

Issues in determining the load of bridges during their standard operation
Many bridges in the world are equipped with various measuring systems. Some of them are of permanent nature, while others are used temporarily. Extensive research has been carried out for several decades in order to develop a methodology for determining service loads based on the measurements of the selected parameters during bridge operation. This area of bridge engineering is called 'weigh in motion' (WIM) (González et al., 2008;Helmi et al., 2014;Lydon et al., 2016;Machelski & Hildebrand, 2015;Rowley et al., 2009).
Usually, the influence lines of certain selected parameters are used, e.g., the influence lines of strains of structural elements or their displacements. The influence line of a selected static parameter allows determining the position and value of the load, however, it is possible in case of the simplest bridges only. Determining the bridge load during its operation becomes a lot more difficult problem in case of fast movement of successive vehicle axes close to each other. In this case, there are vibrations and overlapping influences from successive axes, which could make it difficult to identify the load correctly. Very fast travels of successive vehicles over a bridge cause vibrations, which do not damp immediately and do not allow for direct separation of loads in the records of the measured parameters. There are, however, algorithms which allow coping with this problem, these algorithms are based on precise modelling of bridge elements and the use of the firstorder Tikhonov regularization method (González et al., 2008). Another research problem related to WIM is assessment of the accuracy of OF ROAD AND BRIDGE ENGINEERING 2 02 1 / 1 6 ( 2 ) determining vehicle weights based on the measurements of strains at selected points of a bridge structure and simultaneous assessment of the movement speed. Research works and analyses that have been carried out are relatively simple in case of bridges with beam scheme and girder structure of the span. It has been pointed out that some methods used to determine the gross vehicle weight (GVW) reliably apply only in the situation where the bridge span is at least twice as long as the vehicle (Helmi et al., 2014). As it is desirable to improve the accuracy of WIM measurement systems, attempts are made to use new types of sensors, but not only (or mainly) the strain gauge sensors including the wireless ones. These actions have different goals: from developing light portable measurement systems, which could be installed in successive selected bridges, to elaborating the systems of high capacity to detect all axes of vehicles driving over the bridge, using also the visual methods (Lydon et al., 2016). It should be mentioned that in a series of cases, while determining vehicle weights on a bridge, real load from a truck is replaced with uniformly distributed load. Determination of loads from each axis is a more difficult task in WIM systems. Likewise, in case of beam-scheme bridges, it is a relatively simple task to develop a reliable measuring system, but it is more difficult to analyse measurement results for complex structure bridges, e.g., for cable bridges, including cable-stayed ones. However, research work has been already underway for this class of bridge structures (Machelski & Hildebrand, 2015). Further in the paper, an exemplary analysis and determination of loads on cable-stayed bridge are presented.
The presented analysis describes a rather special case because: • a permanent measuring system was used, which was not intended to directly determine vehicle weights; • the existing force sensors in stays were used, but not the strain gauges; • the static scheme of the bridge is not of the beam type, it is rather a multi-span cable-stayed one; • the vehicles on the bridge were not in motion, they took almost constant positions over several minutes; • the vehicles on the bridge represented a random sample of the bridge load with vehicles, hence a comparison could be made between the existing load and standard or other loads; • the influence functions and matrix approach were used; • three various degrees of local structure discretizing were analysed. The matrix method was used in analysing the structure due to its compact nature and complexity of results obtained (Papadrakakis & Sapountzakis, 2018). It should be added that despite relatively Maciej Hildebrand long-lasting experience in constructing cable-stayed bridges over the world and in spite of existing complex systems of computer modelling and structural analyses, further analytical work is being carried out to improve the theoretical description of this class of structures. The works of Straupe & Paeglitis (2012) and Arellano et al. (2019) can be mentioned as an example.

Description of the bridge and measurements provided by the monitoring system
The bridge examined was put into service in 2007 (Hajdin et al., 2004). Its total length is c. 1.71 km. It includes three main sections, the middle one being a five span entirely steel structure supported with stay cable set forming a single plane. This stayed section is 615 m long. Theoretical length of the longest span is 375 m. Due to the configuration of the bridge piers and appearance of short side spans, the structure was anchored to the piers using cables. Bridge pylons reaching about 63 m above the deck are of steel column type rigidly connected to the span structure. The superstructure rests on the bearings and reinforced concrete piers. The towers and girders are box type constructions, the tower is a rectangular single-cell box, while the spans are three-cell trapezoidal boxes with side cantilevers. The deck is an orthotropic steel plate. On the bridge, there are two roadways, each more than 8 m wide, and two sidewalks, each 2.5 m wide. The total deck width is c. 27.5 m and the height of the box girder is about 3.5 m. The bridge includes 28 pairs of Freyssinet system stay cables.
In view of the large size of the bridge and the record span of the main girder on the national level, and also considering the important location of the bridge in the national transport network, the bridge has been equipped with the monitoring system with 25 sensors (Hildebrand et al., 2008). Measurements were made for the following values: • forces in 8 cables, • rotational displacement of both pylon tops from the initial position, • strains of the steel structure in 10 points, • velocity and direction of wind in two points (on the deck and on the pylon top), • temperature of steel structure in two points of the main girder and temperature of the air inside the girder.

Forces in cables versus load intensity of bridge roadway
In order to calculate the internal forces in the main structural elements of the cable-stayed bridge under consideration, the flat 2D model can be used (Machelski & Janusz, 2017), as shown in Figure 2 and also in Figure 3 in a more detailed form. This model uses bar/ truss elements to discrete the girder, pylons and cables (Papadrakakis & Sapountzakis, 2018;Zienkiewicz & Taylor, 2000). In this way, it is possible to create the influence lines of the internal forces. The lines are suitable for calculations for the case when loads change their positions along the roadway axis on the bridge. The effectiveness of this model was verified in the paper by Machelski & Hildebrand (2015). The calculation algorithm makes use of the axial force influence functions in cables f(x) in the forms shown in Figures 4 and 5. In these diagrams, the characteristic points are the positions of piers at distances x from bridge beginning: 0 m, 60 m, 120 m, 495 m, 555 m, 615 m. The cables 24 and 28 are symmetrical to elements 5 and 1, respectively. The cables Nos. 1, 5, 7, 12, 13, 14, 24, 28 were selected for analysis as they include the force sensors within the bridge monitoring system (SHM).
In the flat 2D model of the bridge, the roadway load q(x) is related to the axial force in cable W k and influence function f k as given in the formula  When distributed load q(x) is fixed, different values of W exist in particular cables because influence functions f(x) have different shapes. Both functions in Eq. (1) are implicit. Due to the fact that the paper aims at determining the load function q(x) using exclusively the forces in cables W, relation (1) is not particularly useful for this solution. Usually, the internal forces are calculated on the basis of a given load, while this paper aims at completing the reverse task. A challenging element of the task is to map a complex form of the load function q(x), which depends on the random layout of vehicles on the bridge.  The task of resolving the continuous system given by Eq.
(1) was transformed to the task of resolving a discrete system where the span is divided into segments with uniformly distributed loads, as shown in Figure 6.
Hence, Equation (1) is transformed into a system of the following equations F · q = w (2) and in detail Loads distributed on the bridge roadway segments are compiled in vector q = col{q 1 q 2 q 3 q 4 q 5 q 6 q 7 q 8 } (4) and the accompanying forces in cables as in the equation In the system of Eq. (3), the elements of matrix F are the surface areas of the influence function f(x) sectors calculated like in the case F k j from the Simpson's rule as follows: Graphic interpretation of Expression (6) is given in Figure 6 where points jand j + are the extreme points of the segment of length c. Eight load segments were assumed (Figures 6 and 7) due to the number of cables (eight) with measurements of axial forces. The division to segments over the length of the selected span is basically arbitrary, which is used further in the paper.
Function q(x) is sought for, assumed in the form of segments with constant values (modelling segments of uniformly distributed load), included in vector q. It is impossible to find q j directly from the system of Eq. (3) when the forces taken during the measurements on site are assumed as values w, that is: 1 5 7 12 13 14 24 28 .
This conclusion results from the differences between the real load on site, i.e., the complex function q(x), as in Eq. (1), and postulated uniformly distributed loads (constant over the segments with length c), in vector q, as in Figure 6. When low number of segments is used (rare dividing grid of the girder), as it happens in this analysis, the values of vector w are obtained on the basis of vector q instead of the values taken on site w o as in Eq. (7).
For the special case, when it is assumed that uniformly distributed loads of constant value marked as q(k) (specific for each stay k) exist over the full length of the span, then from Eq. (3) eight independent equations for each cable of the form are obtained: Hence, F k in Eq. (8) is the surface area under the line of influence function of axial force in cable k, as in the formula When the forces in cables from measurements W k o and those calculated with Formula (8)  In the iterative algorithm used to calculate the effect of traffic congestion on the bridge, a special assumption that load q(x) is close to the uniform scheme along the whole bridge is made. For this reason, the values calculated from Formula (10) are used here as initial values of the iteration process. They create vector The terms of vector q k cannot be directly referred to the elements of vector q as the first ones are connected with cable numbers, while the second ones -with division to discretization segments. However, the elements of vector q k allow for initial estimation of terms in vector q.

Iterative algorithm of the solution
Due to small number of cables (eight) covered by the measurements and a complex scheme of the construction geometry, it is necessary to use iterative procedure consisting in narrowing the calculation area. Figure 7 illustrates the segments (eight) of span division as used in the calculation schemes, which are marked in this paper as follows: •  is calculated based on the differences between forces in the cables. When the differences between terms of vectors w and w o are small enough, which is reflected in coefficient r k , the iteration in scheme G is terminated.
Transition to the calculation scheme M consists in modifying vector w o , hence in the change of terms on the right side of Eq. (3) to the value: where j = 1, 2, ... 7, 8 in scheme G. Because of shorter division grid segments, the terms of matrix F in Eq. (3) are modified. Coefficient r k is also verified in this calculation scheme. Successive transition to the next scheme C consists in modifying vector w S by changing the terms on the right side of Eq. (3) to the value where j = 1, 2, ... 7, 8 in scheme S. The terms of matrix F in Eq. (3) are also modified in this scheme. Due to the differences between the accurate solution given by Formula (1) and approximate solutions (3) consisting in approximating continuous function q(x) to constant values q j on calculation segments, as in Figures 6 and 7, the total zeroing of coefficient r k is not attained. Even if the density of division into discrete elements of the model would be greater, the error would still remain. It should be remembered that calculations are made using the 2D model, while the objective is to calculate the forces in the bridge cables and to compare them with the real values measured in the physical 3D structure loaded with complex loads from vehicles q(x).

Accuracy of the solution
Resulting from Formula (1), both functions, q(x) and f(x), under the integral implicate each other. To get solution of the equation system (3), it was assumed in the iterative algorithm that the loads are constant over the segments of length c. According to such assumption, successive simplifications of Formula (1) for the calculation of component W k (over the section j) can be used, as in terms of the equation The final portion in Eq. (15) includes the value calculated from Formula (6) where the term f k (x) is as shown in Figure 6. Hence, in such formulation, F k j can be a load-independent term of matrix F, like in Eq. (3).
In a general situation, when approximation of function q(x) has the form of the second order parabola, like function f(x), relationships of complex forms are obtained. Assuming a parametric approach to coordinates ξ = x/c, these functions acquire the forms as follows and On this basis, the general formula is obtained  (16) and (17). In the exemplary calculation further in the text, the convex form of function f(x) in the central part of the main span, as in Figure 5, was assumed.
Geometric parameters of function f k with ordinates assumed as in Figure 6 are summarized in Table 1. On this premise, a simplified form of Formula (18) In a specific case, when load q j (x) is a linear function, so in Eq. (17) a = 0, Formula (6) can be still conversed to the entangled relation of the influence function and load: where q j-and q j+ are the values at both ends of segment j, similarly as in case of function f k shown in Figure 6. An example of such load function is discussed in the next section of this paper, where the load is considered as an open polygon instead of being segment-by-segment variable as in Figure 6.

Parameters from Eq. (16) Ordinates
In case of function q(x), four characteristic forms are considered: • that of constant value, as in Formula (15); • that of linear variability, as in Formula (20); • that of convex parabola, analogous to function f(x); • that of concave parabola reverse to the convex one. Table 2 summarizes parameters of the functions considered. The last column of Table 2 provides calculation results from Formula (19). When these values (related to table row no. 1) are compared, the influence of the shape of q(x) diagram on the calculation result is visible. Obviously, the results will be dependent on the quotient of the ordinate values of function q(x). This example shows that assuming the segments of constant load q j in the paper, as determined by necessity to run calculations with iterative algorithm, can significantly disturb the shape of function q(x). The diagrams provided further in this paper show that the form of function q(x), by its very nature, is very complex. This effect -variation of the shapes of q(x) function -is not considered in the results of the calculation algorithm applied. This leads to system errors, which are located in the terms of matrix F. This note is the basis for justifying the conclusion on the inefficiency of the direct solution of the system of equations (3), and also on the shortage of zeroing the iterative coefficient from Formula (12).

Exemplary calculations
The results of calculations related to the measured forces included in the following vector are the example of the implemented iterative procedure: w o = col{3.5 2.7 1.7 3.8 3.9 3.2 2.5 3.8}, MN.
(21)  The values of the forces correspond to the load caused by the traffic congestion shown in the photos (Figure 8). Both photos were taken by the video monitoring system. At that time, like each day, the measurements were taken by the bridge structure monitoring system. As it can be seen, the right roadway was free of any essential load. Photo b) was taken about 1 minute after Photo a) was taken.
This was the basis for calculating the elements of vector q k using Formula (11); they are provided in Table 3.  Table 3 should be read as follows (for instance): in the column marked 1 (number of the stay) there are values of loads calculated on the basis of the influence function of stay No. 1 and the real force measured in stay No. 1. And so on.
When the values in columns of vector q k are compared, a small change in terms of the right side of Eq. (3) for schemes G and M can be observed. In case of scheme C, there is a higher modification of free terms. When the aforementioned algorithm was used, the resultant values of q after several iterations were as given in Table 4. Iteration was terminated when coefficients r k from (12) amounted to several percents. Obviously, there are no relationships between the figures from the columns of Table 4 referring to schemes G, M, C, as according to Figure 7, these are segments located in other places. However, when the values are compared, it may be noticed that the results from Table 3 used in the iteration process are helpful. For instance, analyzing scheme G, term q = 60.8 kN/m for segment No. 4 can be found. The same value can be found in Figure 9, curve G, above the abscissa 262.5 m, i.e., in the mid-point of segment No. 4 in scheme G. Similarly, analyzing scheme C, term q = 76.8 kN/m for segment No. 4 can be found. The same value can be found in Figure 9, curve C, above the abscissa 296.25 m, i.e., in the mid-point of segment No. 4 in scheme C (Figure 7), etc. Figure 9 illustrates the results from Table 4. Compatibility of calculation results for schemes G, M, C is a positive feature of these diagrams. Obviously, the best representation of function q(x) is for the range 210 m < x < 390 m when the segments in scheme C are of minimum length c = 22.5 m. Striving for further division of segments in scheme C seems unjustified due to iterative mode of running to the results. On the other side, a change in locating the modelling scheme C, e.g., 120 m < x < 210 m or 390 m < x < 495 m is recommendable. An attempt to include the range 0 m < x < 120 m is ineffective due to small values of ordinates in diagrams f(x) given in Figures 4 and 5.
Due to relatively small values of ordinates of the influence line near the supports, the results obtained for the sections of the span located very close to the piers are less reliable, see Figures 4 and 5.

Numerical analysis of the central part of the river span
The example above referred to the initial phase of traffic congestion on the bridge. During the later phase of congestion (from several minutes to two hours later), the cable forces changed due to the changes in vehicle positions. Figure 10 illustrates changes of forces in particular cables in the next moments of time.
Successive schemes, i.e., C1, C2, C3, C4, occurred in the subsequent time intervals. Table 5 provides selected results of the measurements. Such moments were selected when the highest forces took place in cables 1, 5, 7, 12 (marked with bold). Rows of Table 5 include the terms of vector w o in Formula (15).   Figure 10 illustrates the diagrams of loads q(x) for the central part of the main span 210 m < x < 390 m from the calculations in scheme C. Significant changes of load intensity during successive congestion phases considered can be seen here. In case when the values of q are close to 100 kN/m, the situation is interpreted as the effect of heavy vehicles concentration in some parts of the bridge. Intensity of this load is similar to the values used during acceptance testing of the bridge. The results shown in Figure 10 can be compared with the dead weight of the span. The dead weight of the span is about 150 kN/m.
It should be noted that there were no substantial changes of the temperature of the structure, as traffic congestion appeared in the afternoon in December. Following the analyses of the records of the forces in stays during the whole day, when the described traffic congestion event occurred, the following conclusion was derived: The influence of loads due to traffic congestion on the forces in stays was at least 3 (usually up to 5) times larger than the influence of temperature changes (during the observed period of congestion). So the problem of temperature influence in this task was assumed as a not significant factor and was not included. However, in a general situation (load event lasting for many hours, especially in the warm season), the influence of temperature changes in such tasks cannot be disregarded. Table 2 is an example of linear function applied on the analyzed segment with length c. In this approach, Relation (6) was formulated as the entanglement of functions q(x) with f(x).

Row 2 in
and W c f q f q f q q k k k k j j j j j+ j j j Assuming q j as the value being analysed from these two formulas, the following term can be created: It results from Figure 11 that for this calculation model, the value obtained in (25) covers two segments with length 2c with the common node j. For this approach, it is very important to separate the loads existing in Formulas (22) and (23). It enables to formulate a system of equations such as Eq. (3) with vector of load q and that of free terms w. However, the interpretation of the essence of terms q j is quite differentterms q j are the values at points of division to segments ( Figure 11) but not the constant values over the whole segment with the central point, as in Figure 6. The terms of matrix F are also subject to change -previously, they were calculated from Formula (6), but now from Eq. (25).
The data from vector (21) were used as exemplary calculations using the presented algorithm. The results of the calculations are given in Table 6. Similar values can be observed comparing the results for scheme C from both Table 4 and 6. The same grids dividing the span into Figure 11. Function of the influence and load in modelling scheme M segments are used in both calculations. The differences refer to the form of the approximated loads as in Figures 6 and 11. Similar conclusion was drawn comparing the values in rows 1 and 2 in Table 2. As in this paper the iterative algorithm is used instead of resolving the system of Eq. (3), it is possible to use larger (or smaller) number of displacement vector terms (e.g., nine as in Figure 11), than the number of measurements in the cables (eight). Hence, the calculation scheme assumed in Figure 11 can be presented in the following form This example leads to the conclusion that basically any division grid over the span length can be used, as well as the segments of various lengths. Therefore, the schemes given in Figure 7 can be combined into one calculation system. The result of iteration is an open polygon of load, as shown in Figure 11. However, it should be noted that the real load function is a complex one and the calculations just approximate it.

Conclusions
According to the discussion presented above, it is possible to map the forms and values of loads existing on the bridge with a relatively complex structural layout on the basis of the forces in several cables. Obviously, this representation demonstrates only limited accuracy. The procedure consisting in making use of the influence function and matrix approach is known for a long time (Papadrakakis & Sapountzakis, 2018;Zienkiewicz & Taylor, 2000). However, some innovation lies in the use of various levels of discretization of the structure or its fragmenting into smaller elements.
An iterative (three-stage) algorithm was used; it consists in shortening the segments under analysis over the bridge length as the area to determine the roadway load in the successive iteration step. Efficiency of such solution was demonstrated in exemplary analyses. The results of numerical analyses proved high intensity of load recorded during traffic congestion as compared with the loads specified in design criteria or codes (for instance, Eurocode). Hence, the analyses carried out provide a substantial result -an opportunity to compare the load of vehicles that create traffic congestion with the standard load, both with the one that is currently known, included in Eurocode, and the load from the previous standard, which had been used to design the bridge in Płock.
Comparison of the loads (see Introduction of this paper) is summarized in Table 7. The data shown in Table 7 can be compared with self-weight of the span, which is equal to approximately 150 kN/m. The iterative algorithm assumes that q(x) is nearly uniform. The results given in the example prove the efficiency of the calculation procedure. It was attained using the iterative algorithm. The best accuracy was reached in the central part of the bridge, which results from relatively high values of influence line ordinates f(k). The r coefficients of several percents were reached in the example analysed.
The areas where loads q were close to 100 kN/m were interpreted as the effect of heavy vehicles concentrating in the bridge area under consideration (Figure 8). Local increase of this load, which was reaching 80…100 kN/m (Figure 10) as demonstrated in the above analysis, is comparable with the load generated by military transport (87.4 kN/m) or by trucks used in the acceptance testing (117.2 kN/m) of the bridge in question.
The following conclusions can be made considering the results of the calculation example presented above: • the matrix method is an effective way to determine the causes (operational loads) examining the effects (forces in the selected stays); • it is possible to determine not only the average load for the whole span, but also to point out span sections of various load values; • solution of the task is possible also for the limited data taken from properly located sensors of suitably chosen parameters. It should be stressed that the analysis aiming at determining bridge load during its operation (WIM) may be applied not only on the beam structures, but also on more complex objects, including cable-stayed bridges; examples of such analyses for cable-stayed bridges are hardly available in the current literature on the subject. The above-described case, including the proposed methodology of the procedure, bridges this gap to a certain extent.