Efficient Design Optimization of Cable-Stayed Bridges: A Two-Layer Framework with Surrogate-Model-Assisted Prediction of Optimum Cable Forces

Abstract

Cable-stayed bridges have commonly been built for crossing large-span obstacles, such as rivers, valleys, and existing structures. Obtaining an optimum design for a cable-stayed bridge is challenging, due to the large number of design variables and design constraints that are typically nonlinear and usually conflict with each other. Therefore, it is a reasonable alternative to turn the large and complex optimization problem into two sub-problems, i.e., optimizing the internal force distribution by adjusting the cable prestressing forces, and optimizing the other sizing or geometrical parameters. However, conventional methods are lacking in efficiency when dealing with the problem of optimization of cable forces in the first sub-problem, under the circumstance that iteration between the two sub-problems is required. To address this, this paper presents a surrogate-model-assisted method to construct a cable forces predictor ahead of the structural optimization process, so that cable forces can be effectively predicted rather than optimized in each iterative round. Additionally, B-spline interpolation curve is adopted for variable condensation when sampling for the surrogate model. Finally, the structure optimization in the second sub-problem is performed by leveraging an optimization program based on particle swarm optimization method. The performance of the proposed framework is tested with a practical engineering application. Results show that the proposed method showcases good efficiency and accuracy. The theoretical raw material consumption of the towers and the cables is 32% lower than the original design.

1. Introduction

Cable-stayed bridges have been widely constructed to span roads and rivers. As a highly redundant structure, cable-stayed bridges have advantages in terms of the stiffness, wind-load resistance, maintenance, and span-crossing ability. As a result, cable-stayed bridges have exceeded a span of 1 km in just 60 years since the first cable-stayed bridge was built. The Chang-Tai Yangtze River Bridge in China, with a total length of 5.3 km and a main span of 1176 m, is the world’s longest cable-stayed bridge for both highway and railway. The construction of such a huge engineering structure is extremely challenging and costly. Therefore, it is becoming increasingly critical to carry out detailed structural optimization and comparison at the design stage.

For a practical engineering project, the implementation of bridge design optimization can face a couple of challenges. First, the involved design parameters for optimization of a cable-stayed bridge are typically high-dimensional. As a highly redundant structure, there are numerous design parameters to be optimized in the bridge design process. Variables involved can be classified into mechanical, sizing, geometrical, and topological [1]. The optimization problem becomes highly dimensional due to the increasing number of variables, which hinders efficiency and consumes more time. The interdependency among design variables exacerbates nonconvexity and nonlinearity, making the problem difficult and challenging. In addition, various constraints need to be considered in the optimization design of cable-stayed bridges. Due to the complexity of real-world environments and bridge-operating conditions, engineers should consider the varieties of site, structure, material, member, section type, load cases, and the complex verification items of specification. These considerations, whether linear or nonlinear, make the feasible domain very limited within the design space. Moreover, the mechanical behavior of cable-stayed bridges must be properly simulated. Cable-stayed bridges transmit loads with cables. They are relatively flexible structures compared to girder and arch bridges. Therefore, the geometric nonlinearity caused by sag effect of cables, large displacement effect, and p-Δ effect, should be considered. Though current commercial software can handle these non-linearity effects by finite element analysis, the simulation can be time-consuming and undermine the feasibility of optimization.

To address these challenges, a lot of research has been conducted to advance the computational methods or strategies for the design optimization of cable-stayed bridges. Early studies conducted by Feder [2] introduced an optimality-criteria-based method to determine the prestressing forces of cables in steel bridges. Similar early studies involved plenty of assumptions and addressed the problem using simplified mathematical formulas. Sung et al. [3] minimized the total strain energy expressed as a quadratic function of the post-tensioning cable forces with an influence matrix. Baldomir et al. [4] optimized cable areas for a long-span steel bridge with the finite differences sensitivity analysis method and solved the problem through a gradient-based sequential quadratic programming algorithm. A three-stage algorithm was presented by Ha et al. [5] to optimize the cable prestressing tensions with a nonlinear inelastic analysis. Besides the above work on optimizing single type of variables, research has also been carried out to pursue a more comprehensive “optimum design” of the bridge. The structural design problem is formulated with various variables including not only mechanical but also sizing, geometrical, and topological variables. Lute et al. [6] proposed an optimization method for cable-stayed bridges that utilized a genetic algorithm to minimize costs while considering geometrical parameters and cross-sectional dimensions as design variables. A support vector machine was utilized for constraint verification, and the presented method was proven to be accurate and computationally efficient for prediction purposes. Gao et al. [7] obtained the optimum design of prestressed concrete bridges. Design variables included the number of prestressing tendons, cable forces, cable areas, and girders’ and towers’ sectional dimensions. Cid et al. [8] examined multi-span cable-stayed bridges while considering geometric nonlinearity effects. The variables included anchorage positions, cable forces, and cable section areas. The SQP algorithm was utilized to minimize the total cost of steel, and sensitivity analysis was conducted using the finite difference method.

