GPU-based parallel programming for FEM analysis in the optimization of steel frames

ABSTRACT

Optimization of large-scale frame structures consumes a vast amount of time since the analysis of such complex systems contains several iterative processes. Mitigating computational burden and reducing this time to a reasonable level is possible by running GPU (Graphical Processing Unit) processors, which can be found on standard computers. This study presents an algorithm for the acceleration of size optimization of steel frames by using the BBO (Biogeography-Based Optimization) method that is suitable for GPU architecture. The GPU-based parallel algorithm, designed for FEM (Finite Element Method) analysis, is applied to three hypothetical steel-frame case structures with different numbers of members and nodes; and processed on four different computers which are available on the market. The presented case studies revealed that the proposed solution’s efficiency increases as the number of members increases and confirmed the ability of the acceleration algorithm for optimization of large-scale frame structures and provided time efficiency.

GRAPHICAL ABSTRACT

KEYWORDS:

• Parallel processing
• GPU-based algorithm
• biogeography-based optimization
• acceleration
• FEM

1. Introduction

Frame construction systems stand as a cornerstone in civil engineering, frequently utilized for various applications. Modern approaches have turned to computer-based methods for optimizing these systems, as evidenced by recent (I. Aydoğdu 2017; Gandomi et al. 2023; Hong, Nguyen, and Nguyen 2022; Kazemzadeh Azad 2020; Kazemzadeh Azad, Hasancebi, and Kazemzadeh Azad 2013; Kociecki and Adeli 2014, 2015; Mikes and Kappos 2023; Vu et al. 2023). Among these methods, metaheuristic techniques have emerged as key players, particularly for large-scale steel frame structures (I. Aydoğdu, Çarbaş, and Akın 2017; Gholizadeh and Mohammadi 2017; Kaveh et al. 2020; Kazemzadeh Azad 2017), thanks to their ability to tackle complex discrete variable problems.

The realm of metaheuristic optimization methods, often referred to as “Nature-inspired optimization”, showcases remarkable flexibility in approaching intricate challenges. Although these techniques have exhibited prowess, it is crucial to acknowledge their inherent drawbacks. Notably, their convergence can be a demanding endeavor, necessitating a substantial number of iterative operations. This characteristic becomes especially prominent when addressing intricate and extensive problems, often involving long simulation times through methods like the Finite Element Method (FEM).

In response to the challenges posed by time-consuming nature-inspired algorithms in optimizing complex systems, researchers have devised a range of strategies tailored to specific aspects of the optimization process. One approach involves the use of surrogate models through FEM (Gholizadeh 2015; Gholizadeh and Milany 2016; Kaveh, Gholipour, and Rahami 2008) and constraint handling (Kaveh, Laknejadi, and Alinejad 2012; Papadrakakis, Lagaros, and Tsompanakis 1998) approximation methods. Surrogate models provide an efficient alternative to exact and computationally demanding finite element analyses. These methods revolve around creating simplified approximations of structural behavior, significantly alleviating the computational burden. These surrogate models can accurately estimate system responses by leveraging a diverse and comprehensive training set. However, it is important to acknowledge that the accuracy of these approximations can vary based on problem complexity and training data quality. Another avenue for improvement is through Improving Optimization Methodology using Hybrid Algorithms. Hybrid algorithms integrate multiple nature-inspired techniques, harnessing their collective strengths (Cheng et al. 2016; Dillen, Lombaert, and Schevenels 2021; Kaveh and Ilchi 2018; Kaveh and Mahdavi 2015; Kaveh, Bakhshpoori, and Afshari 2014; Kazemzadeh Azad 2017; Tayfur, Yilmaz, and Daloglu 2021; Ting et al. 2015). By combining algorithms with complementary characteristics, researchers seek to bolster the robustness and efficiency of the optimization process. It is essential to recognize that the performance of these hybrid algorithms can be problem-dependent, meaning their efficacy may vary based on the specific optimization challenge at hand. A distinct category is comprised of Incorporating Design-Driven Methods. These techniques tailor the optimization process to the unique characteristics of the problem, aiming to reduce the maximum number of iterations required. Drawing from domain-specific knowledge, design-driven methods guide the algorithm toward optimal solutions more efficiently. While highly effective for their intended problems, these methods are confined by their problem-dependent nature, thus limiting their applicability across a spectrum of optimization scenarios. A particularly promising strategy entails Parallel or Distributed Computing with Additional Programming Challenges. Unlike other strategies, this method taps into the power of parallel or distributed computing to significantly expedite the optimization process. By executing tasks concurrently across multiple processors or machines, this approach yields substantial time savings. Remarkably, parallel and distributed computing does not necessarily necessitate methodological alterations. It capitalizes on hardware architecture to partition intricate tasks into manageable components that can be processed simultaneously. This versatility renders it applicable to problems of varying structures, rendering it a compelling solution for addressing time-intensive optimization challenges.

Parallel computing, as traditionally depicted in the literature, has largely been synonymous with deploying clusters of interconnected computers or CPUs to simultaneously analyze complex structures (Hasançebi et al. 2011; Martínez-Frutos and Herrero-Pérez 2017; Papadrakakis, Lagaros, and Fragakis 2003; Sarma and Adeli 2001; Zegard and Paulino 2013; Zhu 2010). While this approach presents considerable potential, it invariably requires significant financial investments to establish and maintain the necessary hardware infrastructure. However, a transformative progression has emerged due to recent advancements in GPU (Graphics Processing Unit) processors. This evolution enables the realization of parallel programming by harnessing the capabilities of standard computer graphics cards. These readily available components empower researchers and engineers to tap into parallel computing without the need for specialized and expensive hardware setups. Although the GPU method’s applicability has been tested in various fields in recent years, ranging from mechanical engineering to medical sciences (Couty et al. 2021; Lei et al. 2019; Luo, Ye, and Chen 2020; Shkurti et al. 2013; Vargas et al. 2022; ZendehAli, Emdad, and Abouali 2023), the studies in the field of structural analysis and structural optimization are limited (Georgescu, Chow, and Okuda 2013; Träff et al. 2023; Xu et al. 2014; Zhu 2010). Therefore, the contribution of GPU-based parallel programming to structural optimization problems stands as a compelling research question.

This study aims to unveil the inherent advantages of employing a GPU-based approach to reduce computation time during the optimization of steel frame systems. To accomplish this goal, a GPU-based parallel algorithm has been meticulously crafted utilizing CUDA (Compute Unified Device Architecture). The Biogeography-Based Optimization (BBO) method is selected for its effectiveness in frame optimization, and its widespread use in Finite Element Method (FEM) analysis of frame structures, making it a popular choice in the literature (Artar and Çarbaş 2021; I. Aydoğdu 2017; Çarbaş 2016; Ghatte 2021; Guo et al. 2017; Hosseini, Ghasemi, and Dizangian 2022). This method is supported by numerous theses and articles in the literature, many of which feature benchmark problems, providing ample opportunities to assess the impact of GPU acceleration on computational performance. In this study, examples from these benchmarks were utilized, ensuring consistency with established practices regarding parameters such as the required number of iterations and load cases/types (Akın and Aydoğdu 2015; Aydoğdu 2010; Çarbaş 2016; Örmecioğlu 2019). The examples considered wind load as the dynamic load, aligning with previous literature. In a bid to assess the program’s efficacy, a version devoid of parallel processing, aligned with the CPU (Central Processing Unit) architecture, was developed for comparative purposes.

Both versions are applied to two hypothetical steel-frame case structures with increasing numbers of members and nodes. A comparison between the two versions is conducted using various load combinations. The constraint functions of the design criteria consist of strength constraints, horizontal and inter-story displacement constraints, as well as column-to-column and column-to-beam constraints of the “Load and Resistance Factor Design Specifications for Steel Buildings” by American Institute of Steel Construction (LRFD-ASIC 2000; AISC, 1986). For the research, a Geforce GTX 1080 titan processor (which consists of 3584 CUDA cores and 28 SMs) was adopted as the main hardware platform. In addition, the algorithm was tested on three other computer configurations for comparison. The results indicate that GPU processors are more successful in FEM operations than CPU processors and verified the GPU-based parallel algorithm’s efficiency in the optimization of large-scale frame structures (Örmecioğlu 2019).

2. Methodology

2.1. Frame optimization

Frame optimization is a process used in structural engineering, aiming to improve structural efficiency and overall stability most efficiently. The primary objectives of this optimization process are minimizing material usage, reducing costs, enhancing structural performance, and meeting safety and building code requirements. One of the fundamental attributes of frame optimization is its systematic approach to enhancing structural frame designs. This process commences with the initial design and Finite Element Method (FEM) analysis. It involves defining objectives and constraints, adjusting design variables, and applying optimization algorithms. In a broader context, optimization problems can generally be delineated by three primary components: objective functions, design variables, and constraints.

For most of the frame optimizations, the mathematical model of the design optimization problem is designed to minimize the overall weight of the frame system which directly affects the material cost. Hence, the objective function of frame optimization can be formulated as;

(1) $\mathit{Minimum}\mathit{\hspace{0.17em}}\mathit{W}\left(\mathit{x}\right)=\sum _{\mathit{r}=1}^{\mathit{ng}}\left(m_r\sum _{\mathit{s}=1}^{t_r}{l}_s\right)\mathit{\hspace{0.17em}}$(1)

where $\mathit{W}$ represents the total weight of the frame system and $m_r$ is the weight per unit length of the steel profile assigned to group $\mathit{r}$. While $t_r$ indicates the total number of elements within group $\mathit{r}$, $\mathit{ng}$ represents the total number of groups in the entire frame system and $l_s$ represents the length of member s. $\mathit{x}\left(\mathit{x}=\left[x_1,x_2,\cdots x_{\mathit{ng}}\right]\right)$ is the vector of integer values representing the sequence numbers of steel sections assigned to member groups, which are referred to as design variables.

Next, the constraint functions are defined. In this study constraint functions are defined as (1) Strength constraints, (2) Displacement constraints, and (3) Geometric constraints based on LRFD-ASIC standards.

Strength constraints for the beam-column members of the structures (LRFD-ASIC 2000) are;

(2) $g_{\mathit{s},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{P_u}{\mathit{\varphi }{P}_n}\right)}_{\mathit{il}}+\frac{8}{9}{\left(\frac{{\mathit{M}}_{\mathit{ux}}}{{\varphi }_b{M}_{\mathit{nx}}}+\frac{{\mathit{M}}_{\mathit{uy}}}{{\varphi }_b{M}_{\mathit{ny}}}\right)}_{\mathit{i},\mathit{l}}-1.0\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\text{breakfor}\mathit{\hspace{0.17em}}\frac{P_u}{\mathrm{\varnothing }{P}_n}\ge 0.2\mathit{\hspace{0.17em}}$(2)

(3) $g_{\mathit{s},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{P_u}{2\mathit{\varphi }{P}_n}\right)}_{\mathit{il}}+{\left(\frac{{\mathit{M}}_{\mathit{ux}}}{{\varphi }_b{M}_{\mathit{nx}}}+\frac{{\mathit{M}}_{\mathit{uy}}}{{\varphi }_b{M}_{\mathit{ny}}}\right)}_{\mathit{i},\mathit{l}}-1.0\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\text{breakfor}\mathit{\hspace{0.17em}}\frac{P_u}{\mathrm{\varnothing }{P}_n}<0.2\mathit{\hspace{0.17em}}$(3)

$\mathit{i}=1,\mathit{\hspace{0.17em}}2,\mathit{\hspace{0.17em}}\dots ,\mathit{nm}$

where, the symbol $M_{\mathit{nx}}$ represents the nominal flexural strength along the strong axis (x-axis), $M_{\mathit{ny}}$represents the nominal flexural strength along the weak axis (y-axis), $M_{\mathit{uy}}$represents the required flexural strength on the weak axis, $M_{\mathit{ux}}$ represents the required flexural strength on the strong axis, $P_n$ represents the nominal axial strength (tension and compression), $P_u$ represents the required axial strength (tension and compression), $\mathit{i}$ stands for the member l denotes the load case, and $\mathit{nm}$ is the total number of members in the frame.

The formula for deformation constraints is;

(4) $g_{\mathit{d},\mathit{j}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}\frac{{\mathit{\delta }}_{\mathit{jl}}}{{\delta }_j^u}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{j}=1,\dots ,n_{\mathit{sm}}\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{l}=1,\dots ,\mathit{\hspace{0.17em}}{n}_{\mathit{lc}}$(4)

where, ${\delta }_{\mathit{jl}}$ represents the maximum deflection of element number j under load case-l, ${\delta }_j^u$ represents the upper limit of deflection that the steel member can endure, $n_{\mathit{sm}}$ represents all elements to which the deflection limit is applied, and $n_{\mathit{lc}}$ represents the number of load cases.

Also, the top-story displacement and inter-story displacement constraints are formulated respectively as;