However, for the implementation of optimizing a practical cable-stayed bridge, handling all types of variables simultaneously is typically not the most effective strategy, because this formula can significantly lead to the increase in the dimensionality, nonconvexity, and nonlinearity of the problem. Moreover, it should be noted that optimizing the distribution of internal forces is not intrinsically contradictory to the optimization of the cost. Taking advantage of this feature, the optimization can be separated into a hierarchical layout, i.e., sequentially optimizing the mechanics-related variables and the other ones. Following this idea, [9] is one of the earliest works that introduced surrogate functions that reveal the potential connection between mechanical variables and other sizing and geometrical variables, therefore decoupling the cable forces optimization from the structural design optimization. To develop the functions, polynomial regression with the ordinary least square method was adopted. The necessary data for the regression were collected from a large parametric study conducted by repeating the finite element technique, while varying three parameters. The surrogate functions, expressed as quadratic polynomials, explicitly related cable forces to three variables concerning the span length, the total length, and the upper structure height. With these functions, the cable forces can be easily determined for different variable values to achieve the optimum post-tensioning distribution, which can minimize the deflection in both the deck and pylon. In later research [10], these functions are used to facilitate the optimum design of a cable-stayed bridge, considering variables such as the cables’ section areas, prestressing forces, and cross-sectional dimensions of the girder.

Although the post-tensioning functions were assessed as accurate in the previous investigation [10], there are still drawbacks when it comes to optimum design problems with more complex design conditions. Firstly, the variables included in the functions are limited in terms of the number and the types, which may lead to inaccurate predictions for the bridges if some of the other design conditions are changed. Secondly, the variables of the surrogate functions are supposed to be consistent with those of the optimum design problem. Therefore, there is a demand for reconstructing the surrogate function if some new variables that may have a visible influence on the post-tensioning distribution are introduced. Thirdly, as the number of design variables increases, the nonlinearity of the problem significantly grows. The polynomial regression will be insufficient to give accurate predictions, due to its poor anti-interference and local fitting ability, compared with other regression methods.

To this end, machine learning techniques have recently garnered significant attention. Researchers have explored various tools, including Random Forest [11], Support Vector Machines [12], and Physics-informed or Data-driven Neural Networks [13], to predict structural performance or responses. In a broad sense, as these tools are strategically designed as a cheaper-to-evaluate substitute for the original sophisticated computational model, they are also typically called surrogate models or meta-models [14]. The conceptual illustration of the construction and the utilization of a surrogate model is shown in Figure 1. Multidisciplinary applications have been conducted based on surrogate models to replace the original time-consuming processes or high-cost experiments. These methods have been successfully applied in a variety of research fields, such as the hydro-environment [15], rock and soil mechanics [16], and bridge engineering [17].

Motivated by the above advancements, the surrogate model from machine learning is used in this research to replace the surrogate functions and is expected to exhibit greater adaptability to a large number and diverse types of variables. Four regression models are studied in this work: Polynomial Regression (PR), Gaussian Process Regression (GPR), Regression Tree (RT), and Support Vector Regression (SVR). To construct the surrogate model, samples are collected by means of a full factorial experiment, concerning the modest number of target sample count and acceptable computational cost. To test the surrogate model, the predicted results and optimized results of a design point distinguished from the samples are compared. Thereafter, the surrogate model is combined with heuristic algorithms to demonstrate its efficiency in the optimum design. By introducing such a predictor, the force variables can be decoupled from other sizing and geometrical variables in the optimization routine, thus enabling a substantial reduction in the problem complexity and boosting the efficiency of optimization. Additionally, a practical application of a cable-stayed bridge project with a main span of 818 m has been used to validate the performance of the proposed framework.

The rest of the paper is organized as follows: Section 2 displays our improved formula of the two sub-problems. Section 3 presents the proposed methods and their integration in the overall procedure. Section 4 shows the implementation of the methods on a (358 + 818 + 358)-meter-long cable-stayed bridge example. Section 5 displays the optimization results and complementary checks. Conclusions are drawn in Section 6.

2. Problem Formula

In a general cable-stayed bridge design, engineers need to carry out comprehensive design and optimization from the structure configuration to the detail members. In terms of the configuration, crucial parameters including the length of the main and the side span, the width of the deck, the height of the towers, and the anchorage position of the cables should be well determined. In terms of the members, cross-sectional dimensions of the main members including the towers, the cables, and the deck are supposed to be well designed. The typical optimum design formula with all design variables optimized altogether can be stated as

$$\begin{array}{cc}\hfill \underset{\mathbf{X},\mathbf{Q}}{\mathrm{m}\mathrm{i}\mathrm{n}}\mathsf{\Phi }& ={\displaystyle \sum _{j=1}^m}{C}_j{V}_j\hfill \\ \hfill s.t.x_{il}⩽x_i& ⩽x_{iu} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} i\in N_X\hfill \\ \hfill q_{kl}⩽q_k& ⩽q_{ku} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} k\in N_Q\hfill \\ \hfill {\sigma }_{jl}⩽{\sigma }_j^{\left(p\right)}& ⩽{\sigma }_{ju} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} j\in M, \mathrm{a}\mathrm{l}\mathrm{l} p\in L\hfill \\ \hfill f_j^{\left(p\right)}& ⩽f_{ju} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} j\in M, \mathrm{a}\mathrm{l}\mathrm{l} p\in L\hfill \\ \hfill K_j^{\left(p\right)}& ⩽K_{ju} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} j\in M, \mathrm{a}\mathrm{l}\mathrm{l} p\in L\hfill \end{array}$$

(1)

where $\mathsf{\Phi }$ represents the total cost of members, $\mathbf{X}=\left\{x_1,x_2,x_3,\dots ,x_{n_X}\right\}$ stands for sizing (e.g., cross-sectional dimensions) and geometrical (e.g., side-span ratio, height-span ratio) design variables, and $\mathbf{Q}=\left\{q_1,q_2,q_3,\dots ,q_{n_Q}\right\}$ stands for mechanical (e.g., cable prestressing forces, tendon prestressing forces in reinforced concrete beams) design variables. $C_k$ and $V_k$ represent the cost coefficient and the volume of the $j$-th member, respectively. $M=\left\{\mathrm{1,2},3,\dots ,m\right\}$ is the member number set. The next two constraints represent the lower and upper boundaries of the $i$-th variable in $\mathbf{X}$ and the $k$-th variable in $\mathbf{Q}$, where $N_X=\left\{\mathrm{1,2},3,\dots ,n_X\right\}$ and $N_Q=\left\{\mathrm{1,2},3,\dots ,n_Q\right\}$. The following three constraints represent, respectively, the strength, stiffness, and stability constraints of the structure. The first constraint of these represents the upper bound of the $j$-th member’s stress response, where $p$ is the load case identifier and $L=\left\{\mathrm{1,2},3,\dots ,l\right\}$ is the load case number set. The second represents the boundary of the deflection response, and the third represents the boundary of the stability coefficient.

It is worth noting that wind or seismic resistance is also critical, especially in the design of long-span bridges; however, due to the limit of the paper length, such dynamic performance is not discussed in this paper. Therefore, the girder sizes are excluded from $\mathbf{X}$ and not optimized in the later example. Meanwhile, considering that cable prestressing forces are more dominant than tendon prestressing forces in preliminary design, $\mathbf{Q}$ consists only of cable forces in the later example.

In this optimization, $\mathbf{X}$ and $\mathbf{Q}$ are variables optimized together, where $n_q$ is generally much larger than $n_X$ due to the high density of cables. However, an optimization process dealing with $\mathbf{X}$ and $\mathbf{Q}$ simultaneously is highly complex and computationally costly, due to their nonlinear and conflicting (i.e., coupling) feature [1]. Thus, in this work, a reasonable two-layer framework is adopted to decouple them into two sub-problems in one optimization round.

Sub-problem one (Sp1) is an internal force distribution optimization problem where $\mathbf{Q}$ is optimized with fixed $\widehat{\mathbf{X}}$, as most of the “cable forces optimization” problems [18,19]. The most common and practical objective function of Sp1 is the weighted sum of members’ strain energy. The formula of Sp1 can be stated as

$$\begin{array}{cc}\underset{\mathbf{Q}}{\mathrm{m}\mathrm{i}\mathrm{n}}U=& \underset{\mathbf{Q}}{\mathrm{m}\mathrm{i}\mathrm{n}}{\int }_S\frac{M^2(S,\mathbf{Q})}{2EI}dS\hfill \\ q_{kl}⩽q_k& ⩽q_{ku}\mathrm{f}\mathrm{o}\mathrm{r}\mathrm{a}\mathrm{l}\mathrm{l}k\in N_Q\hfill \\ f^{\left(deadload\right)}& ⩽f_u\mathrm{f}\mathrm{o}\mathrm{r}\mathrm{a}\mathrm{l}\mathrm{l}p\in L\hfill \end{array}$$

(2)

where $U$ represents the total bending strain energy of members, usually including both the main girder and the towers, $E$ is Young’s modulus, and $I$ is the bending moment of inertia. The first constraint is lower and upper bounds set to ensure that the cable force magnitude does not turn negative or exceed its design strength. In the second constraint, $f^{\left(deadload\right)}$ is the deflection under the dead load case. It is set to ensure that the configuration of the bridge meets the design upon its completion.

Sub-problem two (Sp2) is a sizing and geometry optimization problem where $\mathbf{X}$ is optimized with $\widehat{\mathbf{Q}}$ determined in Sp1. The formula of Sp2 can be stated as

$$\begin{array}{cc}\hfill \underset{\mathbf{X}}{\mathrm{m}\mathrm{i}\mathrm{n}}\mathsf{\Phi }& ={\displaystyle \sum _{j=1}^m}{C}_j{V}_j\hfill \\ \hfill s.t.x_{il}⩽x_i& ⩽x_{iu} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} i\in N_X\hfill \\ \hfill {\sigma }_{jl}⩽{\sigma }_j^{\left(p\right)}& ⩽{\sigma }_{ju} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} j\in M, \mathrm{a}\mathrm{l}\mathrm{l} p\in L\hfill \\ \hfill f_j^{\left(p\right)}& ⩽f_{ju} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} j\in M, \mathrm{a}\mathrm{l}\mathrm{l} p\in L\hfill \\ \hfill K_j^{\left(p\right)}& ⩽K_{ju} \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} j\in M, \mathrm{a}\mathrm{l}\mathrm{l} p\in L\hfill \end{array}$$

(3)

Conceptual differences are shown in Figure 2.

Formula (1) has a drawback in dealing with multiple interdependent variables together, whereas the routine methods dealing with Formula (2) are overly time-consuming for the proposed iterative approach. Thus, this paper focuses on how to solve Sp1 in Formula (2) with greater efficiency and to make the proposed methods more adaptable to different and complex circumstances. Rather than directly optimizing the cable forces in each iteration round, this research seeks to solve Sp1 with a surrogate model from machine learning, by predetermining $\mathbf{Q}$ with a cable forces predictor trained in advance of the optimization.

3. The Proposed Method

3.1. Cable Forces Optimization

In the beginning of this section, the method for determining the optimum cable forces, i.e., Formula (2), is discussed, which builds a foundation for establishing the training database of the force predictor. It should be noted that a rational distribution of cable forces in cable-stayed bridges is generally smoothly distributed without sudden changes. Therefore, B-spline curve presented by French engineer Pierre Bézier, is handy for fitting the cable force distribution in cable-stayed bridges [20,21,22]. Compared to the precise optimization of each original cable force, introducing B-spline greatly condenses design variables into a few control points, while causing only a limited impact on the overall structural behavior.