(5) $g_{\mathit{td},\mathit{j}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}\frac{{\left({\mathrm{\Delta }}_{\mathit{top}}\right)}_{\mathit{jl}}}{\mathit{Limit}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{j}=1,\dots ,n_{\mathit{jtop}}\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{l}=1,\dots ,\mathit{\hspace{0.17em}}{n}_k$(5)

and,

(6) $g_{\mathit{id},\mathit{j}}=\mathit{\hspace{0.17em}}\frac{{\left({\mathrm{\Delta }}_{\mathit{oh}}\right)}_{\mathit{jl}}}{\mathit{Limit}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{j}=1,\dots ,n_{\mathit{st}}\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{l}=1,\dots ,\mathit{\hspace{0.17em}}{n}_{\mathit{lc}}$(6)

where, $n_{\mathit{jtop}}$ represents the number of node points on the top floor, and ${\left({\mathrm{\Delta }}_{\mathit{top}}\right)}_{\mathit{jl}}$ represents the displacement condition of node j on the top floor with load case-l. In Eq. 6(6) $g_{\mathit{id},\mathit{j}}=\mathit{\hspace{0.17em}}\frac{{\left({\mathrm{\Delta }}_{\mathit{oh}}\right)}_{\mathit{jl}}}{\mathit{Limit}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{j}=1,\dots ,n_{\mathit{st}}\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{l}=1,\dots ,\mathit{\hspace{0.17em}}{n}_{\mathit{lc}}$(6) $n_{\mathit{st}}$ denotes the number of floors, $n_{\mathit{lc}}$ denotes the number of load conditions, and ${\left({\mathrm{\Delta }}_{\mathit{oh}}\right)}_{\mathit{jl}}$ denotes the value of displacement at floor number j, in load case-l. In Eqs. (2)(2) $g_{\mathit{s},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{P_u}{\mathit{\varphi }{P}_n}\right)}_{\mathit{il}}+\frac{8}{9}{\left(\frac{{\mathit{M}}_{\mathit{ux}}}{{\varphi }_b{M}_{\mathit{nx}}}+\frac{{\mathit{M}}_{\mathit{uy}}}{{\varphi }_b{M}_{\mathit{ny}}}\right)}_{\mathit{i},\mathit{l}}-1.0\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\text{breakfor}\mathit{\hspace{0.17em}}\frac{P_u}{\mathrm{\varnothing }{P}_n}\ge 0.2\mathit{\hspace{0.17em}}$(2) -(6(6) $g_{\mathit{id},\mathit{j}}=\mathit{\hspace{0.17em}}\frac{{\left({\mathrm{\Delta }}_{\mathit{oh}}\right)}_{\mathit{jl}}}{\mathit{Limit}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{j}=1,\dots ,n_{\mathit{st}}\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{l}=1,\dots ,\mathit{\hspace{0.17em}}{n}_{\mathit{lc}}$(6) ), the allowable deflection, displacement, and drift values are computed per the ASCE Ad Hoc Committee report (Ellingwood 1986)

Finally, geometric constraints for the column-to-column (CtoC) and beam-to-column (B to C) connections of the frame structures are formulated as;

(7) $g_{\mathit{cd},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{d_u}{d_l}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{cc}}\mathit{\hspace{0.17em}}$(7)

(8) $g_{\mathit{cm},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{m_u}{m_l}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{cc}}\mathit{\hspace{0.17em}}$(8)

(9) $g_{\mathit{bc},\mathit{i}}\left(\mathit{x}\right)\mathit{\hspace{0.17em}}=\mathit{\hspace{0.17em}}{\left(\frac{b{f}_b}{d_c-2\mathit{t}{f}_c}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{bc}}\mathit{\hspace{0.17em}}$(9)

(10) $g_{\mathit{bb},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{b{f}_b}{\mathit{b}{f}_c}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{bc}}\mathit{\hspace{0.17em}}$(10)

In Eq. (7)(7) $g_{\mathit{cd},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{d_u}{d_l}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{cc}}\mathit{\hspace{0.17em}}$(7) and (8(8) $g_{\mathit{cm},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{m_u}{m_l}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{cc}}\mathit{\hspace{0.17em}}$(8) ), subscripts $\mathit{i}$, $\mathit{u}$, and $\mathit{l}$ respectively donate the CtoC connection ID, upper column ID on connection $\mathit{i}$, and upper column ID on connection $\mathit{i}$. Here, $\mathit{d}$ represents the depth of the column, while $\mathit{m}$ indicates the weight of the column per unit length. The variable $n_{\mathit{cc}}$ signifies the number of CtoC connections considered in the optimization problem.

In Eq. (9)(9) $g_{\mathit{bc},\mathit{i}}\left(\mathit{x}\right)\mathit{\hspace{0.17em}}=\mathit{\hspace{0.17em}}{\left(\frac{b{f}_b}{d_c-2\mathit{t}{f}_c}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{bc}}\mathit{\hspace{0.17em}}$(9) and (10(10) $g_{\mathit{bb},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{b{f}_b}{\mathit{b}{f}_c}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{bc}}\mathit{\hspace{0.17em}}$(10) ), the subscripts $\mathit{i}$, $\mathit{b}$, and $\mathit{c}$ respectively stand for the BtoC connection ID, beam ID on connection $\mathit{i}$, and column ID on connection $\mathit{i}$. Meanwhile, $n_{\mathit{bc}}$ represents the number of BtoC connections within the frame. Additionally, $\mathit{bf}$ and denote the flange width and flange thickness of the section. An illustration of the BtoC connection is given in Figure 1.

Figure 1. The beam-column connection in a steel frame.

Eq. (9)(9) $g_{\mathit{bc},\mathit{i}}\left(\mathit{x}\right)\mathit{\hspace{0.17em}}=\mathit{\hspace{0.17em}}{\left(\frac{b{f}_b}{d_c-2\mathit{t}{f}_c}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{bc}}\mathit{\hspace{0.17em}}$(9) is applied when the beam is connected depth of the column, while Eq. (10)(10) $g_{\mathit{bb},\mathit{i}}\left(\mathit{x}\right)=\mathit{\hspace{0.17em}}{\left(\frac{b{f}_b}{\mathit{b}{f}_c}\right)}_{\mathit{i}}-1\mathit{\hspace{0.17em}}\le 0\mathit{\hspace{0.17em}},\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,n_{\mathit{bc}}\mathit{\hspace{0.17em}}$(10) applies when the beam is connected flange of the column. For more comprehensive information, please refer to the following reference studies: Aydoğdu (2010) and Örmecioğlu (2019).

The values derived from these constraint functions are then converted to a single value within the framework of the structure’s weight using a cost function ($f_{\mathit{cost}}$). In doing so, a discriminative equation (Eq.11(11) ${\mathit{f}}_{\mathit{cost}}=\mathit{W}\left(\mathit{x}\right){\left[1+\sum {\mathit{c}}_{\mathit{i}}\right]}^{\mathit{P}};{\mathit{C}}_{\mathit{i}}=\left\{\begin{array}{c}0\mathit{for}\mathit{\hspace{0.17em}}{\mathit{g}}_{\mathit{i}}\left(\mathit{x}\right)\le 0\\ {\mathit{g}}_{\mathit{i}}\left(\mathit{x}\right)\mathit{for}\mathit{\hspace{0.17em}}{\mathit{g}}_{\mathit{i}}\left(\mathit{x}\right)>0\end{array}\right\}\mathit{i}=1,2,\cdots ,\mathit{NC}$(11) ) is used to distinguish successful structures from failed ones, highlighting the difference.

(11) ${\mathit{f}}_{\mathit{cost}}=\mathit{W}\left(\mathit{x}\right){\left[1+\sum {\mathit{c}}_{\mathit{i}}\right]}^{\mathit{P}};{\mathit{C}}_{\mathit{i}}=\left\{\begin{array}{c}0\mathit{for}\mathit{\hspace{0.17em}}{\mathit{g}}_{\mathit{i}}\left(\mathit{x}\right)\le 0\\ {\mathit{g}}_{\mathit{i}}\left(\mathit{x}\right)\mathit{for}\mathit{\hspace{0.17em}}{\mathit{g}}_{\mathit{i}}\left(\mathit{x}\right)>0\end{array}\right\}\mathit{i}=1,2,\cdots ,\mathit{NC}$(11)

Here, $\mathit{P}$ is the violation factor recommended as 2, $\mathit{NC}$ is the number of contains in the optimization problem. This ensures that steel frame structures that cannot fulfill the constraints are given significantly higher values than those that pass, placing them at the lower end of the population set.

2.2. GPU (graphical processing unit) based parallel processing and CUDA (compute unified device architecture)

To optimize frame structures, civil engineers utilize CAD software, specialized structural analysis tools, and optimization algorithms driven by recent advancements in computer technology (Meon et al. 2012). These tools enable modeling and analysis of various design alternatives, considering factors like load types, site conditions, materials, and project-specific parameters, aiming for a safe, efficient, and cost-effective structural solution. Finite element analysis (FEM) allows for detailed simulations of frame behavior under different loads, aiding in the evaluation of design changes. However, the Gaussian elimination phase in FEM is time-consuming, especially with the increasing number of structural members and stiffness matrices in complex frame structures, exceeding the capacity of existing computers.

This study focuses on accelerating the Gaussian elimination process, crucial in FEM analysis, by utilizing GPU architecture due to its ability to handle parallel tasks (Duff et al. 1988). The phase, which solves the system stiffness matrix, becomes more time-consuming as the number of matrices increases in large-scale frame structures. The iterative nature of this operation in structure optimization requires repeated FEM analysis for each emerging cluster, further extending the computational time. Metaheuristic optimization, such as the Biogeography-Based Optimization (BBO) method, is then applied to approach the optimal solution (Onan 2013). This computational intensity makes GPU architecture highly suitable for optimizing steel frame systems.

GPU (Graphical Processing Unit) processors were originally designed specifically for the computation of a serial sequence of repetitive graphical operations called GP (Graphics Pipeline). However, over the years, their design has evolved towards more flexible programmability, making them suitable for very intensive repetitive and data-intensive computations called “High-Performance Computing” (HPC) (Figure 2). Especially the release of CUDA (Compute Unified Device Architecture) by NVIDIA in 2007 became a milestone in this regard (McClanahan 2011; Örmecioğlu 2019).

Figure 2. Graphic showing the performance of GPU processors in numerical computations over the years (2000–2023).

In GPU parallel programming and the CUDA framework, “kernels” and “threads” are foundational concepts. A kernel is essentially a self-contained, parallelizable function written in CUDA C/C++ that is executed on the GPU. It serves as a unit of work that can be processed concurrently by multiple threads. Kernels are designed to perform specific tasks or computations on GPU data, often operating on different portions of the data in parallel. When a kernel is launched, it runs in parallel across multiple threads on the GPU, with each thread executing the same code but typically working on distinct data elements or performing slightly different computations.

Threads, on the other hand, refer to individual execution units within a kernel (Figure 3). Threads are the smallest units of work on the GPU and are organized into a hierarchical structure: threads within a block, and blocks within a grid. Threads within a block can communicate and synchronize with each other. This parallelism at the thread level is what enables the GPU to achieve high-performance parallel processing. Modern GPUs can execute many threads simultaneously, making them well-suited for computationally intensive tasks such as scientific simulations, deep learning, and rendering.

Figure 3. Parallel thread hierarchy in CUDA.

Nevertheless, typically, conventional optimization algorithms are often unsuitable for GPU-based systems due to their design, which relies on the sequential nature of each cycle affecting the next. Considering such chains of interactions, it is, therefore, necessary to develop a specialized GPU-based approach that focuses on the most time-consuming stage in structural metaheuristic optimization algorithms: the structural analysis phase. Within the context of this article, the methods to accelerate this pivotal phase by utilizing the parallel processing capabilities of GPUs have been examined.

2.3. The BBO algorithm

The Biogeography-Based Optimization method is a nature-inspired optimization algorithm based on the principles of biogeography, which delves into the geographical distribution of biological organisms. In BBO, potential solutions to an optimization problem are represented as habitats, and the algorithm simulates the migration of species (solution components) between these habitats based on their suitability values. The goal is to improve the suitability of the solutions over iterations by allowing them to adapt and evolve in a way like how species evolve in real ecosystems. In this optimization technique where individuals, affected by regional influences, strive to achieve optimal outcomes within the context of their movements under these influences. The algorithm mimics the process of species migration and adaptation across diverse habitats to find optimal solutions to optimization problems. In this method, individuals tend to either migrate or remain in their current regions based on their Habitat Suitability Index (HSI), while also exhibiting a simultaneous tendency to migrate to other regions under the guidance of external dispersal influences due to Suitability Index Variable (SIV) by undergoing mutations (Kucukkulahli, Erdogmus, and Polat 2017; Simon 2008).

Through these movements, rather than concentrating solely on a single solution point, the BBO approach shifts towards focusing on different solution regions guided by specified external influences. In BBO, these external influences represent diverse environmental factors that impact the suitability of different solution regions. Instead of fixating on one point, individuals (solutions) explore multiple regions, akin to species exploring various habitats. Certain regions might offer more favorable conditions for optimization, while others might present challenges. This broader exploration diversifies the optimization search, increasing the likelihood of finding a more optimal solution. By responding to external cues, individuals are prompted to explore new regions, analogous to how species migrate in response to environmental changes. This dynamic exploration enables a better understanding of the optimization landscape and helps avoid getting trapped in local optima. In summary, the essence of BBO lies in its ability to adaptively explore multiple solution regions guided by external influences. This approach fosters the discovery of more optimal solutions by avoiding fixation on a single point and allowing the optimization process to traverse diverse solution spaces.

The mathematical modeling of the BBO process comprises two primary phases: migration and mutation. During the migration phase, individuals within a habitat move from one habitat to another within a specific population set. During their movement, each individual within the population changes as a result of external influences (design variable). As a result of this process, the migrating population is referred to as the changing solution, and the habitat receiving the migration is called the reference solution. The reference solution is selected from among the solutions within the population using the roulette method. The probability of each habitat (solution) being the reference solution depends on the probability ratio of the habitat receiving the migration, calculated by the formula below:

(12) $\mathit{P}\left(x_j\right)=\frac{{\mu }_j}{{\sum }_{i=1}^N{\mu }_i};\mathit{\hspace{0.17em}}\mathit{j}=1,\dots ,\mathit{\hspace{0.17em}}\mathit{PS}$(12)

where N represents the number of organisms in the habitat, ${\mathit{\mu }}_{\mathrm{j}}$ represents the migration coefficient for habitat, and PS represents the population size.

During the mutation phase, individuals within the migrating population undergo mutations due to external influences and are then added to other habitats. In the BBO algorithm, the mutation process occurs if the generated random number falls below a previously defined mutation rate. If mutation occurs, the design variable that underwent mutation during the migration phase is randomly determined.

(13) $x_i=x_{\mathit{li}}+\mathit{rand}\left(0,1\right)\left(x_{\mathit{ui}}-x_{\mathit{li}}\right);\mathit{\hspace{0.17em}}\mathit{i}=1,\dots ,\mathit{\hspace{0.17em}}\mathit{PS}$(13)

Where, ${\mathrm{x}}_{\mathrm{i}}$ represents the design variable of the solution undergoing mutation, ${\mathrm{x}}_{\mathrm{ui}}$ and ${\mathrm{x}}_{\mathrm{li}}$ represent the upper and lower design limits of the variable, respectively, and PS represents the population size. “rand(0,1)” is a function that assigns a random value between 0 and 1.

In the context of the study, the FEM analysis of steel frame systems detailed above has been integrated with BBO optimization, as visualized by algorithmic flow steps in Figure 4. Subsequently, the cyclical nature of the process has been fundamentally streamlined into a series of sequential phases by using Algorithm 1. Beginning with the “sorting of the population” to the “selection of the best individual”, this approach has enabled the entire process to be compartmentalized into distinct phases, allowing the identification of the smallest processing steps. These include vital BBO procedures like Migration and mutation, as well as Structural Analysis and Evaluation phases. Among these, selected phases have been reconfigured to run exclusively on the GPU processor, based on their compatibility with parallel programming paradigms. GPU-based Parallel Algorithm for FEM is detailed in section 2.4.

Figure 4. Programmatic diagram visualizing the algorithmic flow specialized for BBO optimization of steel frame systems. Grey indicated the GPU-accelerated structural analysis phase.

Algorithm 1. The bbo optization algorithm.

2.4. GPU-based parallel algorithm for FEM

2.4.1. FEM method

The basic idea of the FEM method is dividing the complex structure into a finite number of interconnected members, to determine the loads acting on each node, to calculate the displacements in the direction of these loads, and thus to reach a result from unit to whole. The joint displacements and end forces are defined in the local axis of a member is visually depicted in Figure 5.

Figure 5. The forces, moments, and displacements that occur at two adjacent nodes of a frame.

By using this method, differential equations are formulated to calculate the end forces and displacements of each discrete member which will be represented by its stiffness matrix, numerically representing the member’s behavior. Since there are 6-degree of freedom defined for each node of the space frame element, the size of the [$k_i$] is 12 × 12.

The member stiffness matrix determines the values that affect the cross-sectional forces for each member in the local coordinate system, which needs to be converted to a global coordinate using the following formula.

(14) $\left[K_i\right]={\left[B_i\right]}^{\mathit{T}}\left[k_i\right]\left[B_i\right];\mathit{i}=1,2,\dots ,\mathit{m}$(14)

Here, $\left[K_i\right]$ is the stiffness matrix of member i for global coordinate, $\left[B_i\right]$ is the transformation matrix, and $\mathit{m}$ is the number of members defined in the structure.

After coordinate transformation, the system stiffness matrix ($K_{\mathit{sys}}$) that enables an understanding of how the system’s overall behavior is affected by the members’ displacements is computed by using Eq. 15

(15) $\left[K_{\mathit{sys}}\right]=\sum _{\mathit{i}=1}^m\left[K_i\right]$(15)

After obtaining ($K_{\mathit{sys}}$), the following Eq. 16(16) $\left\{P_{\mathit{sys}}\right\}=\mathit{\hspace{0.17em}}\left[K_{\mathit{sys}}\right]\left\{d_{\mathit{sys}}\right\}$(16) is derived by transferring the external loads ($P_{\mathit{sys}}$) acting on the structural system to find the system displacement vector $\left\{d_{\mathit{sys}}\right\}$ which provides the displacements and rotations for each degree of freedom in the structural system.

(16) $\left\{P_{\mathit{sys}}\right\}=\mathit{\hspace{0.17em}}\left[K_{\mathit{sys}}\right]\left\{d_{\mathit{sys}}\right\}$(16)

The following matrix, resulting from Eq. 16(16) $\left\{P_{\mathit{sys}}\right\}=\mathit{\hspace{0.17em}}\left[K_{\mathit{sys}}\right]\left\{d_{\mathit{sys}}\right\}$(16) is a complex system of equations, hence, linear numerical analysis techniques such as the Gaussian elimination method are employed in solving it and obtaining the system displacement vector $\left\{d_{\mathit{sys}}\right\}$.

Consequently, utilizing the value of the system displacement vector $\left\{d_{\mathit{sys}}\right\}$, the internal forces within the members and section requirements in the system are determined.

After computation of $\left\{d_{\mathit{sys}}\right\},$ the member forces and moments ($\left\{f_i\right\}$) are computed using the following formula;

(17) $\left\{f_i\right\}=\left[k_i\right]\left[B_i\right]\left\{u_i\right\};\mathit{i}=1,2,\dots ,\mathit{m}$(17)

where $\left\{u_i\right\}$ is the global displacement vector of member i.

2.4.2. Acceleration of FEM for GPU

Within the framework of the BBO algorithm’s sequential execution, the Structural Analysis stage, distinguished by its intensive time demands, is broken down into the steps shown in Figure 6, with parallel processing on the GPU enabled by Algorithm 2. The phase’s initial task, “Creation of $K_{\mathit{sys}}$”, is accomplished following the steps outlined in Eq. 15(15) $\left[K_{\mathit{sys}}\right]=\sum _{\mathit{i}=1}^m\left[K_i\right]$(15) , through a process of triplet matrix multiplication that integrates the outcomes via the interconnected structural elements. The procedural formulations for each component in Equation 14(14) $\left[K_i\right]={\left[B_i\right]}^{\mathit{T}}\left[k_i\right]\left[B_i\right];\mathit{i}=1,2,\dots ,\mathit{m}$(14) are extracted, arranged into a matrix with dimensions of 12 × 12, and made operational through CUDA. Following this, the CUDA function allows for the parallel formation of the $K_i$ result for each component, thereby completing the requisite operations for Eq. 15(15) $\left[K_{\mathit{sys}}\right]=\sum _{\mathit{i}=1}^m\left[K_i\right]$(15) in a parallel manner.

Figure 6. Programmatic diagram visualizing the algorithmic flow of the structural analysis module.

Algorithm 2. Structural optimization module.

The next step in the analysis is determining the “U displacement Vector.” This vector provides a comprehensive view of the displacements at the structure’s nodes, illustrating the three-dimensional movements at each node’s coordinates, which is indispensable for evaluating the system’s behavior under different conditions. The K_sys matrix, derived from Eq. 15(15) $\left[K_{\mathit{sys}}\right]=\sum _{\mathit{i}=1}^m\left[K_i\right]$(15) , is applied in Eq. 16(16) $\left\{P_{\mathit{sys}}\right\}=\mathit{\hspace{0.17em}}\left[K_{\mathit{sys}}\right]\left\{d_{\mathit{sys}}\right\}$(16) . At this juncture, a parallel Gauss elimination process is conducted using the Parallel CuBLAS library functions, which are optimized for the GPU processor. This process is carried out distinctly for each structural loading condition, leading to the accurate calculation of the $d_{\mathit{sys}}$ displacement vector. In the parallel Gauss elimination process, the inherent advantages of stiffness matrices always being sparse matrices have been effectively utilized. In the context of structural FEM analysis, the stiffness matrices consistently take the form of a sparse matrix and some non-zero elements adhere to a diagonal arrangement. As a result, these matrices exhibit symmetry along this diagonal axis, thereby providing a computational advantage through this discernible pattern.

In the ultimate phase, the “Creation of the f vector,” the calculations for loads and moments on each member are performed. The operations for each element outlined in Eq. 17(17) $\left\{f_i\right\}=\left[k_i\right]\left[B_i\right]\left\{u_i\right\};\mathit{i}=1,2,\dots ,\mathit{m}$(17) have been systematically derived and assembled into a 12 × 1vector, which is then operationalized with CUDA. Through this CUDA function, the parallel computation of the $f_i$ result for each element is achieved, ensuring that the function is executed in parallel for all structural elements, significantly improving the computational process and the precision of the structural analysis.

All these features have been seamlessly integrated into the design of the algorithm developed for this purpose. Subsequently, formulas were developed on CUDA (Compute Unified Device Architecture) for all loads acting in each direction, and the algorithm was customized to effectively harness the processing power of GPUs through CUDA integration. Thus, the developed formulas were assigned to GPU processor threads for all nodes, resulting in a parallel process that gains processing speed.

3. Design examples and test computers

In this study, three distinct hypothetical structures, each with varying numbers of elements, were implemented under three progressively increasing loading scenarios. The initial structure comprises 105 members, followed by the second structure with 460 and the third with 1024 members. The load cases for each structural category are quantified at 1, 4, and 7, respectively. All design examples used were selected from benchmark problems established in the literature, certifying consistency with established examples regarding parameters such as load types and load cases (Akın and Aydoğdu 2015; İ. Aydoğdu and Saka 2012; Çarbaş 2016; Örmecioğlu 2019). Hence, the whole load combinations of the study for both design examples are defined as; 1.2DL + 1.6LL + 0.5SL, 1.2DL + 0.5LL + 1.6SL, 1.2DL + 1.6WXL+LL + 0.5SL, 1.2DL + 1.6WYL+LL + 0.5SL, 1.2DL + 0.8WXL +1.6SL, 1.2DL + 0.8WYL +1.6SL and 1.4DL where DL, LL, SL and WL represents dead load, live load, snow load and wind load respectively. The design loads, basic wind speed, drift, and deflection limits of the frames are based on ASCE 7–05 (ASCE 7–05 2005) and Ad Hoc Committee on Serviceability, and are given in Table 1.

Table 1. Comparison of loading details and displacement limitations for the first and second design examples.

Type of load	105-member frame	460-member frame	1024-member frame
Live load	$2.39\mathit{kN}/m^2$	$2.39\mathit{kN}/m^2$	$2.39\mathit{kN}/m^2$
Dead load	$2.88\mathit{kN}/m^2$	$2.88\mathit{kN}/m^2$	$2.88\mathit{kN}/m^2$
Snow load	$0.755\mathit{kN}/m^2$	$0.755\mathit{kN}/m^2$	$0.755\mathit{kN}/m^2$
Wind speed	$65.00\mathit{m}/s^{\mathit{\hspace{0.17em}}}$	$38.00\mathit{m}/s^{\mathit{\hspace{0.17em}}}$	$38.00\mathit{m}/s^{\mathit{\hspace{0.17em}}}$
Top story Dep.	$6.67\mathit{cm}$	$3.50\mathit{cm}$	$3.50\mathit{cm}$
Inter-story Dep.	$1.33\mathit{cm}$	$0.875\mathit{cm}$	$0.875\mathit{cm}$

This study aimed to prepare an algorithm that would facilitate the optimization analyses of steel frame structures using commonly available hardware in the market, rather than specialized hardware. Using the algorithm developed for this purpose, the study attempted to evaluate the suitability of different hardware configurations commonly found in the market for structural optimization tasks. Therefore, hardware selection was based on market availability and categorized by age (old, medium, new). Subsequently, the algorithm developed to accelerate the structural analysis phase of the BBO process was executed independently on four distinct test computers. These test computers were selected to represent a range of hardware ages and capabilities, including both desktop and laptop configurations. The goal was to evaluate how different hardware setups influenced the optimization process, with a focus on achieving efficient and effective structural analyses for steel frame structures.

Test computer-A, Geforce GTX 1080 titan processor is a standard desktop graphics processor that includes features such as video decoding, OpenGL, Vulkan, and DirectX graphics processing pipelines, and it consists of a total of 3584 CUDA cores and 28 streaming multiprocessors (SMs). Each SM contains 128 CUDA cores, and within a single SM, all cores execute the same commands, forming an integrated parallel structure known as the MP (Massively Parallel) type. Alongside this processor, a desktop computer with an i7-8700K CPU, featuring 6 cores, served as the primary hardware platform (main test computer) for the study, while the other test computers were randomly selected for comparison.

Additionally, independent from main test computer, others were utilized, including a laptop with a GeForce GT 730 M GPU processor and a 4-core i7-3630QM CPU (test computer-B), a laptop with a GeForce GTX 950 M graphics processor and a 4-core i7-6700HQ CPU (test computer-C), and another laptop with a GeForce GTX 1050 Ti GPU processor and an i7-7700HQ CPU (test computer-D), each equipped with the specifications listed in Table 2. Each test computer functioned separately, and the performance metrics were collected accordingly.

The selection of benchmark examples was crucial for ensuring a fair comparison between the different hardware configurations. These examples were chosen based on the graphics card with the lowest memory capacity among the various hardware options selected for this study. Additionally, the design of the algorithm was tailored to accommodate this lowest capacity for each hardware configuration, rather than being designed specifically for each one. This approach helped standardize the testing conditions and allowed for a more meaningful comparison of the hardware’s performance in optimizing steel frame structures.

Table 2. Specifications of the test computers.

Name of the test computer	Processor	No. of cores	Core freq.	Double precision perf.	Single precision perf.
A	I7-8700K	6	3700 Mhz	32.11 GFLOPS	61.41 GFLOPS
A	Geforce GTX 1080 Ti	3584	1481 Mhz	11.34 TFLOPS	354.4 GFLOPS
B	I7-7700HQ	4	3800 Mhz	43.07 GFLOPS	33.34 GFLOPS
B	Geforce GTX 1050 Ti	768	1392 MHz	2.138 TFLOPS	66.82 GFLOPS
C	I7-6700HQ	4	3500 Mhz	17.66 GFLOPS	8.39 GFLOPS
C	Geforce GTX 950M	640	1124 MHz	1.439 GFLOPS	44.96 GFLOPS
D	I7-3630QM	4	2400 Mhz	8.83 GFLOPS	27.30 GFLOPS
D	Geforce GT 730M	384	900 MHz	556.8 GFLOPS	23.20 GFLOPS
Ref	Tesla P100	3584	1190 Mhz	4.763 TFLOPS	19.05 TFLOPS

In the study, 8-byte double precision was utilized for both GPU and CPU computations to ensure high precision in the optimization analyses of steel frame structures. With the aim of comparison, the identical analytical procedure was repeated twice on each test computer, employing both CPU and GPU processors. This led to a total of 48 analysis runs, encompassing two distinct case structures – one with 105 members and another with 1024 members. These analyses were conducted on four different computer systems, each subjected to three distinct loading conditions, utilizing both CPU and GPU processors.

Figure 8 and 9 provide time values for all the test computers during analyses performed under 1, 4, and 7 loading conditions.

3.1. Design example-1: five-story, 105-member steel frame

Figure 7. Design example-1. 3D, side, isometric, and plan views of frame structure with 105 members.

The first design example is a five-story steel space frame structure with 105 members, and 54 joints, which is based on previous studies in the literature (Akın and Aydoğdu 2015; I. Aydoğdu 2010; Çarbaş 2016; Örmecioğlu 2019). The 3D, side, and plan views of the frame structure are shown in Figure 7. In the design example-1, the frame members are grouped into 11 independent design variables which are illustrated in Table 3, and the group definitions of the structure were taken from Aydoğdu (2010). The frame members are chosen from a range of 272 W-sections, spanning from W100 × 19.3 to W1100 × 499 mm, as specified in LRFD-AISC (LRFD 2000). In addition to dead and live loads, the structure is exposed to snow and wind loads.

Table 3. Group definition of the design example-1:The 105-member frame structure.

Story	Side beam	Inner beam	Corner column	Side column	Inner column
5	1	2	3	3	5
4	1	2	3	3	5
3	1	2	6	7	8
2	1	2	6	7	8
1	1	2	9	10	11

In an instance comprising 105 elements, optimization was carried out with CPU and GPU, employing the load combinations provided above and the computers specified in Table. This process was repeated a total of 2 times, and the cumulative duration of results for an average iteration is displayed in Figure 8 Furthermore, Figure 9 depicts the distribution of CPU and GPU core utilization and memory during this iteration.

Figure 8. GPU/CPU kernel times of the 105-member structure under 1,4 and 7 load cases. LC: the number of load cases.

Figure 9. Process/Memory, time percent graph of internal GPU kernel for the 105-member structure under 1, 4, and 7 load cases. LC: the number of load cases.

In Figure 8, it is evident that the GPU processing time is significantly higher than that of the CPU, ranging from 2.7 to 70 times faster. The discrepancy in GPU-CPU processing time is most pronounced in the “D” computer configuration, with an average 49-fold difference, and least pronounced in the “A” configuration, with an average fourfold difference. As the number of load combinations increases, both CPU and GPU solution times on the sample also increase.

In Figure 9, the percentage of time spent on memory operations falls within the range of 65% to 80%. Notably, in the case of single-loading scenarios, there exists a more inconsistent distribution between memory and processing times, whereas, under conditions of 4 and 7 loading combinations, the relationship between memory and processing times appears to be more closely aligned.

Within Table 4, the minimum weight of the 105-member structure, as computed by our developed algorithm, is contrasted against findings from existing literature. A solution is yielded by our BBO-GPU algorithm that is %4.79 heavier than the lightest weight, %40.6 lighter than the heaviest weight, and %6.36 lighter than the average weight. These deviations suggest that the BBO-GPU solution falls within an acceptable range.

Table 4. Design sections and limit values of the optimum designs for the 105-member space frame.

Member group	Type	ACO (I. Aydoğdu 2010)	HS (I. Aydoğdu 2010)	hTLBO-HS (Akın and Aydoğdu 2015)	SSO	SJ_SSSO	EFF (Çarbaş 2016)	CS (Çarbaş 2016)	PSO (Çarbaş 2016)	FF (Çarbaş 2016)	BBO-GPU
1	Beam	W460X52	W530X66	W460X52	W460X52	W460X52	W410X46.1	W460X52	W460X52	W310X52	W460X52
2	Beam	W200X35.9	W310X38.7	W200X35.9	W360X44	W200X35.9	W200X41.7	W250X44.8	W200X41.7	W410X53	W200X35.9
3	Column	W200X35.9	W200X35.9	W310X38.7	W310X38.7	W310X38.7	W360X44	W610X82	W310X67	W360X101	W250X49.1
4	Column	W310X38.7	W200X35.9	W200X46.1	W200X41.7	W200X35.9	W250X49.1	W200X35.9	W460X113	W610X101	W200X35.9
5	Column	W360X57.8	W360X44	W360X44	W460X52	W410X53	W410X60	W610X125	W760X161	W310X117	W530X123
6	Column	W460X52	W310X38.7	W310X74	W310X38.7	W410X53	W410X60	W610X82	W310X67	W530X101	W250X49.1
7	Column	W310X86	W360X72	W250X73	W360X72	W250X73	W410X100	W250X67	W460X113	W760X134	W530X92
8	Column	W610X101	W610X92	W610X101	W760X134	W690X125	W530X123	W760X185	W760X161	W690X125	W690X152
9	Column	W530X66	W410X53	W460X74	W410X53	W410X60	W530X74	W690X170	W310X67	W690X140	W250X49.1
10	Column	W460X89	W360X72	W250X73	W460X74	W310X74	W460X106	W760X147	W460X113	W760X134	W610X101
11	Column	W690X170	W760X147	W760X173	W840X176	W760X134	W840X176	W920X488	W760X161	W690X170	W760X173
Max. Str. Ratio		0.886	0.979	0.921	0.964	0.994	0.946	0.842	0.831	0.873	0.872
TSD (cm)		4.983	4.837	4.708	4.945	5.083	4.822	3.95	3.663	4.588	4.859
ISD (cm)		0.569	1.333	1.325	1.332	1.326	1.332	1.139	1.075	1.333	1.314
Max. Iter.		50000	50000	50000	50000	50000	50000	50000	50000	50000	75000
Weight (kN)		265.38	278.196	269.184	270.976	258.416	284.04	336.05	338.38	391.06	278.1

Figure 10. Search history of 105-member space frame.

Figure 10 shows the search history behavior of the BBO-GPU concerning the other optimization algorithm previously used in the 105-member example. This figure reveals that the BBO-GPU’s convergence remains in the interval of the optimization algorithms.

3.2. Design example-2: twenty-story, 460-member steel frame

Figure 11. Design example-2. 3D, side, isometric, and plan views of frame structure with 460 members.

The second design example involves a twenty-story steel space frame structure comprising 460 members and 210 joints. This design is based on previous research outlined in Aydoğdu (2010), and the group definitions of the structure were taken from Aydoğdu (2010) article; they were conducted as depicted in Table 5. Figure 11 displays 3D, side, isometric, and plan views of the frame structure. For this particular design example, we have categorized the frame members into 13 distinct design variables, which are listed in Table 5. These frame members have been selected from a range of 272 W-sections, with sizes ranging from W100 × 19.3 to W1100 × 499 mm, following LRFD-AISC guidelines (LRFD 2000).

Table 5. Group definition of the design example-2:The 460-member frame structure. See Figure 10 for a-, b-, c-, d, and e-beam locations.

Story	Beams					Corner column	Side column	Inner column
	(a)	(b)	(c)	(d)	(e)
20–18	1	2	3	4	5	6	6	6
17–15	1	2	3	4	5	7	7	7
14–12	1	2	3	4	5	8	8	8
11–9	1	2	3	4	5	9	9	9
8–7	1	2	3	4	5	10	10	10
6–5	1	2	3	4	5	11	11	11
4–3	1	2	3	4	5	12	12	12
2–1	1	2	3	4	5	13	13	13

Optimization was conducted on a system comprising 460 elements, utilizing both CPU and GPU processing. This optimization process incorporated the load combinations detailed earlier and was executed on the computer systems specified in Table. This procedure was iterated twice, and Figure 12 provides an overview of the average duration of these iterations. Additionally, Figure 13 illustrates the utilization of CPU and GPU cores, as well as memory distribution throughout this iterative process.

In Figure 12, we observe that the GPU processing time outperformed the CPU in only 4 out of 12 combinations, achieving a speedup factor of 1.1–1.4 times. Conversely, in the remaining 8 scenarios, the GPU lagged behind the CPU, with processing times ranging from 1.1 to 4.0 times slower. The most significant discrepancy in GPU-CPU processing times was observed in computer “A,” with an average improvement of 1.2 times, while the smallest difference occurred in computer “C,” where the GPU was 3.4 times slower on average. Notably, as the number of load combinations increased, the CPU and GPU solution times also exhibited an upward trend. Across computers A, B, and C, the average time differences amounted to 6.5 milliseconds, whereas computer D showed a considerably larger average difference of 395.5 milliseconds.

Figure 12. GPU/CPU kernel times of design example-2:The 460-member structure under 1,4 and 7 load cases.

Figure 13 illustrates that the time allocated to memory operations ranges from 65% to 93% of the total processing time. Remarkably, computer “D” consistently registered the highest processing times for each load type, while computer “A” consistently demonstrated the shortest processing times. Furthermore, it is noteworthy that, for computer “A,” the proportion of memory time relative to the overall time remained above 90% across all three load types.

Figure 13. Process/Memory, Time Percent Graph of Internal GPU kernel for the 460-member structure under 1, 4, and 7 load cases.

Previously, the optimization of the 460-member structure was conducted using the Ant Colony Optimization (ACO) algorithm (İ. Aydoğdu and Saka 2012). When comparing the outcomes of BBO-GPU against those of ACO, a 0,35% difference becomes apparent between the two sets of results (Table 6). This variation remains within an acceptable range for evaluation.

Table 6. Design sections and limit values of the optimum designs for the 460-member space frame.

Group number	Group type	ACO (İ. Aydoğdu and Saka 2012)	BBO-GPU
1	Beam	W610X92	W840X176
2	Beam	W310X28.3	W360X91
3	Beam	W760X196	W610X101
4	Column	W460X68	W610X92
5	Column	W530X66	W360X91
6	Column	W310X202	W840X176
7	Column	W360X237	W840X176
8	Column	W360X237	W840X176
9	Column	W610X262	W840X176
10	Column	W760X314	W840X176
11	Column	W760X314	W840X176
12	Column	W760X314	W840X210
13	Column	W840X329	W840X226
Minimum weight (kN)		3191.15	3180
Maximum top story drift (cm)		17.9	18.27
Maximum inter- story drift (cm)		1.133	1.18
Maximum strength constraint ratio		0.983	0.996
Maximum number of Iterations		75000	75000

Figure 14. Search history of 460-member space frame.

Figure 14 presents a comparison of the search history patterns between the BBO-GPU and the optimization algorithm that was used earlier in the 460-member example. The data from this figure indicates that the BBO-GPU’s convergence performance falls within the boundaries set by the earlier optimization algorithms.

3.3. Design example-3: eight-story, 1024-member steel frame

The third design example is an eight-story steel space frame structure with 460 members and 384 joints in 40 independent design groups, taken from previous studies (I. Aydoğdu and Akın 2014; İ. Aydoğdu and Saka 2012; Örmecioğlu 2019). Figure 15 displays the isometric, side, and plan views for design example 3, while Table 7 outlines the grouping of frame members. The group definitions of the structure were taken from Aydoğdu and Saka (2012). Like the design example-1, the frame members are chosen from a selection of 272 W-sections, ranging from W100 × 19.3 to W1100 × 499 mm, as outlined in LRFD-AISC (LRFD 2000).

Table 7. Group definition of the design example-3:The 1024-member frame structure.

Story	Side beam	Inner beam	Corner column	Side column	Inner column
8	1	2	17	18	19
7	3	4	20	21	22
6	5	6	23	24	25
5	7	8	26	27	28
4	9	10	29	30	31
3	11	12	32	33	34
2	13	14	35	36	37
1	15	16	38	39	40

Figure 15. Design example-3. 3D, side, isometric, and plan views of frame structure with 1024 members.

In the context of a scenario involving 1024 elements, optimization was conducted using both CPU and GPU, employing the load combinations previously outlined and the computer systems specified in Table. This optimization process was replicated twice, as depicted in Figure 16, which presents the cumulative duration of results for an average iteration. Additionally, Figure 17 offers a visual representation of the distribution of CPU and GPU core utilization, as well as memory usage, throughout these optimization iterations. This graphical representation provides a comprehensive insight into the allocation of computational resources during the optimization process.

Figure 16 presents a clear comparison between GPU and CPU processing times across 12 distinct combinations. Notably, in 7 of these scenarios, the GPU outperformed the CPU, completing tasks 1.1 to 3.8 times faster. Conversely, in the remaining 5 scenarios, the GPU exhibited slower processing times, ranging from 1.0 to 1.7 times longer than the CPU. The largest discrepancy between GPU and CPU processing time was observed in the “A” computer configuration, with an average reduction of 2.6 times, while the “C” configuration showed the least significant difference, with an average increase of 1.1 times. It is important to note that as the number of load combinations increased, both CPU and GPU solution times demonstrated a gradual upward trend in the sample data.

Figure 16. GPU/CPU kernel times of the design example-3:The 1024-member structure under 1,4 and 7 load cases.

Figure 17 presents data indicating that the percentage of time dedicated to memory operations lies in the range of 53% to 74%. It is noteworthy that, across all loading types, the “D” computer configuration consistently records the highest processing time, whereas the “A” configuration consistently records the lowest processing time.

Figure 17. Process/Memory, time percent graph of internal GPU kernel for the 1024-member structure under 1, 4, and 7 load cases.

The minimum weight achieved using the developed algorithm is compared to the results presented in the literature (I. Aydoğdu et al. 2017; Çarbaş 2016), as shown in Table 8. Upon comparison, a difference of 0,44 % is observed between our findings and those of the referenced studies. Given the stochastic nature of metaheuristic algorithms, such a variance can be deemed acceptable.

Table 8. Design details of the best solutions for the 1024-member structure.

#	Group Type	BBO	BBO-GPU	#	Group Type	BBO	BBO-GPU
1	Beam	W310X21	W310X38.7	21	Column	W360X196	W460X128
2	Beam	W310X21	W200X46.1	22	Column	W250X58	W360X134
3	Beam	W410X53	W310X79	23	Column	W360X162	W360X196
4	Beam	W250X49.1	W310X117	24	Column	W360X237	W1000X258
5	Beam	W610X101	W460X82	25	Column	W250X67	W360X147
6	Beam	W250X58	W460X82	26	Column	W360X162	W690X240
7	Beam	W410X100	W360X110	27	Column	W360X347	W1000X321
8	Beam	W250X80	W360X110	28	Column	W310X97	W410X149
9	Beam	W610X155	W610X125	29	Column	W610X174	W760X284
10	Beam	W530X165	W460X52	30	Column	W1100X390	W1000X477
11	Beam	W690X125	W610X82	31	Column	W530X150	W530X165
12	Beam	W310X107	W460X68	32	Column	W760X314	W920X417
13	Beam	W610X125	W840X176	33	Column	W1100X499	W1000X477
14	Beam	W460X52	W310X32.7	34	Column	W530X150	W840X176
15	Beam	W840X176	W760X147	35	Column	W920X446	W920X417
16	Beam	W310X28.3	W310X32.7	36	Column	W1100X499	W1100X499
17	Column	W310X107	W200X35.9	37	Column	W920X201	W840X176
18	Column	W360X196	W250X80	38	Column	W1000X443	W1100X433
19	Column	W200X41.7	W310X60	39	Column	W1100X499	W1100X499
20	Column	W310X117	W310X158	40	Column	W1000X258	W1000X222
BBO				BBO-GPU
Max. Strength Ratio			1	Max. Strength Ratio			1
Top Story Drift (cm)			6.508	Top Story Drift (cm)			6.860
Inter Story Drift (cm)			0.875	Inter Story Drift (cm)			0.875
Maximum Iteration			75000	Maximum Iteration			75000
Weight (kN)			6462.79	Weight (kN)			6491.28

Figure 18. Search history of 1024-member space frame.

Figure 18 illustrates the comparative search history behavior of the BBO-GPU against the previously utilized optimization algorithm in the context of the 1024-member example. It is evident from this figure that the convergence of the BBO-GPU is maintained within the range defined by the other optimization algorithms.

4. Discussion

When analyzing the performance metrics of three design examples, a key distinction arises where the GPU processor struggles to outperform the CPU, especially in scenarios involving relatively few elements within the framework systems. However, a notable performance improvement is observed when analyzing a structure with 1024 members on “test computer-A,” where the GPU processor proves to be over three times faster than the CPU. To comprehend the reasons behind this phenomenon, a software profile analysis is conducted. This analysis involves capturing a comprehensive record of events and their durations during computations on the GPU processor.

As depicted in Figures 9, 13 , and 17, nearly 60% of the total time consumed by the GPU is attributed to data transfer operations between the CPU (main computer) and the GPU processor (device). The remaining 40% represents the actual time spent by the GPU processor on mathematical computations. Comparing processing times without memory-related tasks to CPU processing time reveals that a performance increase of up to 9.5 times can be achieved, particularly evident in the analysis of the 1024-member structure in the case of “test computer-A,” while it has less impact on the 105-member structure. Additionally, it is observed that the 460-member structure achieves similar results across all conditions except for computer-D.

When examining the analyses conducted within this study, the results can be categorized into two main performance groups. The first classification is primarily based on architectural differences between CPU and GPU processor groups. The most significant difference is the scalability of performance with the amount of processed data. This becomes particularly evident when comparing the analysis of a 105-member framework, which requires fewer computations and less memory, to that of a 1024-member structure, which is more computation-intensive and memory-demanding. Furthermore, the 460-member structure tends to evenly allocate processor and memory consumption, positioning it between the results of the other two structure examples. Figure 8 demonstrates that under all load conditions, classical CPU processors outperform GPU processors in the analysis of the 105-member structure. However, as shown in Figure 16, the CPU processor falls behind the GPU processor in the analysis of the 1024-member structure. Figure 12 indicates that CPU and GPU processor utilization rates are similar. Therefore, in data-intensive and repetitive operations, the GPU’s parallel processing capability tends to yield significantly better results.

The performance gap between GPU and CPU analyses can be attributed to two key factors. Firstly, differences in hardware configurations contribute significantly to this gap. Secondly, a bottleneck occurs in memory transfer between the host and the device. The study highlights the impact of GPU processor architecture, as different GPU models show varying performance levels due to differences in their command sets. GPUs with advanced command sets demonstrate superior parallel computational capabilities compared to those with older architectures, underscoring the importance of GPU model selection for achieving optimal performance in parallel programming for structural optimization. The hardware selection deliberately included devices from various release years, enhancing the diversity of GPU architectures tested. This variation in processing capabilities among the hardware led to differences in speed and performance. Moreover, a fixed portion of time is required for handshake, initialization, and data transfer during computation, which is necessary for GPU calculations but not for CPU calculations. Consequently, it is not possible to reduce the computation time below a certain threshold when using GPUs. To illustrate this, the study includes a specific example with minimal computation and memory transfer involving 105 elements, presenting this phenomenon through tables and graphs. Figures 8,12 and illustrate the percentage of the memory bottleneck, which is the reason for the occurrence of situations where the GPU spends more time than the CPU in the analysis of structures with a low number of members.

The second classification pertains to different GPU processor models, each characterized by unique processor architectures, based on their production years. As indicated in Table 9, variations in performance among GPU processors stem not only from differences in data pathways but also from disparities in processor architectures, particularly in the command sets. In the context of the SM 3.0 architecture, GPU processors in test machines other than “D” exhibit superior capabilities for parallel mathematical computations due to their more advanced command sets. The efficiency of these commands becomes evident in the obtained results. Conversely, the GPU processor in the “D” computer requires the use of older software solutions due to the absence of specific commands compared to others. Importantly, within the SM 3.0 architecture, no command set allows atomic addition for parallel floating-point calculations, unlike the processors used in other test machines. As a result, it has been observed that commands in the “D” test computer operate significantly slower than others due to this limitation.

The fluctuations In processing times observed In the single-load scenario shown in Figure may be attributed to the necessity of data transfer between the main computer and the GPU during GPU processing. These fluctuations serve as an indicator that using a GPU may not provide significant advantages in scenarios characterized by low computational requirements. Additionally, the similarity in processor/memory usage ratios for the 460-member example in Figure 13 highlights the threshold at which GPU usage may only be advantageous beyond a certain level of computational demand.

Table 9. GPU processor architectures with their production years.

Name of the test computer	Manufacture year	SM ver.(Major.Minor)	GPU architecture name	GPU model name	Busline (GB/sn)
A	2017	SM 6.1	Pascal	1080 titan	484.4
B	2016	SM 6.1	Pascal	1050 titan	112.1
C	2015	SM 5.0	Maxwell	GTX 950M	28.8
D	2013	SM 3.0	Kepler	GT 730M	28.8

To validate the optimality of the developed algorithm, we compared its optimum results with those found in existing literature. The comparison reveals a deviation of 6,36% for the first example and 0,35% for the second and third examples. These observed differences comfortably fall within acceptable ranges of evaluation, affirming the algorithm’s effectiveness.

5. Conclusion

The study explores the application of GPU-based parallel programming to optimize large-scale steel frame structures, addressing the computational complexities associated with iterative structural analysis processes. The research demonstrates that GPU processors provide a significant performance advantage over CPUs. The extent of this advantage depends on several factors, including processor architecture, computational demands, and the specific GPU model used. Notably, in scenarios with high computational requirements, such as analyzing a 1024-member frame structure, GPUs outperform CPUs by a factor of more than three. This highlights the effectiveness of GPU parallel processing in managing time-intensive optimization challenges in structural engineering.

Furthermore, the research emphasizes the importance of data transfer operations between the CPU and GPU, especially in situations with lower computational demands. As computational requirements increase, the benefits of GPU parallel processing become more evident. This underscores the necessity of carefully considering computational demands and GPU architecture when implementing parallel computing methods.

The study also highlights the impact of GPU processor architecture, as different GPU models exhibit varying levels of performance due to differences in their command sets. GPUs with advanced command sets demonstrate superior parallel computational capabilities compared to those with older architectures, emphasizing the role of GPU model selection in achieving optimal performance in parallel programming for structural optimization.

In conclusion, this research emphasizes the promising role of GPU-based parallel programming in accelerating the optimization of large-scale steel frame structures. It delivers significant improvements in time efficiency and performance, especially in scenarios involving complex and computationally intensive iterative processes. While this study focuses on the possibility of expediting the optimization of steel structures using ordinary devices, further research is warranted to fully explore the potential of GPU-based solutions in civil engineering and related fields. Specifically, future studies should focus on evaluating the performance of GPU-based methods in large-scale finite element analysis to gain a deeper understanding of their capabilities and limitations. Nevertheless, it is important to note that increasing the number of members in a group may lead to convergence delays. Therefore, future studies should consider the use of guide-based techniques for larger member numbers to address this omission and further enhance the optimization process for steel frame structures.

Declaration of generative AI and AI-assisted technologies in the writing process

During the preparation of this work, the authors used ChatGPT to improve readability and language. After using this tool/service, the authors reviewed and edited the content as needed and take(s) full responsibility for the content of the publication.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Additional information

Notes on contributors

Tevfik Oğuz Örmecioğlu

MSc Örmecioğlu is a PhD Candidate and a civil engineer with extensive experience in construction project management and infrastructure development. He has a strong academic background with a focus on structural optimization and artificial intelligence. He is also proficient in various programming languages and software applications, making him a versatile and skilled professional in his field.

İbrahim Aydoğdu

Dr. Aydoğdu is a faculty member at Akdeniz University, where he teaches courses on artificial intelligence, machine learning, and deep learning. His passion lies in sharing his knowledge and expertise with students, inspiring them to explore the exciting world of AI. Alongside teaching, he conducts research and publishes in the areas of natural inspired optimization algorithms and computational mechanics.

Hilal Tuğba Örmecioğlu

Dr. Örmecioğlu is an academic specializing in architecture, with a focus on the history of Modern Architecture, Construction, and Technology. Her interests cover earthquake-resistant architectural design to the preservation of modern architecture. Her research also extends to the integration of AI into architectural education and the utilization of new technologies in architectural and structural design.