Suppose a $p$-th degree B-spline curve is introduced, the original dense cable prestressing forces $\mathbf{Q}=\left\{q_1,q_2,q_3,\dots \dots q_{n_q}\right\}$ can be replaced by its control points $\mathbf{P}=\left\{p_1,p_2,p_3,\dots \dots p_{n_p}\right\}$, where $n_p$ is far less than $n_q$. The relation between $\mathbf{Q}$ and $\mathbf{P}$ can be described as

$$q_k=\sum _{i=0}^{n_p}{N}_{i,p}\left(u_k\right)p_i \mathrm{f}\mathrm{o}\mathrm{r} \mathrm{a}\mathrm{l}\mathrm{l} k\in \left\{1,2,3,\dots \dots ,n_q\right\}$$

(4)

where $u_k$ is the relative position of the $k$-th cable ($0⩽u_k⩽1$), and $N_{i,p}$ is the shape function. More detail can be found in [22].

In this work, assisted by the B-spline technique, the following design variables are considered: four fitting points ($d_1$–$d_4$) in the side span, four fitting points ($d_5$–$d_8$) in the main span, the force of the No. 1 cable ($d_9$), the range of the counterweight ($d_{10}$), and the load magnitude ($d_{11}$), as shown in Figure 3. The counterweight and its distribution range are set as variables as a complement to deal with the unbalanced self-weight between the main span and the side.

Since Formula (2) is a typical nonlinear optimization problem and the number of design variables is condensed with B-spline curve, it can be efficiently solved with several optimization methods, such as the differential evolution and surrogate-model-assisted differential evolution [17]. The sequential quadratic programming (SQP) solver [23] is adopted in this research due to its high efficiency (note that the initial point should be carefully selected to avoid local optimum solutions).

During the sampling to construct the predictor, Formula (2) is solved in collaboration with MATLAB R2018b and ANSYS 2022 R1 APDL. For each sample, modeling parameters are derived in the beginning. An initial FEM is built with these modeling parameters in ANSYS. When the cable forces and counterweight are updated during the cable forces optimization in MATLAB, the initial strains of the cables are modified for the FEM and the counterweight is imposed in ANSYS. Then, the static analysis is launched to obtain the bending energy. If convergence is checked, the result of the sample is output to a text file. The sample process is completed when it reaches the target sample count.

3.2. Surrogate Model Assisted Predictor

Because of the high computational cost to solve the cable forces optimization problem, it is appropriate to introduce the surrogate model method to predict the cable forces of the desired completion state. As shown in Figure 4, the following steps are developed in this work to train and select a proper surrogate model.

• Sampling: Generate training data (or samples) from the design of experiment (DOE) [24]. In this work, the samples are collected by means of a classic full factorial experiment. As for the cable-force predictor, the input data are the identified design variables, and the output data are the optimum cable forces taking the form of controlling points mentioned in Section 3.1. Considering that the cable forces optimization problem is well solved with the minimum bending energy method and the proper algorithm, there should be no anomalous sample.
• Training: Configure the output function based on the machine learning algorithm. The property of the chosen surrogate model can be evaluated with root mean square error (RMSE) and training time. The surrogate model with the least RMSE and an acceptable training time is the best solution. In this work, the following models are considered and compared: Polynomial Regression (PR) [25], Gaussian Process Regression (GPR) [26], Regression Tree (RT) [27], and Support Vector Regression (SVR) [28].
• Testing: Create new testing data, conduct the experiment and obtain the actual responses. Expose the new data to the output function and get the predictions. Compare the actual responses and the predictions, then evaluate the accuracy of the surrogate model. If the test fails, it means that the discretization of the variables is either not dense enough or not rational; so, the design of the experiment should be adjusted.

With the predictor constructed, the optimized cable forces in each iteration round can be predetermined in subsequent Sp2. The introduction of the surrogate model significantly speeds up the iteration process, thus enabling efficient design optimization later.

3.3. Modification Strategies

To simplify some parts of the overall optimization procedure, two modification strategies are further adopted to improve the accuracy of simulations without significantly reducing efficiency.

3.3.1. Elastic Modulus of the Cables

The first modification is made on the elastic modulus of the cables. Among the three typical expressions of geometric nonlinearity, namely, sag effect, p-Δ effect, and large geometric deformation [29], sag effect has the most obvious impact on the optimized cable forces under dead load. Instead of performing the time-consuming nonlinear analysis, the Ernst formula modifying the elastic modulus of cable is adopted.

$$E_{eq}=\frac{E}{1+\frac{EA}{12T^3}(qH{)}^2}$$

(5)

where $E_{eq}$ is the effective elastic modulus of the cable, $E$ is the material elastic modulus, $A$ is the section area, $T$ is the cable prestressing force, $q$ is the cable weight per meter, and $H$ is the horizontal projection length of the cable.

3.3.2. Section Area of the Cables

The second modification is made on the section area of the cables. Instead of identifying areas for each cable as a variable, 3 to 5 adjacent cables are grouped and share the same section area. However, there are still too many area variables that hinder the optimization efficiency. Therefore, an iterative method is adopted in this paper to estimate each cable area with a certain designated safety factor and calculate the cable force changes caused by the modification.

The steps are as follows:

• Calculate objective cable stress ${\sigma }_{obj}={\sigma }_s/\eta$, where ${\sigma }_s$ is the design strength of cable material and $\eta$ is the designated safety factor ($\eta =3$ in this paper).
• Designate initial section area $A_{i0}$ for each cable. Initialize $∆A_i^{\left(k\right)}=0$.
• Get the cable force determined with the cable forces predictor $q_i$.
• Calculate the cable stress ${\sigma }_{id}^{\left(k\right)}$ under dead load according to ${\sigma }_{id}^{\left(k\right)}=q_i/\left(A_{i0}+∆A_i^{\left(k\right)}\right)$.
• Impose load cases, launch analysis and combine the results. Get the cable stress response ${\sigma }_{il}^{\left(k\right)}$ under live load.
• Calculate $∆A_i^{\left(k+1\right)}$ with: $\left(A_{i0}+∆A_i^{\left(k+1\right)}\right)\cdot {\sigma }_{obj}=A_{i0}\cdot \left({\sigma }_{id}^{\left(k\right)}+{\sigma }_{il}^{\left(k\right)}\right)+∆A_i^{\left(k+1\right)}\cdot l_i\cdot \gamma$, where $l_i$ is the length of the cable and $\gamma$ is the volumetric weight of the cable material.
• Designate $A_{i0}+∆A_i^{\left(k+1\right)}$ for each cable area.
• Check convergence of $∆A_i^{\left(k+1\right)}$. If convergence is achieved, perform the last load case analysis; otherwise, return to step 4 and update $∆A_i^{\left(k\right)}$ with $∆A_i^{\left(k+1\right)}$.

The cable areas and stresses are supposed to converge after only a few iterations (two times in our later example). This method benefits the optimization by cutting off the design variables of cable section areas, thus reducing the complexity of the optimization problem.

It is worth pointing out that the presented method causes a slight disturbance to the prediction results due to the redistribution effect of stiffness change. This disturbance is closely connected with the initial section area $A_{i0}$. The closer $A_{i0}$ is to its convergent value $A_{i0}+∆A_i^{\left(k\right)}$, the less the disturbance and the more desirable the corresponding completion state will be.

3.4. The Proposed Optimization Framework

The optimum design of cable-stayed bridges concerning varieties of variables is typically a nonconvex problem. Gradient algorithms can often be trapped in local optima. Therefore, it is advisable to use heuristic algorithms with global search strategies to obtain the global optimum. Among the commonly used heuristic algorithms, Particle Swarm Algorithm (PSO) stands out due to its fast convergence, few parameters, and easy-to-implement formula. It is effective for high-dimensional optimization problems and converges quickly to the optimum solution. Thus, PSO is chosen as the optimization method to solve Formula (3) in Sp2. Based on the above considerations, the PSO-integrated overall procedure, presented in Figure 5, is as follows.

• Problem Definition. Basic information of the optimization problem is identified, including variables, objective, and constraints. In this work, the problem is defined as follows:
  • Objective function—the total theoretical material cost of the towers and the cables;
  • Design variables—several parameters determining the volume or the weight of the towers and the cables;
  • Constraints—variable boundaries to achieve practical design and strength, stiffness and stability verifications;
  • Load cases—dead load, vehicle load, wind load, as well as their combinations according to Chinese Design Code.
• Predictor Construction. The cable forces predictor is constructed with the surrogate model method. Detailed construction steps can be found in Section 3.2.
• Design Optimization. An optimization program is made based on PSO algorithm. Readers can refer to Appendix A for the basic theory of PSO and Appendix B for the detailed steps of the program. The program consists of three modules:
  • Initialization module—particles are generated randomly within the search area. Those that pass the constraint verification remain in the swarm. The Pbest, the Gbest, and their fitness are initialized.
  • Iteration module—velocity and position of each particle are updated in each iteration round. Given that the objective function (cost) can be calculated explicitly and directly by design variables, constraint verification is only performed for particles whose cost is lower than the previously lowest. The improvement leads to a reduction in verification times. If the updated particle passes the constraint verification, the Pbest, the Gbest, and their fitness are updated.
  • Structural analysis module—constraints relating to design code verifications are verified by calling the structural analysis module. In the analysis module, a finite element model is established, the desired completion state is obtained utilizing the cable forces predictor, the responses under various load combinations are analyzed, and finally, the constraints are verified according to the code.

4. Implementation

4.1. Description

A (358 + 818 + 358)-meter-long cable-stayed bridge was chosen as the example to implement our formula and methods, as shown in Figure 6. It was a five-span cable-stayed bridge with two towers. The geometry of the vehicle lanes was planarly straight and vertically circular (with a radius of 29,000 m). The towers were both in a typical H-shape. As for the main girder, it was 38.9 m in width and 3.6 m in height.

4.2. Objective Function

The material cost constitutes a large proportion of the total cost of a construction project. The theoretical material cost of bridge towers and cables was taken as the objective function in this paper. The cost of the beam was excluded, because its cross-section should be designed with priority to satisfy the wind resistance performance. The objective function is stated as

$$\mathsf{\Phi }={\mathsf{\Phi }}_t+{\mathsf{\Phi }}_c={\sum }_{i=1}^{N_t}{C}_t{L}_{ti}{A}_{ti}+{\sum }_{j=1}^{N_c}{C}_c{L}_{cj}{A}_{cj}{\gamma }_c$$

(6)

where $C$ is the unit steel weight or concrete volume cost, $L$ is the length of the member, $A$ is the section area of the member, $\gamma$ is the volumetric weight of material, and $N$ is the number of elements.

It is worth noting that the theoretical cost only considers the raw material consumption, and ignores the other costs including rebars, detailing, labor, and machinery.

4.3. Design Variables

Parameters directly determining the volume or the weight were supposed to be chosen as design variables, because the optimization objective was to reduce the total cost of towers and cables. In the implementation, six variables were chosen for the example in

Table 1. Design variables.

Illustration	Variable	Format	Type	Symbol
	Ratio of tower height to span	Continuous	Geometrical	$X_1$
	Cable anchorage distance on girder	Discrete	Topology	$X_2$
	Longitudinal length at tower bottom	Continuous	Sizing	$X_3$
	Ratio of cross-sectional traverse width to longitudinal length	Continuous	Sizing	$X_4$
	Ratio of cross-sectional thickness to traverse width	Continuous	Sizing	$X_5$
	Cross-sectional size reduction ratio of tower top to tower bottom	Continuous	Sizing	$X_6$

. These variables covered different variable formats and types and were a suitable example for the optimum design problem.

After identifying the design variables, other dependent parameters essential for forming the parametric model were obtained according to the derivations in

Table 2. Determination of dependent design parameters.

Illustration
Parameter	Type	Symbol	Derivation
Length of the side span	Set	$l_s$	$=358 \mathrm{m}$
Length of the main span	Set	$l_m$	$=818 \mathrm{m}$
Distance of standard anchorage	Variable	$l_n$	$=X_2$ in Table 1
Ratio of densification length	Set	${\phi }_e$	$=\frac{l_e}{l_s}=0.65$
Half count of mid-span cables	Dependent	$N$	$=\mathrm{r}\mathrm{o}\mathrm{u}\mathrm{n}\mathrm{d}\left(\frac{l_m}{2\cdot l_n}\right)-1$
Count of symmetric cables	Dependent	$N_s$	$=\mathrm{r}\mathrm{o}\mathrm{u}\mathrm{n}\mathrm{d}\left(\frac{l_s\cdot \left(1-{\phi }_e\right)}{l_n}\right)$
Count of densified cables	Dependent	$N_e$	$=N-N_s$
Length of none-cable area near the towers	Dependent	$l_{d1}$	$=\left(\frac{l_m-l_{d2}}{2l_n}-\mathrm{r}\mathrm{o}\mathrm{u}\mathrm{n}\mathrm{d}\left(\frac{l_m-l_{d2}}{2l_n}\right)+1\right)\cdot l_n$
Length of none-cable area in the middle	Dependent	$l_{d2}$	$=l_n$
Distance of densified anchorage	Dependent	$l_{d3}$	$=\frac{l_s\cdot {\phi }_e}{N_e}$

.

4.4. Load Cases

To cover all the required verifications in the Chinese Design Code, corresponding load and their combination cases were included in the example.

• Dead Load (DL). The dead load consisted of the self-weight, the secondary loads, the cable tension forces, and the counterweight in the side span.
• Vehicle Load (VL). The vehicle loads consisted of uniformly distributed loads (q, 40.2 kN/m) and concentrated loads (P, 1377.1 kN). As simplified static loads, four static load cases with different load layouts were adopted, as shown in Figure 7.
• Wind Load (WL). The wind speed for the bridge site was 30.1 m/s. Two wind loads with different return periods were imposed as the equivalent static gust wind loads on the main girder, the towers, and the cables, according to the Chinese Design Code. The return period of the design wind speed of WL1 was 10 years, and that of WL2 was 100 years.
• Load Combination (LCB). A total of 11 different combinations of the dead load, the vehicle load, and the wind load were considered according to the Chinese Design Code.

4.5. Constraints

Constraints regarding the value ranges of the design variables could be selected based on engineering experience. They were necessary to make sure the results did not deviate too much from practical engineering sizes. Note that the ranges here also defined the boundaries of the samples when constructing the predictor. The consistency ensures the accuracy of the prediction results for cable forces. Meanwhile, constraints regarding strength, stiffness, and stability were defined in accordance with the Chinese Design Code. These constraints chosen to be applied here are the most important items in the preliminary design of cable-stayed bridges. They are the most typical representatives of the strength, stiffness, and stability of the bridge. The constraints are shown in

Table 3. Constraints and boundaries.

Name	Type	Load Cases	Boundary
$X_1$	Variable	All	0.20~0.26
$X_2$	Variable	All	12 m~18 m
$X_3$	Variable	All	13 m~14 m
$X_4$	Variable	All	0.65~0.79
$X_5$	Variable	All	0.15~0.29
$X_6$	Variable	All	0.5~0.9
Concrete Design Strength of the towers	Strength	LCB1~LCB5	−22.4 MPa~1.83 MPa
Steel Design Strength of the main girder	Strength	LCB1~LCB5	−270 MPa~270 MPa
Safety Stress Factor of the cables	Strength	LCB1~LCB5	≥2.5
Vertical Deflection of the main girder	Stiffness	VL1~VL4	≤2.045 m
Elastic Buckling Stability Factor	Stability	LCB6~LCB11	≥4

.

4.6. Cable Forces Predictor

4.6.1. Sampling

A modest number of 972 samples were generated from uniformly discretizing the design variables listed in Table 1, as shown in

Table 4. A total of 972 samples generated from design variables.

Variable	Boundary	Discretization	Detail
$X_1$	0.20~0.26	4	0.20, 0.22, 0.24, 0.26
$X_2$	12 m~18 m	3	12 m, 15 m, 18 m
$X_3$	13 m~14 m	3	13 m, 13.5 m, 14 m
$X_4$	0.65~0.79	3	0.65, 0.72, 0.79
$X_5$	0.15~0.29	3	0.15, 0.22, 0.29
$X_6$	0.5~0.9	3	0.5, 0.7, 0.9

.

Constraints of the cable forces optimization problem mentioned in Section 3.1 included the following: 1000 kN to 10,000 kN cable force boundary for all the cable forces derived from variables $d_1$ to $d_9$; 0 to half side span length for variable $d_{10}$; 0 to the load magnitude when the main girder was filled with steel grit concrete for variable $d_{11}$. After defining the optimization problem, the fmincon function was called in MATLAB to solve it.

4.6.2. Training

The data (i.e., 972 samples including their optimized cable force and counterweight results) were imported into MATLAB and trained with the Statistics and Machine Learning Toolbox. The RMSE and training time of variable $d_1$ are shown in

Table 5. Training model comparison.

Training Model	RMSE	Training Time (s)
PR	1964.4	5.8
RT	253.6	2.9
SVR	761.7	3.4
GPR	168.7	11.6

. Results show that Gaussian Process Regression should be the most suitable training model for our problem, because the RMSE is the least among the four concerned models. Though the training time of GPR is the longest, it was still acceptable in terms of absolute timespan. Polynomial Regression was inadequate in terms of accuracy to deal with optimum design problems, which further highlights the motivation of this work. Regression Tree was both relatively accurate in prediction results and efficient in training time. The fast training advantage of RT is rooted in its efficient binary tree data hierarchy, but it did not show up much because of the modest number of training data. Support Vector Regression was also inaccurate in predicting the results, demonstrated by its relatively large RMSE.

4.6.3. Testing

Selecting a testing design distinguished from our samples ($X_1$ = 0.23, $X_2$ = 15 m, $X_3$ = 13.2 m, $X_4$ = 0.69, $X_5$ = 0.25, $X_6$ = 0.85), a test was performed by comparing its predicted results to individually optimized results. Results in

Table 6. Prediction error of design variables for the testing design.

Design Variables	Predicted Result	Individually Optimized Result	Error (%)
$d_1$	2,690,799.02 N	2,693,537.33 N	−0.102
$d_2$	4,233,299.06 N	4,235,208.77 N	−0.045
$d_3$	6,032,463.85 N	6,032,898.59 N	−0.007
$d_4$	7,723,535.68 N	7,725,499.78 N	−0.025
$d_5$	2,033,779.01 N	2,033,014.80 N	0.038
$d_6$	3,360,536.53 N	3,359,055.61 N	0.044
$d_7$	5,259,640.37 N	5,259,081.00 N	0.011
$d_8$	6,625,880.87 N	6,624,414.47 N	0.022
$d_9$	3,557,614.41 N	3,557,905.04 N	−0.008
$d_{10}$	102.91 m	102.91 m	−0.036
$d_{11}$	162,867.39 N/m	162,925.93 N/m	0.005

and Figure 8 show that the predictor was accurate enough to determine the optimized cable forces and counterweight of our example. The predicted bending moment was close to the optimized one over the entire length. The maximum error of the moment did not exceed 5%, and the relatively large errors were located at the positions where the moment was small in terms of absolute value.

4.7. Optimization Settings

Problem-defining parameter settings for PSO (as shown in

Table 7. Problem-defining parameters settings for PSO algorithm.

Parameter	Symbol	Value
Population	$N_p$	30
Iteration Count	$T_{\mathrm{m}\mathrm{a}\mathrm{x}}$	500
Learning Factor	$c$	$c_1=2.0, c_2=2.0$
Inertia Weight	$\omega$	${\omega }_{\mathrm{m}\mathrm{a}\mathrm{x}}=0.9, {\omega }_{\mathrm{m}\mathrm{i}\mathrm{n}}=0.4$

) were referenced from Cao [30], who dealt with the cable-supported bridge optimization problem with PSO as well.

The results recorded the Gbest (refer to Appendix A) during the 500 iteration rounds. The values of variables were stable after around 320 iterations. Meanwhile, there were a total of 20 price decreases, as shown in Figure 9. After seven times of decrease, the theoretical material cost nearly stopped decreasing. The results indicate the fast convergence speed of PSO algorithm in optimum design problems.

5. Results

The optimized results of the design variables are shown in

Table 8. The optimum design.

Design Variable	Boundary	Initial	Optimized
$X_1$	0.20~0.26	0.253	0.244
$X_2$	12 m, 15 m, 18 m	12 m	18 m
$X_3$	13 m~14 m	13.5 m	13.012 m
$X_4$	0.65~0.79	0.746	0.670
$X_5$	0.15~0.29	0.256	0.150
$X_6$	0.5~0.9	0.804	0.511
Objective Function $\mathsf{\Phi }$ (Relative)	-	1	0.68

. Compared to the initial values, there was a significant decrease (32%) in the theoretical material cost of the optimum design, demonstrating the efficiency of our improved formula.

The internal force distribution of the bridge in the obtained optimum solution was checked to validate the rationality of the design. As shown in Figure 10, the predicted moment of the main girder was in a jagged shape without any sudden changes. Therefore, the completion stage of the optimum design was in the desired state.

Checks were performed to identify if strength, stiffness, and stability constraints were met in the optimum design. As shown in

Table 9. Check of strength, stiffness, and stability constraints.

Constraints	Boundary	Load Cases	Optimum Design
Concrete Design Strength of the towers	−22.4 MPa~1.83 MPa	LCB1~LCB5	−18.6 MPa~1.83 MPa
Steel Design Strength of the main girder	−270 MPa~270 MPa	LCB1~LCB5	−144.51 MPa~52.1 MPa
Safety Stress Factor of the cables	≥2.5	LCB1~LCB5	2.7
Vertical Deflection of the main girder	≤2.045 m	VL	0.953 m
Elastic Buckling Stability Factor	≥4	LCB6~LCB11	7.22

, all constraints were met. Meanwhile, the stress of the towers served as the controlling structural factor in our optimization.

Detailed checks of strength, stiffness, and stability constraints under each load combination are listed in

Table 10. Check of constraints under each load combination.

Type	Verification Item	Load Cases
Strength		LCB1	LCB2	LCB3	LCB4	LCB5
	Tensile stress of the main girder (MPa)	32.26	49.2	52.1	36.45	17.11
	Compressive stress of the main girder (MPa)	−126.73	−129.43	−122.97	−144.51	−118.39
	Tensile stress of the towers (MPa)	0.66	0	1.83	0.42	0.69
	Compressive stress of the towers (MPa)	−18.2	−18.13	−18.65	−18.66	−17.76
	Maximum stress of the cables (MPa)	571.33	619.74	572.6	566.61	526.1
Stiffness		VL1	VL2	VL3	VL4
	Maximum deflection of the main girder (m)	0.611	0.952	0.429	0.551
Stability		LCB6	LCB7	LCB8	LCB9	LCB10	LCB11
	Safe factor	8.13	7.63	7.44	7.79	7.22	8.13

The stresses and deflections were obtained as the maximum over the whole girder or towers. The safe factor of stability was obtained as the elastic buckling stability factor that corresponds to the first modality. Elements of were generated by meshing the line between geometry key points with a discretization number of five for the girder and towers and one for the cables. (In Ansys 2022 R1) The element type of the girder and the towers was Beam4. The element type of the cables was Link10.

. The results show that the stress of towers under LCB3 was close to the boundary (1.83 MPa), which means LCB3 (LCB3 = 1.1 × (1.2 × DL + 1.4 × VL3 + 0.75 × 1.1 × WL1)) served as the controlling load factor of our optimization. This was because the LCB3 was composed of all kinds of load types and the VL3 had a side-span layout, which was unfavorable for the towers.

Finally, the complexity of the optimization was evaluated and compared for the improved formula. In this example, the particle swarm optimization algorithm consisted of 30 particles and 500 iterations. If the routine method was adopted, the number of cable forces optimization would be 15,000 times. By means of a surrogate-model-assisted cable forces predictor, only 972 times were needed for the samples when constructing the predictor before structural optimization. During the structural optimization process, the optimum cable forces and counterweight of the particles were determined with the predictor, eliminating the need for cable forces optimization.

6. Conclusions

In the preliminary design of cable-stayed bridges, it is strategically significant to properly determine the design parameters. However, it is challenging to obtain the optimum design of a cable-stayed bridge because of numerous variables, multiple load cases, and diverse constraints. Integrated methods are adopted in this paper to enable efficient design optimization of cable-stayed bridges:

• to simplify the complexity of the problem, an improved two-layer framework is presented for cable-stayed bridge optimum design problems. The formula consists of two mutually iterative sub-problems: optimizing the internal force distribution by adjusting the cable prestressing forces and optimizing the sizing and geometrical parameters. The sub-problems exhibit fewer variable coupling features, making them easier to solve.
• to decrease the dimension of the design variables, B-spline interpolation curve is adopted to condense the variables, instead of setting all cable forces as variables. B-spline curve stands out when confronting large-span cable-stayed bridges with dense cables, because it can fit cable force distribution in cable-stayed bridges with only a few controlling points or fitting points.
• to improve optimization efficiency, a surrogate model-assisted predictor for optimum cable forces is constructed. The predictor addresses the time-consuming problem of determining the optimum cable prestressing forces in each of the iteration rounds in the optimization problem. This predictor is expected to be the highlight of this paper.
• to deal with the nonconvex optimization problem, an optimization program consisting of the initialization module, iteration module, and structural analysis module is made based on PSO algorithm. PSO is well-known for its global searching ability. The global searching ability is enhanced by the following measures: a moderate population number of 30 which is five times to the number of the design variables; well-set defining parameters of the algorithm; and the dual strategy (refer to Appendix B) adopted when initializing particles.

Finally, the optimum design of a (358 + 818 + 358)-meter-long cable-stayed bridge was chosen as the implementation of the proposed methods. First and foremost, the theoretical material cost was set as the objective function. Six design parameters closely connected with cost of the towers and the cables were identified as design variables. Constraints were set consistently with the Chinese Design Code. Then, 972 samples were generated by uniformly discretizing the design variables. After training and comparison, Gaussian Process Regression was demonstrated as the best surrogate model for prediction. The cable forces predictor was successfully constructed when GPR passed the testing. Afterward, the optimization program was launched to obtain the optimum design of the variables. The cable forces predictor was utilized to determine the optimum cable tensions and counterweight in the program. The predictor eliminated the need to solve the cable forces optimization problem in each of the iteration rounds, resulting in improved efficiency. The results show that the theoretical material cost of the optimum design is 32% lower than the original design. The feasibility and reliability of the structural optimization process were verified by several checks.

Despite our proposed method demonstrated accuracy and efficiency in our example, there are still some possible improvements when confronting larger scale problems. For instance, the surrogate model is naturally a regression trained with given samples, which indicates that if some new design variables are introduced, or if a wider searching range is to be explored, the model needs to be updated for the sake of accuracy. Therefore, it is beneficial to construct a larger version of the cable forces predictor covering a broader range of potential design variables for cable forces optimization problems. In addition, the sampling strategy and different surrogate models can also be more comprehensively compared. These directions can be potentially investigated in future work.